Göbel's Sequence

Author

Eric W. Weisstein
April 10, 2024

This notebook downloaded from http://mathworld.wolfram.com/notebooks/IntegerSequences/GoebelsSequence.nb.

For more information, see Eric's MathWorld entry http://mathworld.wolfram.com/GoebelsSequence.html.

©2024 Wolfram Research, Inc. except for portions noted otherwise

Exponent

A115632

N[CC=(10^(-104))FromDigits[{1,0,4,7,8,3,1,4,4,7,5,7,6,4,1,1,2,2,9,5,5,9,9,0,9,4,6,2,7,4,3,1,3,7,5,5,4,5,9,0,5,8,7,6,1,2,8,6,0,2,3,3,0,9,6,9,5,1,0,4,0,6,4,8,5,3,5,3,6,0,5,9,0,4,9,7,2,6,2,3,1,7,9,7,5,1,3,0,9,7,9,0,0,0,7,0,9,9,4,7,9,5,1,1,3}],20]

Exact

f[0]:=1
f[n_]:=f[n]=(1+Sum[f[k]^2,{k,0,n-1}])/n

f/@Range[0,10]

N[With[{n=25},(((((f[n])/(n+2-(n^(-1))+4(n^(-2))))))^(2^(-n)))],100]//Timing

N[With[{n=30},(((((f[n])/(n+2-(n^(-1))+4(n^(-2))-21(n^(-3))))))^(2^(-n)))],100]//Timing

seq=N[Table[(((((f[n])/(n+2-(n^(-1))+4(n^(-2))-21(n^(-3))))))^(2^(-n))),{n,31}],100]

ListPlot[seq,PlotRange->CC+{-1,1}(10^(-6))]

ResourceFunction["SequenceLimit"][seq]

CC1=ResourceFunction["SequenceLimit"][Drop[seq,15]]

seq[[-1]]-CC//N

CC1-CC//N

Numeric

f[0]:=1`100
f[n_]:=f[n]=(1+Sum[f[k]^2,{k,0,n-1}])/n

N[With[{n=25},(((((f[n])/(n+2-(n^(-1))+4(n^(-2))))))^(2^(-n)))],100]//Timing

N[With[{n=30},(((((f[n])/(n+2-(n^(-1))+4(n^(-2))-21(n^(-3))))))^(2^(-n)))],100]//Timing

seq=Table[(((((f[n])/(n+2-(n^(-1))+4(n^(-2))-21(n^(-3))))))^(2^(-n))),{n,55}];

ListPlot[seq,PlotRange->CC+{-1,1}(10^(-6))]

ResourceFunction["SequenceLimit"][seq]

CC2=ResourceFunction["SequenceLimit"][Drop[seq,15]]

seq[[-1]]-CC//N

CC2-CC//N

Somos

From Michael Somos in A115632:

{a(n)=local(t=log(2)/2); for(k=2, 14, t+= (log(1+(k-1)/exp(2^(k-1)*t))-log(k))/2^k); t=exp(t-suminf(k=15, log(k)/2^k)); floor(t*10^(n-1))%10}

t=((Log[2])/2);

Do[t+=((Log[1+((k-1)/(Exp[2^(k-1)t]))]-Log[k])/(2^k)),{k,2,14}]

tinf=Sum[((Log[k])/(2^k)),{k,15,∞}]//FullSimplify

CC3=Exp[t-tinf]

N[CC3,120]

N[CC3-CC,150]//N

Transformed sequence

A061322

s[1]=2;
s[n_]:=s[n-1]+(((s[n-1])^2)/(((n-1))^2))

Table[s[n],{n,10}]

FactorInteger/@%//Column

Do[Print[{n,IntegerQ[s[n]]//Timing,IntegerLength[s[n]]}],{n,40}]

Table[Mod[s[n],k],{k,Prime[Range[7]]},{n,k}]//Column

a[1,m_]=2;
a[n_,m_]:=a[n,m]=PolynomialMod[a[n-1,m]+(((a[n-1,m])^2)/(((n-1))^2)),m]
a[n_]:=a[n,n]

Table[a[n,k],{k,Prime[Range[14]]},{n,k}]//Column

Prime[14]

nonintegral=Table[a[p],{p,Prime[Range[100]]}]

Prime[Position[nonintegral,Except[0,_Integer]]//Flatten]