login
Quet transform of A006519 (see A101387 for definition). Also, least k such that n+k has at most k ones in its binary representation.
6

%I #51 Oct 31 2020 12:51:46

%S 1,2,1,2,3,2,1,2,3,2,3,4,3,2,1,2,3,2,3,4,3,2,3,4,3,4,5,4,3,2,1,2,3,2,

%T 3,4,3,2,3,4,3,4,5,4,3,2,3,4,3,4,5,4,3,4,5,4,5,6,5,4,3,2,1,2,3,2,3,4,

%U 3,2,3,4,3,4,5,4,3,2,3,4,3,4,5,4,3,4,5,4,5,6,5,4,3,2,3,4,3,4,5,4,3,4,5,4,5

%N Quet transform of A006519 (see A101387 for definition). Also, least k such that n+k has at most k ones in its binary representation.

%C If n+a(n) has exactly a(n) 1's in binary, then a(n+1) = a(n)+1, but if n+a(n) has less than a(n) 1's, then a(n+1) = a(n)-1. a(n) is the number of terms needed to represent n as a sum of numbers of the form 2^k-1. [_Jeffrey Shallit_]

%C Is a(n) = A080468(n+1)+1?

%C Compute a(n) by repeatedly subtracting the largest number 2^k-1<=n until zero is reached. The number of times a term was subtracted gives a(n). Examples: 5 = 3 + 1 + 1 ==> a(5) = 3 6 = 3 + 3 ==> a(6) = 2. Replace all zeros in A079559 by -1, then the a(n) are obtained as cumulative sums (equivalent to the generating function given); see fxtbook link. [_Joerg Arndt_, Jun 12 2006]

%H Alois P. Heinz, <a href="/A100661/b100661.txt">Table of n, a(n) for n = 1..10000</a>

%H Joerg Arndt, <a href="http://www.jjj.de/fxt/#fxtbook">Matters Computational (The Fxtbook)</a>, section 1.26.3, p. 72.

%H David Wasserman, <a href="/A100661/a100661.txt">Quet transform + PARI code</a> [Cached copy]

%F a(2^k-1) = 1. For 2^k <= n <= 2^(k+1)-2, a(n) = a(n-2^k+1)+1.

%F G.f.: x*(2*(1-x)*prod(n>=1, (1+x^(2^n-1))) - 1)/((1-x)^2) = x*(1 + 2*x + 1*x^2 + 2*x^3 + 3*x^4 + 2*x^5 + 1*x^6 + 2*x^7 + 3*x^8 + 2*x^9 + ...) [_Joerg Arndt_, Jun 12 2006]

%e a(4) = 2 because 4+2 (110) has two 1's, but 4+1 (101) has more than one 1.

%e Conjecture (Joerg Arndt):

%e a(n) is the number of bits in the binary words of sequence A108918

%e ......A108918.A108918..n..=..n.=.(sum.of.term.2^k-1)

%e ........00001.1.....00001.=..1.=..1

%e ........00011.2.....00010.=..2.=..1.+.1

%e ........00010.1.....00011.=..3.=..3

%e ........00101.2.....00100.=..4.=..3.+.1

%e ........00111.3.....00101.=..5.=..3.+.1.+.1

%e ........00110.2.....00110.=..6.=..3.+.3

%e ........00100.1.....00111.=..7.=..7

%e ........01001.2.....01000.=..8.=..7.+.1

%e ........01011.3.....01001.=..9.=..7.+.1.+.1

%e ........01010.2.....01010.=.10.=..7.+.3

%e ........01101.3.....01011.=.11.=..7.+.3.+.1

%e ........01111.4.....01100.=.12.=..7.+.3.+.1.+.1

%e ........01110.3.....01101.=.13.=..7.+.3.+.3

%e ........01100.2.....01110.=.14.=..7.+.7

%e ........01000.1.....01111.=.15.=.15

%p hb:= n-> `if`(n=1, 0, 1+hb(iquo (n, 2))):

%p a:= proc(n) local m, t;

%p m:= n;

%p for t from 0 while m>0 do

%p m:= m - (2^(hb(m+1))-1)

%p od; t

%p end:

%p seq(a(n), n=1..100); # _Alois P. Heinz_, Jan 22 2011

%t hb[n_] := If[n==1, 0, 1+hb[Quotient[n, 2]]];

%t a[n_] := Module[{m=n, t}, For[t=0, m>0, t++, m = m - (2^(hb[m+1])-1)]; t];

%t Array[a, 100] (* _Jean-François Alcover_, Oct 31 2020, after _Alois P. Heinz_ *)

%o (Sage) A100661 = lambda n: next(k for k in PositiveIntegers() if (n+k).digits(base=2).count(1) <= k) # _D. S. McNeil_, Jan 23 2011

%o (PARI)

%o A100661(n)=

%o { /* method: repeatedly subtract Mersenne numbers */

%o local(m, ct);

%o if ( n<=1, return(n) );

%o m = 1;

%o while ( n>m, m<<=1 );

%o m -= 1;

%o while ( m>n, m>>=1 );

%o /* here m=2^k-1 and m<=n */

%o ct = 0;

%o while ( n, while (m<=n, n-=m; ct+=1); m>>=1 );

%o return( ct );

%o }

%o vector(100,n,A100661(n)) /* show terms */

%o /* _Joerg Arndt_, Jan 22 2011 */

%o (PARI)

%o TInverse(v)=

%o {

%o local(l, w, used, start, x);

%o l = length(v); w = vector(l); used = vector(l); start = 1;

%o for (i = 1, l,

%o while (start <= l && used[start], start++);

%o x = start;

%o for (j = 2, v[i], x++; while (x <= l && used[x], x++));

%o if (x > l,

%o return (vector(i - 1, k, w[k]))

%o , /* else */

%o w[i] = x; used[x] = 1

%o )

%o );

%o return(w);

%o }

%o PInverse(v)=

%o {

%o local(l, w);

%o l = length(v); w = vector(l);

%o for (i = 1, l, if (v[i] <= l, w[v[i]] = i));

%o return(w);

%o }

%o T(v)=

%o {

%o local(l, w, c);

%o l = length(v); w = vector(l);

%o for (n = 1, l,

%o if (v[n],

%o c = 0;

%o for (m = 1, n - 1, if (v[m] < v[n], c++));

%o w[n] = v[n] - c

%o , /* else */

%o return (vector(n - 1, i, w[i]))

%o )

%o );

%o return(w);

%o }

%o Q(v)=T(PInverse(TInverse(v)));

%o /* compute terms: */

%o v = vector(150);

%o for (n = 1, 150, m = n; x = 1; while (!(m%2), m\=2; x *= 2); v[n] = x); Q(v)

%Y Cf. A006519, A080468, A101387, A100808.

%Y Records at indices in (essentially) A000325.

%K easy,nonn

%O 1,2

%A _David Wasserman_, Jan 14 2005