%I #44 Apr 19 2019 18:35:24
%S 0,0,0,0,1,0,0,1,0,0,1,0,1,1,0,1,1,0,1,1,1,1,1,1,0,1,1,1,3,0,1,2,0,1,
%T 2,1,2,2,1,2,1,0,3,2,1,2,1,2,1,2,2,1,4,1,2,3,0,2,4,1,3,2,1,4,1,1,3,3,
%U 2,2,4,2,1,3,2,3,4,2,3,3,1,2,5,2,4,3,2,4,1,1,6,4,3,4,2,3,0,4,4,3,5,1,5,5,1,4,5,2
%N Number of partitions of n into 4 nonzero squares.
%C Records occur at n= 4, 28, 52, 82, 90, 130, 162, 198, 202, 210,.... - _R. J. Mathar_, Sep 15 2015
%H T. D. Noe, <a href="/A025428/b025428.txt">Table of n, a(n) for n = 0..10000</a>
%H <a href="/index/Su#ssq">Index entries for sequences related to sums of squares</a>
%F For n>0, a(n) = ( A063730(n) + 6*A213024(n) + 3*A063725(n/2) + 8*A092573(n) + 6*A010052(n/4) ) / 24. - _Max Alekseyev_, Sep 30 2012
%F a(n) = ( A000118(n) - 4*A005875(n) - 6*A004018(n) - 12*A000122(n) - 15*A000007(n) + 12*A014455(n) - 24*A033715(n) - 12*A000122(n/2) + 12*A004018(n/2) + 32*A033716(n) - 32*A000122(n/3) + 48*A000122(n/4) ) / 384. - _Max Alekseyev_, Sep 30 2012
%F a(n) = [x^n y^4] Product_{k>=1} 1/(1 - y*x^(k^2)). - _Ilya Gutkovskiy_, Apr 19 2019
%F a(n) = Sum_{k=1..floor(n/4)} Sum_{j=k..floor((n-k)/3)} Sum_{i=j..floor((n-j-k)/2)} A010052(i) * A010052(j) * A010052(k) * A010052(n-i-j-k). - _Wesley Ivan Hurt_, Apr 19 2019
%p A025428 := proc(n)
%p local a,i,j,k,lsq ;
%p a := 0 ;
%p for i from 1 do
%p if 4*i^2 > n then
%p return a;
%p end if;
%p for j from i do
%p if i^2+3*j^2 > n then
%p break;
%p end if;
%p for k from j do
%p if i^2+j^2+2*k^2 > n then
%p break;
%p end if;
%p lsq := n-i^2-j^2-k^2 ;
%p if lsq >= k^2 and issqr(lsq) then
%p a := a+1 ;
%p end if;
%p end do:
%p end do:
%p end do:
%p end proc:
%p seq(A025428(n),n=1..40) ; # _R. J. Mathar_, Jun 15 2018
%p # second Maple program:
%p b:= proc(n, i, t) option remember; `if`(n=0, `if`(t=0, 1, 0),
%p `if`(i<1 or t<1, 0, b(n, i-1, t)+`if`(i^2>n, 0, b(n-i^2, i, t-1))))
%p end:
%p a:= n-> b(n, isqrt(n), 4):
%p seq(a(n), n=0..100); # _Alois P. Heinz_, Apr 14 2019
%t nn = 100; lim = Sqrt[nn]; t = Table[0, {nn}]; Do[n = a^2 + b^2 + c^2 + d^2; If[n <= nn, t[[n]]++], {a, lim}, {b, a, lim}, {c, b, lim}, {d, c, lim}]; t (* _T. D. Noe_, Sep 28 2012 *)
%t f[n_] := Length@ IntegerPartitions[n, {4}, Range[ Floor[ Sqrt[n - 1]]]^2]; Array[f, 105] (* _Robert G. Wilson v_, Sep 28 2012 *)
%o (PARI) A025428(n)=sum(a=1,n,sum(b=1,a,sum(c=1,b,sum(d=1,c,a^2+b^2+c^2+d^2==n))))
%o (PARI) A025428(n)=sum(a=1,sqrtint(max(n-3,0)), sum(b=1,min(sqrtint(n-a^2-2),a), sum(c=1,min(sqrtint(n-a^2-b^2-1),b),issquare(n-a^2-b^2-c^2,&d) & d <= c )))
%o (PARI) A025428(n)=sum(a=sqrtint(max(n,4)\4),sqrtint(max(n-3,0)), sum(b=sqrtint((n-a^2)\3-1)+1,min(sqrtint(n-a^2-2),a), sum(c=sqrtint((t=n-a^2-b^2)\2-1)+1, min(sqrtint(t-1),b), issquare(t-c^2) ))) \\ - _M. F. Hasler_, Sep 17 2012
%o for(n=1,100,print1(A025428(n),","))
%o (PARI) T(n)={a=matrix(n,4,i,j,0);for(d=1,sqrtint(n),forstep(i=n,d*d+1,-1,for(j=2,4,a[i,j]+=sum(k=1,j,if(k<j&&i-k*d*d>0,a[i-k*d*d,j-k],if(k==j&&i-k*d*d==0,1)))));a[d*d,1]=1);for(i=1,n,print(i" "a[i,4]))} /* _Robert Gerbicz_, Sep 28 2012 */
%Y Cf. A000414, A000534, A025357-A025375, A216374, A025416 (greedy inverse).
%Y Column k=4 of A243148.
%K nonn,easy
%O 0,29
%A _David W. Wilson_
%E Values of a(0..10^4) double-checked by _M. F. Hasler_, Sep 17 2012