%I #40 Feb 24 2017 02:57:09
%S 17,29,37,73,89,97,109,149,157,193,229,241,269,277,349,409,433,541,
%T 601,661,709,769,829,853,929,937,1009,1021,1069,1109,1117,1129,1249,
%U 1321,1409,1429,1489,1549,1609,1669,1753,1789,1801,1873,2029,2089,2161,2221
%N Primes which are half the sum of 2 squares of primes.
%C Primes of the form x^2 + y^2, where x > y > 0, such that x-y = p and x+y = q are primes. Proof: (p^2+q^2)/2 = ((x-y)^2+(x+y)^2)/2 = x^2+y^2 so we have x = (p+q)/2 and y = (q-p)/2. - _Thomas Ordowski_, Sep 24 2012
%C All terms == 1 or 5 (mod 12). - _Thomas Ordowski_, Jun 28 2013
%C Or, primes in A143850. - _Zak Seidov_, Jun 06 2015
%H T. D. Noe, <a href="/A103739/b103739.txt">Table of n, a(n) for n = 1..1000</a>
%e 17 is in the sequence because (3^2 + 5^2) / 2 = 17.
%p Primes:= select(isprime,[seq(2*i+1,i=1..400)]):
%p Psq:= map(`^`,Primes,2):
%p M:= max(Psq):
%p S:= select(t -> t <= M/2 and isprime(t),{seq(seq((Psq[i]+Psq[j])/2, j=1..i-1),i=1..nops(Psq))}):
%p sort(convert(S,list)); # _Robert Israel_, Jun 08 2015
%o (PARI) list(lim)=my(v=List(), p2, t); lim\=1; if(lim<9, lim=9); forprime(p=3, sqrtint(2*lim-9), p2=p^2; forprime(q=3, min(sqrtint(2*lim-p2), p), if(isprime(t=(p2+q^2)/2), listput(v, t)))); Set(v) \\ _Charles R Greathouse IV_, Feb 14 2017
%Y Intersection of A143850 and A000040.
%Y Cf. A001248, A002313.
%K easy,nonn
%O 1,1
%A _Giovanni Teofilatto_, Mar 28 2005
%E Corrected and extended by _Walter Nissen_, Jul 19 2005