A001037 long by Dingledooper (xnor,explanation)

# 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).

To protect the system from spam, please input your favorite sport (hint: I believe its name must start with 'g', case insensitive)

download

return to the top page