A133580 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).
1, 1, 3, 1, 3, 7, 15, 31, 63, 7, 15, 31, 63, 127, 255, 17, 35, 71, 143, 287, 575, 1151, 2303, 4607, 9215, 1843, 3687, 1229, 2459, 4919, 9839, 19679, 39359, 78719, 157439, 314879, 629759, 1259519, 2519039, 5038079, 10076159, 20152319, 40304639
Offset: 0
Keywords
Examples
Write the GCD of a(n-1) and n under a(n-1): n = : 0 1 2 3 4 5 6 7 8 9 ... a(n): 1 1 3 1 3 7 15 31 63 7 ... gcd : 1 1 3 1 1 1 1 1 9 1 ...
Links
- Harvey P. Dale, Table of n, a(n) for n = 0..1000
Programs
-
Mathematica
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 *) 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 *)
-
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
-
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
Extensions
Corrected and extended by Stefan Steinerberger, Dec 31 2007
Offset changed to 0 by N. J. A. Sloane, Feb 13 2015
The Mathematica programs are correct; b-file corrected by Harvey P. Dale, Feb 14 2015
Comments