%I #43 Jun 29 2022 02:18:34
%S 0,4,8,8,12,20,20,20,24,28,36,36,36,44,44,44,48,56,60,60,68,68,68,68,
%T 68,80,88,88,88,96,96,96,100,100,108,108,112,120,120,120,128,136,136,
%U 136,136,144,144,144,144,148,160,160,168,176,176,176,176,176,184,184
%N Number of integer solutions to x^2 + y^2 <= n excluding (0,0).
%C a(32)/32 = 100/32 = 3.125; lim_{n->infinity} a(n)/n = Pi.
%C The terms of this sequence are four times the running total of the excess of the 4k + 1 divisors of the natural numbers (from 1 through to n) over their 4k + 3 divisors. - _Ant King_, Mar 12 2013
%D G. H. Hardy and E. M. Wright, An Introduction to the Theory of Numbers, 5th ed., Oxford Univ. Press, 1979, th. 339
%H Seiichi Manyama, <a href="/A014198/b014198.txt">Table of n, a(n) for n = 0..10000</a>
%H Eric Weisstein's World of Mathematics, <a href="http://mathworld.wolfram.com/SumofSquaresFunction.html">Sum of Squares Function</a>
%F a(n) = 4*A014200(n).
%F a(n) = A057655(n)-1.
%e For n=2 the 8 solutions are (x,y) = (+-1,0), (0,+-1), (+-1,+-1).
%p A014198 := proc(n)
%p nops([ numtheory[thue]( abs( x^2+y^2) <= n, [ x, y ] ) ]);
%p end proc:
%p seq(A014198(n),n=0..60) ;
%t Prepend[SquaresR[2,#] &/@Range[59],0]//Accumulate (* _Ant King_, Mar 12 2013 *)
%o (PARI) a(n)=local(j); j=sqrtint(n); sum(x=-j,j,sum(y=-j,j,x^2+y^2<=n))-1
%o (Python)
%o from math import prod
%o from itertools import count, accumulate, islice
%o from sympy import factorint
%o def A014198_gen(): # generator of terms
%o return accumulate(map(lambda n:prod(e+1 if p & 3 == 1 else (e+1) & 1 for p, e in factorint(n).items() if p > 2) << 2, count(1)),initial=0)
%o A014198_list = list(islice(A014198_gen(),30)) # _Chai Wah Wu_, Jun 28 2022
%Y Cf. A014200, A057655.
%K nonn,easy
%O 0,2
%A _N. J. A. Sloane_