cp's OEIS Frontend

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

%I A081738 #20 Jan 31 2025 03:30:31
%S A081738 0,4,13,13,38,38,87,87,87,87,208,208,377,377,377,377,666,666,1027,
%T A081738 1027,1027,1027,1556,1556,1556,1556,1556,1556,2397,2397,3358,3358,
%U A081738 3358,3358,3358,3358,4727,4727,4727,4727,6408,6408,8257,8257,8257,8257,10466,10466
%N A081738 a(n) = Sum_{2 <= p <= n, p prime} p^2.
%H A081738 Daniel Suteu, <a href="/A081738/b081738.txt">Table of n, a(n) for n = 1..10000</a>
%t A081738 Table[Total[Prime[Range[PrimePi[n]]]^2],{n,48}] (* _Stefano Spezia_, Aug 22 2022 *)
%o A081738 (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 A081738 (PARI) a(n) = norml2(primes(primepi(n))); \\ _Michel Marcus_, Aug 22 2022
%o A081738 (Magma)
%o A081738 A081738:= func< n | n eq 1 select 0 else (&+[k^2: k in PrimesInInterval(1, n)]) >;
%o A081738 [A081738(n): n in [1..60]]; // _G. C. Greubel_, Jan 31 2025
%o A081738 (Python)
%o A081738 from sympy import prime, primerange
%o A081738 def A081738(n): return sum(p**2 for p in primerange(2,n+1))
%o A081738 print([A081738(n) for n in range(1,61)]) # _G. C. Greubel_, Jan 31 2025
%Y A081738 Cf. A024450, A024525, A133391.
%K A081738 nonn
%O A081738 1,2
%A A081738 _N. J. A. Sloane_, Apr 07 2003