This site is supported by donations to The OEIS Foundation.

User:Ruud H.G. van Tol/Collatz

From OeisWiki
Jump to: navigation, search

A060322 add PROG

first(N) = {
  my( a=Vec([1, 1, 1], N), x=[4] );
  for( n=4, N
  , my( w=List() );
    for( i=1, #x
    , my( q=2*x[i] );
      if( 1 == q%3, listput(w, (q-1)/3) );
      listput(w, q);
    );
    x= Vec(w);
    a[n]= #x;
  );
  a;
}

? first(20)
%84 = [1, 1, 1, 1, 2, 3, 4, 5, 6, 8, 12, 18, 24, 31, 39, 50, 68, 91, 120, 159]



Faster:

? first(N) = {
  my( a= Vec([1,0,0,0,1], N), Y= [5], Z=[] );
  for( n= 6, N
  , my( X= List() );
    foreach( Z, x
    , listput(X, 4*x+1);
      if( 1 == x%6, listput(X, (4*x-1)\3) )
    );
    Z= Y;
    foreach( Z, x
    , if( 5 == x%6, listput(X, (2*x-1)\3) );
    );
    Y= Vec(X);
    a[n]= #Y;
  );
  for( i= 2, #a, a[i]+= a[i-1] );
  a;
}

? first(20)
[1, 1, 1, 1, 2, 3, 4, 5, 6, 8, 12, 18, 24, 31, 39, 50, 68, 91, 120, 159]
 

On A254070 (first 1 (mod 4) successor of 4n-3)

a(n) = my(x=3*n-2); x=(x>>valuation(x, 2))+1; (3/2)^(valuation(x, 2)-1)*x-1;
 

On A350082 (Collatz trajectory)

Sequence of numbers that are not in this sequence:

[13, 21, 31, 35, 45, 49, 53, 61, 65, 67, 69, 77, 83, 85, 91, 93, 99, 101, 103, ...]
 

To do: Come up with a variant without the zeroes.

(m,v)

((2*m+1)*2^v - 1) => ((2*m+1)*3^v - 1) / 2

 
I'm currently exploring how presenting numbers x as

m >= 0: x = m * 2^(v+1) + 2^v - 1.

Even numbers: v = 0.
Odd numbers: v >= 1.

exposes more of how they behave under (3x+1)/2.

With x in base-2, this is a split at the rightmost 0-bit,
so v is the number of trailing 1-bits ("valuation"),
and m is the "multiplier" in front.

Example: 27 = 0b11_0_11 = 3 * 2^3 + (2^2 - 1).

I'll also use an (m,v)-tuple-notation.
Example: 27 = (3,2).

See also A025480 (Grundy or nim-values).

For v >= 1:
m * 2^(v+1) + 2^v - 1 ->
(3 * (m * 2^(v+1) + 2^v - 1) + 1) / 2 =
3 * m * 2^v + 3 * 2^(v-1) - 1 (is odd) -> ... ->
(3^v * (m * 2 + 1) - 1) / 2.

Alternative view: (2*m + 1) * 2^v - 1.

Example for (m,v) = (1,2) :
11 =  1 * 2^3 + 2^2 - 1 = ( 1,2) ->
17 =  4 * 2^2 + 2^1 - 1 = ( 4,1) ->
26 = 13 * 2^1 + 2^0 - 1 = (13,0) -> 13 = (3,1).
Or in a single step:
11 = (1,2) -> (3^2*(1*2+1)-1)/2 = 13 = (3,1).
Next: (3,1) -> (3*7-1)/2 = 10 -> 5 = (1,1).
Next: 4 -> 1 = (0,1).


Example for (m,v) = (3,2) :
27 =  3 * 2^3 + 2^2 - 1 = ( 3,2) ->
41 = 10 * 2^2 + 2^1 - 1 = (10,1) ->
62 = 31 * 2^1 + 2^0 - 1 = (31,0) -> 31 = (0,5).
Or in a single step:
27 = (3,2) -> (3^2*(3*2+1)-1)/2 = 31 = (0,5).
Next: (0,5) -> (0+3^5*1-1)/2 = 121 = (30,1).
Next: (3*61-1)/2 = 91 = (11,2).
 

(m,v)_2 -> (m,v)_3

m * 2^(v+1) + 2^v - 1 -> m * 3^v + (3^v - 1)/2

Examples:
 27 = (  3,2)_2 -> (  3,2)_3 =   3*  9+(  9-1)/2 =  31
 31 = (  0,5)_2 -> (  0,5)_3 =   0*243+(243-1)/2 = 121
121 = ( 30,1)_2 -> ( 30,1)_3 =  30*  3+(  3-1)/2 =  91
 91 = ( 11,2)_2 -> ( 11,2)_3 =  11*  9+(  9-1)/2 = 103
103 = (  6,3)_2 -> (  6,3)_3 =   6* 27+( 27-1)/2 = 175
175 = (  5,4)_2 -> (  5,4)_3 =   5* 81+( 81-1)/2 = 445
445 = (111,1)_2 -> (111,1)_3 = 111*  3+(  3-1)/2 = 334 -> 167
167 = ( 10,3)_2 -> ( 10,3)_3 =  10* 27+( 27-1)/2 = 283
(etc.)
 

See also: A259663, and
Francis Laclé; 2-adic parity explorations of the 3n+1 problem; 2021;
https://hal.science/hal-03201180v2/document

(m,v) notation of powers of 3

n-even: 3^n -> ((3^n-1)/4, 1)_2 -> (3^n-1)/4 * 3^1 + (3^1-1)/2
                                 = 3*(3^n-1)/4 + 1
Examples:
  1 ->   1
  9 ->   7
 81 ->  61
729 -> 547


n-odd:  3^n -> ((3^n-3)/8, 2)_2 -> (3^n-3)/8 * 3^2 + (3^2-1)/2
                                 = 9*(3^n-3)/8 + 4
Examples:
    3 ->     4,   1
   27 ->    31
  243 ->   274, 137
 2187 ->  2461
19683 -> 22144, 173
 

On A100982

a_iter(N) = {
  my(a=List([[0,0,2,1,1,0]]), p0=#a, p1=#a);
  for(p3=0,N-1
  , my(p2=logint(3^p3,2)+1, m=if(logint(3^(p3+1),2)>p2,2,1));
    for(p=p0,p1
    , my(s=a[p][1]);
      for(i=0,if(s,valuation(s,2)-1,0)
      , my(t=[(a[p][1]+(1<<i))<<m,0,0,p2+m,0,p3+1]);
        listput(a,t)
      )
    );
    p0=p1+1;
    p1=#a
  );
  for(i=2,#a
  , my(t=a[i][1], p2=a[i][4], p3=a[i][6], r=0);
    for(j=1,p2, r=if(bittest(t, p2-j), 3*r+1, r) / 2);
    a[i][2]=numerator(r);
    a[i][3]=-r * (2^p2) / (3^p3) % (2^p2);
    a[i][5]=a[i][3] * 3^p3 / 2^p2 + r
  );
  Vec(a);
}

? printp( Mat(a_iter(5)~) )
[  0   0   2 1   1 0]
[  2   1   1 2   1 1]
[ 12   5   3 4   2 2]
[ 26  23  11 5  10 3]
[ 28  19  23 5  20 3]
[108  85  59 7  38 4]
[116  73   7 7   5 4]
[120  65  15 7  10 4]
[218 319 123 8 118 5]
[220 287 219 8 209 5]
[234 283 199 8 190 5]
[236 251  39 8  38 5]
[242 259  79 8  76 5]
[244 227 175 8 167 5]
[248 211  95 8  91 5]

Examples:
  2 + i*2^ 1 ->   1 + i*3^ 0
  1 + i*2^ 2 ->   1 + i*3^ 1
  3 + i*2^ 4 ->   2 + i*3^ 2
 11 + i*2^ 5 ->  10 + i*3^ 3
 23 + i*2^ 5 ->  20 + i*3^ 3
 59 + i*2^ 7 ->  38 + i*3^ 4
  7 + i*2^ 7 ->   5 + i*3^ 4
 15 + i*2^ 7 ->  10 + i*3^ 4
123 + i*2^ 8 -> 118 + i*3^ 5
219 + i*2^ 8 -> 209 + i*3^ 5
199 + i*2^ 8 -> 190 + i*3^ 5
 39 + i*2^ 8 ->  38 + i*3^ 5
 79 + i*2^ 8 ->  76 + i*3^ 5
175 + i*2^ 8 -> 167 + i*3^ 5
 95 + i*2^ 8 ->  91 + i*3^ 5
 

Permutation of positive integers

? my(N=2^5, v=Vec(0,N), seen=Map(), i=0, m=1, n=1); while(i<N, v[i++]=n; mapput(~seen,n,0); n=if(n%2,n*3+1,n)/2; if(mapisdefined(seen,n),while(mapisdefined(seen,m+=2),);n=m));v
% [1, 2, 3, 5, 8, 4, 7, 11, 17, 26, 13, 20, 10, 9, 14, 15, 23, 35, 53, 80, 40, 19, 29, 44, 22, 21, 32, 16, 25, 38, 27, 41]