Project Euler Number 338

Question

I'm stuck on Project Euler problem 338. Here is what I've done so far...

Let's denote a rectangle with width and height x and y respectively (x,y). To form new rectangles you can consider cutting a kind of stairway down the diagonal (as is shown in in the problem description) with d steps. But to form a new rectangle the following must hold: d|x and either (d-1)|y or (d+1)|y. The new rectangle then becomes either (x/d*(d-1), y/(d-1)*d) or (x/d*(d+1), y/(d+1)*d). Obviously the new rectangles area is the same as the old rectangle's.

That was enough to confirm that G(10)=55 and G(1000)=971745 by looping through all relevant d and adding all new rectangles to a set being careful to count (x,y) and (y,x) only once.

The main issue with this method is that it's possible to make a new rectangle in two different ways. For example, (9,8) can transform into both (6,12) and (12,6) with d=3 and either d-1 or d+1 dividing y. Or another example of (4,4) transforms into both (2,8) and (8,2) with d=2 and d=1 respectively.

I was then lucky enough to read this blog post. It removes the need to check for duplicates by searching for one of the sides instead.

def F(w, h):
    if w&1 and h&1: return 0
    if w<h: w,h = h,w

    r = 0
    x = 1
    while x**2 <= w*h:
        if (w*h)%x!=0 or x==h:
            x += 1
            continue

        if w%(w-x)==0 or x%(x-h)==0:
            r += 1

        x += 1

    return r

def G(N):
    s = 0
    for w in range(1, N+1):
        for h in range(1, w+1):
            s += F(w,h)

    return s

G(10¹²) would require far too long to solve regardless of how fast F is though. I think it's necessary to either use some kind of sieving algorithm where we loop through all x < 10¹² counting how many (w,h) satisfy h <= w <= 10¹², x|(w*h), x != h and (w-x)|w or (x-h)|x.

I think an O(n^2/3) algorithm must be possible... but I'm stuck here!

---

Edit: I don't have access to the forum since I am unable to solve it. That's why I'm asking for help. I have completed most other questions and want to tackle this question now!

Edit 2: I think considering the areas by prime factors is a dead end. That's because there are 10²⁴ different areas. Rectangles with prime areas have 0 solutions, rectangles with semiprime areas have 1 solution if one of the primes is 2 otherwise they have 0 solutions. But counting all the semiprime solutions alone would take too long since we'd need to count all the primes p such that 2*p < 10²⁴ which is not feasible.

Edit 3: I've stripped down the code:

def G(N):
    s = 0
    for x in range(1, N):
        for h in range(1, N+1):
            if x==h: continue
            for w in range(max(h, x**2//h), N+1):
                if (w*h)%x==0 and x%(w-x)==0 and x%(x-h)==0:
                    s -= 1

    for x in range(1, N):
        for h in range(1, N+1):
            if x==h: continue
            for w in range(max(h, x**2//h), N+1):
                if (w*h)%x==0 and w%(w-x)==0:
                    s += 1

    for x in range(1, N):
        for h in range(1, N+1):
            if x==h: continue
            for w in range(max(h, x**2//h), N+1):
                if (w*h)%x==0 and h%(x-h)==0:
                    s += 1

    return s

I don't think breaking the brute-force code down will work though. Remember it's enough that we just count the solutions (x, w, h) to each of these three subproblems. The last such summation would have the constraints 0 < x < N, 0 < h < N+1, x!=h, max(h, x²/h) < w < N+1, x|wh and x-h|h.

I think we should start with the assumption that some prime p divides x, w, h or even x-h and then see what we can deduce about the other variables. If that works well, maybe consider p^k for arbitrary k.

Answer 1

I don't have a solution yet, but something interesting for Python. I realized that Python can be used as a convenient tool for notation of algorithms! Basically I wrote down a program similar to yours and started transforming the program logically which leave the results unchanged. I came up with

def order(x,y):
    if x>=y:
        return (x,y)
    else:
        return (y,x)

N=1000
num=set()
for n in range(1, N+1):
    for a in range(1,N//n+1):
        for b in range(1,N//(n+1)+1):
            if a==b: continue
            num.add((order(a*n,b*(n+1)), order(b*n,a*(n+1))))

print(N, len(num))

Obviously brute force or even a simple loop over 10^12 is not feasible, but maybe with this algorithm one can find a closed form expression at some point. If it wasn't for the set character of num, it would be doable. Maybe one can find duplicate point this way.

This could be a dead end, but still it's pretty cool that Python can be used for notation and work with algorithms :)

Any progress on your side?