A081738 a(n) = Sum_{2 <= p <= n, p prime} p^2.

0, 4, 13, 13, 38, 38, 87, 87, 87, 87, 208, 208, 377, 377, 377, 377, 666, 666, 1027, 1027, 1027, 1027, 1556, 1556, 1556, 1556, 1556, 1556, 2397, 2397, 3358, 3358, 3358, 3358, 3358, 3358, 4727, 4727, 4727, 4727, 6408, 6408, 8257, 8257, 8257, 8257, 10466, 10466

Offset: 1

Author

N. J. A. Sloane, Apr 07 2003

Keywords

• nonn

Links

• Daniel Suteu, Table of n, a(n) for n = 1..10000

Crossrefs

Cf. A024450, A024525, A133391.

Programs

Magma

A081738:= func< n | n eq 1 select 0 else (&+[k^2: k in PrimesInInterval(1, n)]) >;
[A081738(n): n in [1..60]]; // G. C. Greubel, Jan 31 2025

Mathematica

Table[Total[Prime[Range[PrimePi[n]]]^2],{n,48}] (* Stefano Spezia, Aug 22 2022 *)

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

PARI

a(n) = norml2(primes(primepi(n))); \\ Michel Marcus, Aug 22 2022

Python

from sympy import prime, primerange
def A081738(n): return sum(p**2 for p in primerange(2,n+1))
print([A081738(n) for n in range(1,61)]) # G. C. Greubel, Jan 31 2025