######################################################### # # # Coded by Indranil Ghosh (indranilg49@gmail.com) # # # ######################################################### #Python 2.7.11, OEIS sequence: A286598 from sympy import prime, factorint, log, floor def A(n): return n - 2**int(floor(log(n, 2))) def b(n): return n + 1 if n<2 else prime(1 + (len(bin(n)[2:]) - bin(n)[2:].count("1"))) * b(A(n)) def P(n): f = factorint(n) return sorted([f[i] for i in f]) def a046523(n): x=1 while True: if P(n) == P(x): return x else: x+=1 def a278222(n): return a046523(b(n)) def ok(n): return n&(n - 1)==0 def a153152(n): return n if n<2 else (n + 1)/2 if ok(n + 1) else n + 1 def x(n): return (int(bin(n)[2:][::-1], 2) - 1)/2 def msb(n): return n if n<3 else msb(n/2)*2 def a059893(n): return x(n) + msb(n) def a153142(n): return 0 if n==0 else a059893(a153152(a059893(n))) def a(n): return n + 1 if n<2 else a278222(a153142(n)) print [a153152(n) for n in xrange(0, 101)]