%I #20 Sep 08 2022 08:45:11
%S 41,809,1321,2729,4649,5801,11689,15401,17449,21929,26921,41641,52009,
%T 55721,59561,71849,80681,94889,99881,126761,156841,169769,190121,
%U 197161,204329,226601,234281,266281,327209,345769,394409,457001,467881,524201
%N Primes of the form 16*m^2 + 25, m=1,3,5,...
%C This is a special case of the theorem that all prime numbers of the form 4k+1 can be expressed as the sum of two squares. Let p = a^2 + b^2 then a=4k+1 and b = 4m. From this it follows that p = 16(m^2 + k^2) + 8k + 1. When k=1 we have p = 16m^2 + 25. If we let j=16m then the arithmetic progression j*m + 25 has an infinite number of primes by Dirichlet's theorem.
%C Primes of the form 64*k^2 + 64*k + 41. - _Vincenzo Librandi_, Dec 11 2011
%D H. Rademacher, Lectures on Elementary Number Theory, 1964, pp. 121-136.
%H Vincenzo Librandi, <a href="/A087856/b087856.txt">Table of n, a(n) for n = 1..5000</a>
%t Select[Table[16m^2 + 25, {m, 1, 201, 2}],PrimeQ] (* _Harvey P. Dale_, Jan 24 2011 *)
%t Select[Table[64n^2+64n+41,{n,0,4000}],PrimeQ] (* _Vincenzo Librandi_, Dec 11 2011 *)
%o (PARI) fourmp1(m,k=1) = { forstep(x=1,m,2, y=16*(x^2+k^2)+8*k+1; if(isprime(y),print1(y", ")) ) }
%o (Magma) [a: n in [0..100] | IsPrime(a) where a is 64*n^2+64*n+41]; // _Vincenzo Librandi_, Dec 11 2011
%o (GAP) Filtered(List([1,3..201],m->16*m^2+25),IsPrime); # _Muniru A Asiru_, Nov 24 2018
%Y Cf. A087857, A087861, A087862.
%K nonn,easy
%O 1,1
%A _Cino Hilliard_, Oct 09 2003