Draw Number Spirals, Part 2

Question

In relation to my post from yesterday:

Draw number spirals

I like to find which of the numbers in the OUTER LAYER of the spiral are Prime numbers. For the range of n=1...1000

In doing so, I found the pattern that the outer layer numbers in terms of n are always in this range: between

1+(n-2)^2 and n^2

ex: n=4 [5-->16] n=5 [10-->25] , etc

Now I like to see which n value generates the MOST primes. I came up with this:

[Expression-1]

Select[Range[1000],Total@Boole@Table[PrimeQ[Range[1 + (# - 2)^2, #^2]]] > 300&]

Which gives me 954

and we can see it is the most and list out those Primes - in fact there are 308 Primes:

[Expression-2]

Intersection[Prime[Range[954^2]], Range[1 + (954 - 2)^2, 954^2]]

However, 2 things I need:

1. Above had a trial and error - I had to tweak that constant from various values. Is there a way to modify Expression-1 so it finds the MAX value wihtout having to specify the criteria >300 ?

2. Can we generalize the expression to find all Prime counts for ALL layers (1...1000) and SORT them descending order with counts?

ex:

954==> 308 primes

999==> 299 primes

962==> 298 primes

.....

(those are real values!)

Answer 1

Define a function that counts the number of primes in the range (using PrimePi), and ReverseSortBy that function (note I change the lower bound to be 1 less so that 1 + (n-2)^2 is included if it's also prime:

primeCount := PrimePi[#^2] - PrimePi[(# - 2)^2] &

primeData = {#, primeCount[#]} & /@ Range[1000];
sorted = ReverseSortBy[primeData, Last];

So the 10 largest are:

sorted[[1 ;; 10]] // MatrixForm

To get the n that maximizes PrimeCount we can just use MaximalBy:

maximalPrimeN := First@MaximalBy[Range@#, primeCount] &

so the maximum less than 1000 is 954:

maximalPrimeN[1000]
(*954*)

FWIW we can also define maximalPrimeN recursively, but we have to make the $RecursionLimit infinite (I really don't like how I defined this, I feel like there is a cleaner way to do this):

maximalPrimeNRecursive[1] = 1;
maximalPrimeNRecursive[n_] := maximalPrimeNRecursive[n] =
  Block[{$RecursionLimit = Infinity},
   
   Module[{ maxN, maxCount },
    
    maxN = maximalPrimeNRecursive[n - 1];
    maxCount = primeCount@maxN;
    
    If[primeCount[n] > maxCount,
     n,
     maxN
     ]
    ]
   ]