|
%I #17 Sep 09 2022 04:24:31
%S 0,4,13,13,38,38,87,87,87,87,208,208,377,377,377,377,666,666,1027,
%T 1027,1027,1027,1556,1556,1556,1556,1556,1556,2397,2397,3358,3358,
%U 3358,3358,3358,3358,4727,4727,4727,4727,6408,6408,8257,8257,8257,8257,10466,10466
%N a(n) = Sum_{2 <= p <= n, p prime} p^2.
%H Daniel Suteu, <a href="/A081738/b081738.txt">Table of n, a(n) for n = 1..10000</a>
%t Table[Total[Prime[Range[PrimePi[n]]]^2],{n,48}] (* _Stefano Spezia_, Aug 22 2022 *)
%o (PARI) a(n, j=2) = if(n <= 1, return(0)); my(r=sqrtint(n)); my(V=vector(r, k, n\k)); my(F(n,j)=(subst(bernpol(j+1),x,n+1) - subst(bernpol(j+1),x,1)) / (j+1)); my(L=n\r-1); V=concat(V, vector(L, k, L-k+1)); my(T=vector(#V, k, F(V[k],j))); my(S=Map(matrix(#V,2,x,y,if(y==1,V[x],T[x])))); forprime(p=2, r, my(sp=mapget(S,p-1), p2=p*p); for(k=1, #V, if(V[k] < p2, break); mapput(S, V[k], mapget(S,V[k]) - p^j*(mapget(S,V[k]\p) - sp)))); mapget(S,n)-1; \\ _Daniel Suteu_, Aug 21 2022
%o (PARI) a(n) = norml2(primes(primepi(n))); \\ _Michel Marcus_, Aug 22 2022
%Y Cf. A024450, A024525, A133391.
%K nonn
%O 1,2
%A _N. J. A. Sloane_, Apr 07 2003
|
