Non-sums of distinct positive powers

Question

There are 31 positive integers that cannot be expressed as the sum of 1 or more distinct positive squares:

2, 3, 6, 7, 8, 11, 12, 15, 18, 19, 22, 23, 24, 27, 28, 31, 32, 33, 43, 44, 47, 48, 60, 67, 72, 76, 92, 96, 108, 112, 128

There are 2788 positive integers that cannot be expressed as the sum of 1 or more distinct positive cubes:

2, 3, 4, 5, 6, 7, 10, 11, 12, 13, 14, 15, 16, 17, 18, 19, 20, 21, 22, 23, 24, 25, 26, 29, 30, 31, 32, 33, 34, 37, 38, 39, 40, 41, 42, 43, 44, 45, 46, 47, 48, 49, 50, 51, 52, 53, 54, 55, 56, 57, 58, 59, 60, 61, 62, 63, 66, 67, 68, 69, 70, 71, 74, 75, 76, ...

It has been proven that this sequence is finite for all powers, and the largest term in each sequence is part of the Anti-Waring numbers, specifically \$N(k,1)-1\$, for a power \$k\$.

You are to take a positive integer \$n > 1\$ as input, and output all positive integers that cannot be expressed as the sum of 1 or more distinct positive integers to the power of \$n\$. You may output the integers in any clear and consistent format (as digits, unary, Church numerals etc.) and you may use any clear and consistent separator between the numbers (newlines, spaces, as a list, etc.).

Your program must terminate, through error or natural termination, in finite time; iterating over all integers and outputting those that cannot be expressed in such a form is not an acceptable method.

This is code-golf, so the shortest code in bytes wins.

Test cases

2 -> 2, 3, 6, 7, 8, 11, 12, 15, 18, 19, 22, 23, 24, 27, 28, 31, 32, 33, 43, 44, 47, 48, 60, 67, 72, 76, 92, 96, 108, 112, 128
3 -> 2, 3, 4, 5, 6, 7, 10, 11, 12, 13, 14, 15, 16, 17, 18, 19, 20, 21, 22, 23, 24, 25, 26, 29, 30, 31, 32, 33, 34, 37, 38, 39, 40, 41, 42, 43, 44, 45, 46, 47, ...
4 -> 2, 3, 4, 5, 6, 7, 8, 9, 10, 11, 12, 13, 14, 15, 18, 19, 20, 21, 22, 23, 24, 25, 26, 27, 28, 29, 30, 31, 32, 33, 34, 35, 36, 37, 38, 39, 40, 41, 42, 43, ...
5 -> 2, 3, 4, 5, 6, 7, 8, 9, 10, 11, 12, 13, 14, 15, 16, 17, 18, 19, 20, 21, 22, 23, 24, 25, 26, 27, 28, 29, 30, 31, 34, 35, 36, 37, 38, 39, 40, 41, 42, 43, ...
6 -> 2, 3, 4, 5, 6, 7, 8, 9, 10, 11, 12, 13, 14, 15, 16, 17, 18, 19, 20, 21, 22, 23, 24, 25, 26, 27, 28, 29, 30, 31, 32, 33, 34, 35, 36, 37, 38, 39, 40, 41, 42, 43, 44, 45, 46, 47, 48, 49, 50, 51, 52, 53, 54, 55, 56, 57, 58, 59, 60, 61, 62, 63, 66, 67, 68, 69, 70, 71, 72, 73, ...
7 -> 2, 3, 4, 5, 6, 7, 8, 9, 10, 11, 12, 13, 14, 15, 16, 17, 18, 19, 20, 21, 22, 23, 24, 25, 26, 27, 28, 29, 30, 31, 32, 33, 34, 35, 36, 37, 38, 39, 40, 41, 42, 43, 44, 45, 46, 47, 48, 49, 50, 51, 52, 53, 54, 55, 56, 57, 58, 59, 60, 61, 62, 63, 64, 65, 66, 67, 68, 69, 70, 71, 72, 73, 74, 75, 76, 77, 78, 79, 80, 81, 82, 83, 84, 85, 86, 87, 88, 89, 90, 91, 92, 93, 94, 95, 96, 97, 98, 99, 100, 101, 102, 103, 104, 105, 106, 107, 108, 109, 110, 111, 112, 113, 114, 115, 116, 117, 118, 119, 120, 121, 122, 123, 124, 125, 126, 127, 130, 131, 132, 133, 134, 135, 136, 137, 138, 139, 140, 141, 142, 143, 144, 145, 146, 147, 148, 149, 150, 151, 152, 153, ...

Answer 1

APL(Dyalog Unicode), 26 bytes ^SBCS

⊢{z~∊(⊢,+)\⍺*⍨z←⍳!⍵}!×2*×⍨

Try it on APLgolf!

This only works theoretically in multiple ways. The upper bound used (factorial of c = n! * 2^(n^2)is way too high to be accurately computed in floats, it would take forever to accurately compute it as a bigint, and it would quickly go out of memory trying to compute its range and all subset sums of that. As a minimal evidence of correctness, the function g in the link instead takes an explicit limit on its right and n on its left, and computes all numbers up to limit that are non-sums of distinct n-th powers.

Deriving a correct and simple upper bound

A001661 contains a formula for the upper bound of the last term for each n. I was able to simplify it in many steps (and greatly overshoot the bound in the process) to just 7 bytes of APL:

d*2^(n-1)*(c*2^n + (2/3)*d*(4^n - 1) + 2*d - 2)^n + c*d
where c = n!*2^(n^2) and d = 2^(n^2 + 2*n)*c^(n-1) - 1

= d*(c + 2^(n-1)*(c*2^n + (2/3)*d*(4^n + 2) - 2)^n)

= d*(2*c + (2*c*2^n + (4/3)*d*(4^n + 2) - 4)^n)/2

by 2c<n^n*(2^n)^n = (n*2^n)^n
and a^n+b^n<(a+b)^n
= d*(n*2^n + 2*c*2^n + (4/3)*d*(4^n+2)-4)^n/2

by n*2^n + 2*c*2^n < c*(2^(n+1)+1) < 2^(n^2+2*n)*c <= d + 1
<= d*(d + 1 + (4/3)*d*(4^n+2)-4)^n/2

= d*(d*(4^(n+1)/3+11/3)-3)^n/2

by 4^(n+1)/3+11/3 < c
and 2^(n^2 + 2*n) = 2^(n^2)*4^n <= 2^(n^2)*n!*64/6
which gives d < c^n*64/6
< c^n*64/6 * (c^n*64/6 * c)^n/2
= c^(n^2+2n) * (32/3)^(n+1)/2

by (32/3)^(n+1)/2 < c^2
(for n=2 = 16384/27 < 607)
(for n=3 < 6475)
< c^(n^2+2n+2)
< c!

The APL formula for `c!`: !!×2*×⍨

The code

⊢{z~∊(⊢,+)\⍺*⍨z←⍳!⍵}!×2*×⍨    Monadic train. ⍵←n
                    !×2*×⍨    c = n!*2^(n^2)
⊢{                 }          Pass to a dyadic dfn, ⍺←n, ⍵←above
              z←⍳!⍵   z ← 1..c!
           ⍺*⍨        Raise each to the power of n
    ∊(⊢,+)\           Compute all sums of distinct powers
  z~                  Remove them from z

Answer 2

Ruby, 80 ... 72 71 bytes

->n,*w{*r=0;a=1;1while w!=w|=[*(z=a**n)..(a+=1)**n]-r|=r.map{|y|y+z};w}

Try it online!

Explanation:

The algorithm stops when no new numbers are found between (x-1)^n and x^n. When that range is covered, then by just adding x^n we can cover the range up to 2*x^n. I assume that's enough to reach (x+1)^n (but alas, I can't prove it, sorry)

Answer 3

Python 2, 97 bytes

n=input()
q=x=1
w=0
a=n**n<<n*n
while q<a**a:
 if w**n<q:w+=1;x|=x<<w**n
 if~x>>q&1:print q
 q+=1

Try it online!

Uses \$\left(n^n2^{n^2}\right)^{n^n2^{n^2}}\$ as an upper bound (Derived from Bubbler's upper bound with \$n^n\ge n!\$). Calculating this takes so much time that there won't be any output for \$n>4\$ on TIO.

---

Here is an alternative approach that doesn't rely on an upper bound but rather terminates when no new number was found in \$\left(w^n, (w+1)^n\right]\$ for some integer \$w\$. I'm not sure if this is valid, but the highest numbers from \$n=2\$ to \$n=6\$ (see below) match up with the paper linked in the challenge, and I think the Ruby answer works in a similar way.

def f(n):
 w=f=q=1;x=3
 while f:
  w+=1;x|=x<<w**n;f=0
  while w**n^q:f|=x>>q&1<1!=print(q);q+=1

Try it online!

Here is a more optimized version of the same algorithm that only prints the last integer: Try it online!

Answer 4

Jelly, 19 bytes

²2*×!!ŒṗQƑƇfƑƇɗÐḟ²€

Try it online!

A monadic link taking an integer n as its argument and returning a list of integers. Implements @Bubbler’s method for determining the upper bound, so be sure to upvote that one too! While this answer would end in a finite time, it would be a very very long time even for n=2, and requires a huge amount of RAM (enough to generate integer partitions of \$3072!\$). I’ve therefore included a footer on TIO which demonstrates that the rest of the code works for integers up to 50.

Answer 5

Python 3.8 (pre-release), 202 bytes, (Compact™)

g=lambda x:x<2or x*g(~-x);from itertools import*
def f(n,q=2):
 while q<g(g(n)*2**n**n):k=[x**n for x in range(1,q)];q not in[sum(h)for l in range(1,-~len(k))for h in combinations(k,l)]and print(q);q+=1

Try it online!

Hmph, the damn upper limit is too big. Won't work for \$n>2\$ because of the big factorial calculation

-1 thanks to @cairdcoinheringaahing

Python 3 (PyPy), 283 bytes, (Performance™)

g=lambda x:x<2or x*g(~-x);from itertools import*
def f(n,q=2):
 c=g(n)*2**n**2;d=2**(n*n+2*n)*c**~-n-1
 while q<d*2**~-n*(c*2**n+.6*d*~-4**n+2*d-2)**n+c*d:k=[x**n for x in range(1,-~int(q**(1/n)))];q not in[sum(h)for l in range(1,-~len(k))for h in combinations(k,l)]and print(q);q+=1

Try it online!

Probably terminates more than billion times faster compared to the above version.

Optimizations:

• Uses the full formula to calculate the upper bound from A001661, the generated upperbound is still big, but not as big as one by the Compact version
• To check if number can be represented in nth powers of 1 or more integers, previous version checks in the range of 1 to the number, but this checks until the nth root of the number, which speeds up the check dramatically.

It works for \$n>2\$ too!