OFFSET
1,2
COMMENTS
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]
Is a(n) = A080468(n+1)+1?
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]
LINKS
Alois P. Heinz, Table of n, a(n) for n = 1..10000
Joerg Arndt, Matters Computational (The Fxtbook), section 1.26.3, p. 72.
David Wasserman, Quet transform + PARI code [Cached copy]
FORMULA
a(2^k-1) = 1. For 2^k <= n <= 2^(k+1)-2, a(n) = a(n-2^k+1)+1.
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]
EXAMPLE
a(4) = 2 because 4+2 (110) has two 1's, but 4+1 (101) has more than one 1.
Conjecture (Joerg Arndt):
a(n) is the number of bits in the binary words of sequence A108918
........00001.1.....00001.=..1.=..1
........00011.2.....00010.=..2.=..1.+.1
........00010.1.....00011.=..3.=..3
........00101.2.....00100.=..4.=..3.+.1
........00111.3.....00101.=..5.=..3.+.1.+.1
........00110.2.....00110.=..6.=..3.+.3
........00100.1.....00111.=..7.=..7
........01001.2.....01000.=..8.=..7.+.1
........01011.3.....01001.=..9.=..7.+.1.+.1
........01010.2.....01010.=.10.=..7.+.3
........01101.3.....01011.=.11.=..7.+.3.+.1
........01111.4.....01100.=.12.=..7.+.3.+.1.+.1
........01110.3.....01101.=.13.=..7.+.3.+.3
........01100.2.....01110.=.14.=..7.+.7
........01000.1.....01111.=.15.=.15
MAPLE
hb:= n-> `if`(n=1, 0, 1+hb(iquo (n, 2))):
a:= proc(n) local m, t;
m:= n;
for t from 0 while m>0 do
m:= m - (2^(hb(m+1))-1)
od; t
end:
seq(a(n), n=1..100); # Alois P. Heinz, Jan 22 2011
MATHEMATICA
hb[n_] := If[n==1, 0, 1+hb[Quotient[n, 2]]];
a[n_] := Module[{m=n, t}, For[t=0, m>0, t++, m = m - (2^(hb[m+1])-1)]; t];
Array[a, 100] (* Jean-François Alcover, Oct 31 2020, after Alois P. Heinz *)
PROG
(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
(PARI)
A100661(n)=
{ /* method: repeatedly subtract Mersenne numbers */
local(m, ct);
if ( n<=1, return(n) );
m = 1;
while ( n>m, m<<=1 );
m -= 1;
while ( m>n, m>>=1 );
/* here m=2^k-1 and m<=n */
ct = 0;
while ( n, while (m<=n, n-=m; ct+=1); m>>=1 );
return( ct );
}
vector(100, n, A100661(n)) /* show terms */
/* Joerg Arndt, Jan 22 2011 */
(PARI)
TInverse(v)=
{
local(l, w, used, start, x);
l = length(v); w = vector(l); used = vector(l); start = 1;
for (i = 1, l,
while (start <= l && used[start], start++);
x = start;
for (j = 2, v[i], x++; while (x <= l && used[x], x++));
if (x > l,
return (vector(i - 1, k, w[k]))
, /* else */
w[i] = x; used[x] = 1
)
);
return(w);
}
PInverse(v)=
{
local(l, w);
l = length(v); w = vector(l);
for (i = 1, l, if (v[i] <= l, w[v[i]] = i));
return(w);
}
T(v)=
{
local(l, w, c);
l = length(v); w = vector(l);
for (n = 1, l,
if (v[n],
c = 0;
for (m = 1, n - 1, if (v[m] < v[n], c++));
w[n] = v[n] - c
, /* else */
return (vector(n - 1, i, w[i]))
)
);
return(w);
}
Q(v)=T(PInverse(TInverse(v)));
/* compute terms: */
v = vector(150);
for (n = 1, 150, m = n; x = 1; while (!(m%2), m\=2; x *= 2); v[n] = x); Q(v)
CROSSREFS
KEYWORD
easy,nonn
AUTHOR
David Wasserman, Jan 14 2005
STATUS
approved