|
%I #34 Feb 16 2015 04:17:11
%S 1,1,3,1,3,7,15,31,63,7,15,31,63,127,255,17,35,71,143,287,575,1151,
%T 2303,4607,9215,1843,3687,1229,2459,4919,9839,19679,39359,78719,
%U 157439,314879,629759,1259519,2519039,5038079,10076159,20152319,40304639
%N a(0)=a(1)=1; for n>1, a(n) = 2*a(n-1) + 1 if a(n-1) and n are coprime, otherwise a(n) = a(n-1)/gcd(a(n-1),n).
%C The initial value a(0)=1 is somehow artificial; using a(0)=0 would yield the same subsequent terms using the recurrence formula already for n=1. - _M. F. Hasler_, Feb 15 2015
%H Harvey P. Dale, <a href="/A133580/b133580.txt">Table of n, a(n) for n = 0..1000</a>
%e Write the GCD of a(n-1) and n under a(n-1):
%e n = : 0 1 2 3 4 5 6 7 8 9 ...
%e a(n): 1 1 3 1 3 7 15 31 63 7 ...
%e gcd : 1 1 3 1 1 1 1 1 9 1 ...
%t a = {1, 1}; Do[If[GCD[n, a[[ -1]]] == 1, b = 2*a[[ -1]] + 1, b = a[[ -1]]/GCD[a[[ -1]], n]]; AppendTo[a, b], {n, 2, 50}]; a (* _Stefan Steinerberger_, Dec 31 2007 *)
%t nxt[{a_,b_}]:={a+1,If[CoprimeQ[b,a+1],2b+1,b/GCD[b,a+1]]}; Join[{1}, Transpose[ NestList[nxt,{1,1},50]][[2]]] (* _Harvey P. Dale_, Sep 16 2012 *)
%o (PARI) A=vector(1000,i,1);for(n=2,#A,A[n]=if(gcd(A[n-1],n)>1,A[n-1]/gcd(A[n-1],n),A[n-1]*2+1)) \\ _M. F. Hasler_, Feb 15 2015
%o (PARI) a=0;#A133580=vector(1000,n,a=if(gcd(a,n)>1,a/gcd(a,n),a*2+1)) \\ _M. F. Hasler_, Feb 15 2015
%Y Cf. A133058, A133579.
%K nonn
%O 0,3
%A _Ctibor O. Zizka_, Dec 26 2007
%E Corrected and extended by _Stefan Steinerberger_, Dec 31 2007
%E Offset changed to 0 by _N. J. A. Sloane_, Feb 13 2015
%E The Mathematica programs are correct; b-file corrected by _Harvey P. Dale_, Feb 14 2015
|
