This is a front-end for the Online Encyclopedia of Integer Sequences, made by Christian Perfect. The idea is to provide OEIS entries in non-ancient HTML, and then to think about how they're presented visually. The source code is on GitHub.
1=1^2 so a(1)=1, 2=1^2+1^2 so a(2)=1 and a(3)=1, 3=1^2+1^2+1^2 so a(4)=1 and a(5)=1 and a(6)=1, 4=2^2 so a(7)=2, 5=1^2+2^2 so a(8)=1 and a(9)=2, 6=1^2+1^2+2^2 so a(10)=1 and a(11)=1 and a(12)=2, 7=1^2+1^2+1^2+2^2 so a(13)=1 and a(14)=1 and a(15)=1 and a(16)=2, 8=2^2+2^2 so a(17)=2 and a(18)=2, 9=3^2 so a(19)=3, etc...
a(1) = 7 because A002828(1*1) = 1, A002828(2*1) = 2, A002828(3*1) = 3, A002828(5*1) = 2, A002828(6*1) = 3, ..., and 7 is the smallest positive multiplier leading to A002828(7*1) = 7.
istwo(n:int)=my(f); if(n<3, return(n>=0); ); f=factor(n>>valuation(n, 2)); for(i=1, #f[, 1], if(bitand(f[i, 2], 1)==1&&bitand(f[i, 1], 3)==3, return(0))); 1; isthree(n:int)=my(tmp=valuation(n, 2)); bitand(tmp, 1)||bitand(n>>tmp, 7)!=7; a002828(n)=if(issquare(n), !!n, if(istwo(n), 2, 4-isthree(n))); \\ A002828 a(n) = {my(m=1); while(a002828(m*n)!=4, m++); m; } \\ Michel Marcus, Apr 17 2018
75 = 1^2 + 5^2 + 7^2 = 5^2 + 5^2 + 5^2, 76 = 2^2 + 6^2 + 6^2, 77 = 2^2 + 3^2 + 8^2 = 4^2 + 5^2 + 6^2, 78 = 2^2 + 5^2 + 7^2. These the first 4 consecutive numbers with the same smallest number of squares needed to represent, so a(4) = 75.
a001481 n = a001481_list !! (n-1) a001481_list = [x | x <- [0..], a000161 x > 0] -- Reinhard Zumkeller, Feb 14 2012, Aug 16 2011
[n: n in [0..160] | NormEquation(1, n) eq true]; // Arkadiusz Wesolowski, May 11 2016
readlib(issqr): for n from 0 to 160 do for k from 0 to floor(sqrt(n)) do if issqr(n-k^2) then printf(`%d,`,n); break fi: od: od:
upTo = 160; With[{max = Ceiling[Sqrt[upTo]]}, Select[Union[Total /@ (Tuples[Range[0, max], {2}]^2)], # <= upTo &]] (* Harvey P. Dale, Apr 22 2011 *)
Select[Range[0, 160], SquaresR[2, #] != 0 &] (* Jean-François Alcover, Jan 04 2013 *)
isA001481(n)=local(x,r);x=0;r=0;while(x<=sqrt(n) && r==0,if(issquare(n-x^2),r=1);x++);r \\ Michael B. Porter, Oct 31 2009
is(n)=my(f=factor(n));for(i=1,#f[,1],if(f[i,2]%2 && f[i,1]%4==3, return(0))); 1 \\ Charles R Greathouse IV, Aug 24 2012
B=bnfinit('z^2+1,1);
is(n)=#bnfisintnorm(B,n) \\ Ralf Stephan, Oct 18 2013, edited by M. F. Hasler, Nov 21 2017
list(lim)=my(v=List(),t); for(m=0,sqrtint(lim\=1), t=m^2; for(n=0, min(sqrtint(lim-t),m), listput(v,t+n^2))); Set(v) \\ Charles R Greathouse IV, Jan 05 2016
is_A001481(n)=!for(i=2-bittest(n,0),#n=factor(n)~, bittest(n[1,i],1)&&bittest(n[2,i],0)&&return) \\ M. F. Hasler, Nov 20 2017
from itertools import count, islice from sympy import factorint def A001481_gen(): # generator of terms return filter(lambda n:(lambda m:all(d & 3 != 3 or m[d] & 1 == 0 for d in m))(factorint(n)),count(0)) A001481_list = list(islice(A001481_gen(),30)) # Chai Wah Wu, Jun 27 2022
15 is in the sequence because it is the sum of four squares, namely, 3^2 + 2^2 + 1^2 + 1^2, and it can't be expressed as the sum of fewer squares. 16 is not in the sequence, because, although it can be expressed as 2^2 + 2^2 + 2^2 + 2^2, it can also be expressed as 4^2.
a004215 n = a004215_list !! (n-1) a004215_list = filter ((== 4) . a002828) [1..] -- Reinhard Zumkeller, Feb 26 2015
N:= 1000: # to get all terms <= N
{seq(seq(4^i * (8*j + 7), j = 0 .. floor((N/4^i - 7)/8)), i = 0 .. floor(log[4](N)))}; # Robert Israel, Sep 02 2014
Sort[Flatten[Table[4^i(8j + 7), {i, 0, 2}, {j, 0, 42}]]] (* Alonso del Arte, Jul 05 2005 *)
Select[Range[120], Mod[ #/4^IntegerExponent[ #, 4], 8] == 7 &] (* Ant King, Oct 14 2010 *)
isA004215(n)={ local(fouri,j) ; fouri=1 ; while( n >=7*fouri, if( n % fouri ==0, j= n/fouri -7 ; if( j % 8 ==0, return(1) ) ; ) ; fouri *= 4 ; ) ; return(0) ; } { for(n=1,400, if(isA004215(n), print1(n,",") ; ) ; ) ; } \\ R. J. Mathar, Nov 22 2006
isA004215(n)= n\4^valuation(n,4)%8==7 \\ M. F. Hasler, Mar 18 2011
def valuation(n, b):
v = 0
while n > 1 and n%b == 0: n //= b; v += 1
return v
def ok(n): return n//4**valuation(n, 4)%8 == 7 # after M. F. Hasler
print(list(filter(ok, range(344)))) # Michael S. Branicky, Jul 15 2021
from itertools import count, islice def A004215_gen(startvalue=1): # generator of terms >= startvalue return filter(lambda n:not (m:=(~n&n-1).bit_length())&1 and (n>>m)&7==7,count(max(startvalue,1))) A004215_list = list(islice(A004215_gen(),30)) # Chai Wah Wu, Jul 09 2022
def A004215(n): def bisection(f,kmin=0,kmax=1): while f(kmax) > kmax: kmax <<= 1 kmin = kmax >> 1 while kmax-kmin > 1: kmid = kmax+kmin>>1 if f(kmid) <= kmid: kmax = kmid else: kmin = kmid return kmax def f(x): return n+x-sum(((x>>(i<<1))-7>>3)+1 for i in range(x.bit_length()>>1)) return bisection(f,n,n) # Chai Wah Wu, Feb 14 2025
a(1) = 0 = 0^2 + 0^2 + 0^2. A005875(0) = 1 = A000164(0). a(9) = 9 = 0^2 + 0^2 + 3^2 = 1^2 + 2^2 + 2^2. A000164(9) = 2. A000164(9) = 30 = 2*3 + 8*3 (counting signs and order). - _Wolfdieter Lang_, Apr 08 2013
isA000378 := proc(n) # return true or false depending on n being in the list local x,y ; for x from 0 do if 3*x^2 > n then return false; end if; for y from x do if x^2+2*y^2 > n then break; else if issqr(n-x^2-y^2) then return true; end if; end if; end do: end do: end proc: A000378 := proc(n) # generate A000378(n) option remember; local a; if n = 1 then 0; else for a from procname(n-1)+1 do if isA000378(a) then return a; end if; end do: end if; end proc: seq(A000378(n),n=1..100) ; # R. J. Mathar, Sep 09 2015
okQ[n_] := If[EvenQ[k = IntegerExponent[n, 2]], m = n/2^k; Mod[m, 8] != 7, True]; Select[Range[0, 100], okQ] (* Jean-François Alcover, Feb 08 2016, adapted from PARI *)
isA000378(n)=my(k=valuation(n, 2)); if(k%2==0, n>>=k; n%8!=7, 1)
list(lim)=my(v=List(),k,t); for(x=0,sqrtint(lim\=1), for(y=0, min(sqrtint(lim-x^2),x), k=x^2+y^2; for(z=0,min(sqrtint(lim-k), y), listput(v,k+z^2)))); Set(v) \\ Charles R Greathouse IV, Sep 14 2015
def valuation(n, b):
v = 0
while n > 1 and n%b == 0: n //= b; v += 1
return v
def ok(n): return n//4**valuation(n, 4)%8 != 7
print(list(filter(ok, range(84)))) # Michael S. Branicky, Jul 15 2021
from itertools import count, islice def A000378_gen(): # generator of terms return filter(lambda n:n>>2*(bin(n)[:1:-1].index('1')//2) & 7 < 7, count(1)) A000378_list = list(islice(A000378_gen(),30)) # Chai Wah Wu, Jun 27 2022
def A000378(n): def f(x): return n-1+sum(((x>>(i<<1))-7>>3)+1 for i in range(x.bit_length()>>1)) m, k = n-1, f(n-1) while m != k: m, k = k, f(k) return m # Chai Wah Wu, Feb 14 2025
Cnt4[n_] := Module[{k = 1}, While[Length[PowersRepresentations[n, k, 4]] == 0, k++]; k]; Array[Cnt4, 100] (* T. D. Noe, Apr 01 2011 *)
seq[n_] := Module[{v = Table[0, {n}], s, p}, s = Sum[x^(k^4), {k, 1, n^(1/4)}] + O[x]^(n+1); p=1; For[k=1, k <= 19, k++, p *= s; For[i=1, i <= n, i++, If[v[[i]]==0 && Coefficient[p, x, i] != 0, v[[i]] = k]]]; v];
seq[100] (* Jean-François Alcover, Sep 28 2019, after Andrew Howroyd *)
seq(n)={my(v=vector(n), s=sum(k=1, sqrtint(sqrtint(n)), x^(k^4)) + O(x*x^n), p=1); for(k=1, 19, p*=s; for(i=1, n, if(!v[i] && polcoeff(p,i), v[i]=k))); v} \\ Andrew Howroyd, Jul 06 2018
from itertools import count from sympy.solvers.diophantine.diophantine import power_representation def A002377(n): if n == 1: return 1 for k in count(1): try: next(power_representation(n,4,k)) except: continue return k # Chai Wah Wu, Jun 25 2024
7=4+1+1+1, so 7 requires 4 squares using the greedy algorithm, so a(7)=4.
a053610 n = s n $ reverse $ takeWhile (<= n) $ tail a000290_list where
s _ [] = 0
s m (x:xs) | x > m = s m xs
| otherwise = m' + s r xs where (m',r) = divMod m x
-- Reinhard Zumkeller, May 08 2011
A053610 := proc(n) local a,x; a := 0 ; x := n ; while x > 0 do x := x-A048760(x) ; a := a+1 ; end do: a ; end proc: # R. J. Mathar, May 13 2016
f[n_] := (n - Floor[Sqrt[n]]^2); g[n_] := (m = n; c = 1; While[a = f[m]; a != 0, c++; m = a]; c); Table[ g[n], {n, 1, 105}]
A053610(n,c=1)=while(n-=sqrtint(n)^2,c++);c \\ M. F. Hasler, Dec 04 2008
from math import isqrt def A053610(n): c = 0 while n: n -= isqrt(n)**2 c += 1 return c # Chai Wah Wu, Aug 01 2023
Comments