# According to the OEIS page, this sequence is Column 2 of A074650. # # One identity from the A074650 page states that: # T(n,k) = (k^n - Sum_{d<n,d|n} d*T(d,k)) / n # # We are interested in the special case where k=2, so this can be rewritten as # a(n) = (2^n - Sum_{d<n,d|n} d*a(d)) / n # # Here is the ungolfed implementation: # # a = [1] # for n in range(1, 51): # x = (2**n - sum(a[d]*d for d in range(1,n) if n%d<1)) / n # print x # a += x, # # Notice that `a[d]*d` stays constant for any `d`, so we can do it while appending `x`: # # a = [1] # for n in range(1, 51): # x = (2**n - sum(a[d] for d in range(1,n) if n%d<1)) / n # print x # a += x*n, # # This can be simplified further because the division and multiplication by `n` cancel out: # # a = [1] # for n in range(1, 51): # x = 2**n - sum(a[d] for d in range(1,n) if n%d<1) # print x/n # a += x, # # Notice that a[n] is used in the sum of every multiple of `n`. Instead # of summing everything at `n`, we can sum "on the fly" like a sieve: # # a = [1]+[0]*50 # for n in range(1, 51): # x = 2**n - a[n] # print x/n # for i in range(0, 51, n): # a[i] += x # # Now for the trick. We change `a` into an integer and perform all array operations # through bit manipulations. The first 64 bits of `a` represents `a[0]`, the next 64 bits # represents `a[1]`, and so on. 64 is an arbitrary number, but we need enough bits in # order to store the maximum possible value that `a[n]` can reach. A naive implementation # goes like this: # # a = C = 1 # for n in range(1, 51): # C <<= 64 # x = 2**n - a/C%2**64 # print x/n # for i in range(0, n*99, n): # a += x<<64*i # # Let's focus on the inner loop. We want to add `x` to every n-th cell up to 51. How to # do it without a loop? Well we can generate the required bitstring beforehand, and simply # perform "vector addition" on `a`. Assume n=3. Then, in base 2^64, we want to generate # this literal: # 'x00x00x...x00x00x00x00x00x00x00x' # # In base 2^(64*n) representation, this is simply: # 'xxx...xxxxxxxx' # # But how do we generate this literal? Observe that the expression 10**n/9 gives us # a string of z ones. So in base C, this expression would be C**z/(C-1). Since we # want a string of z `x`s, not ones, we just need to multiply C**z/(C-1) by x. # Now what about z? How large does it need to be? Technically n needs to be at least 51, # but because of floor division we can sneakily get away with n=9, like in @xnor's solution. # # Here's the updated code: # # a = C = 1 # for n in range(1, 51): # C <<= 64 # x = 2**n - a/C%2**64 # print x/n # a += C**9/(C-1)*x # # One more nice trick. a[n] will always be less than 2**n, so we can replace the # `%2**64` with `%2**n`. This allows us to get rid of the other 2**n because # y-x%y == -x%y, except when x=0. Unfortunately, x=a/C does equal 0 on the first # iteration, but we hack it out by replacing `print x/n` with `print-~x/n`. # # a = C = 1 # for n in range(1, 51): # C <<= 64 # x = -a/C%2**n # print-~x/n # a += C**9/(C-1)*x # # A few more golfs, and we reach @xnor's 66-byte solution, hooray!: a=C=n=1 exec'C<<=64;x=-a/C%2**n;print-~x/n;a+=C**9/~-C*x;n+=1;'*50
Note that non-ascii characters in the above source code will be escaped (such as \x9f).