Проект Эйлера номер 338

Я застрял в Project Euler задаче 338. Вот что я сделал до сих пор...

Обозначим прямоугольник шириной и высотой x и y соответственно (x,y). Чтобы образовать новые прямоугольники, вы можете сократить своего рода лестницу по диагонали (как показано в описании задачи) с d ступенями. Но для формирования нового прямоугольника должны выполняться следующие условия: d|x и либо (d-1)|y, либо (d+1)|y. Затем новый прямоугольник становится либо (x/d*(d-1), y/(d-1)*d), либо (x/d*(d+1), y/(d+1)*d). Очевидно, что площадь нового прямоугольника такая же, как у старого прямоугольника.

Этого было достаточно, чтобы подтвердить, что G(10)=55 и G(1000)=971745, перебрав все соответствующие d и добавив все новые прямоугольники в набор, стараясь подсчитать (x,y) и (y,x) только один раз.

Основная проблема с этим методом заключается в том, что новый прямоугольник можно создать двумя разными способами. Например, (9,8) может превратиться как в (6,12), так и в (12,6) с d=3 и либо d-1, либо d+1, разделив y. Или другой пример преобразования (4,4) в (2,8) и (8,2) с d=2 и d=1 соответственно.

Затем мне посчастливилось прочитать эту запись в блоге. Это устраняет необходимость проверять наличие дубликатов, вместо этого выполняя поиск одной из сторон.

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¹²) потребовало бы слишком много времени для решения, независимо от того, насколько быстра F. Я думаю, необходимо либо использовать какой-то алгоритм просеивания, в котором мы перебираем все x ‹ 10¹², подсчитывая, сколько (w,h) удовлетворяют h ‹= w ‹= 10^{12, x|(w*h), x != h и (wx)|w или (xh)|x.}

Я думаю, что алгоритм O(n^2/3) должен быть возможен... но я застрял здесь!

Изменить: у меня нет доступа к форуму, так как я не могу его решить. Вот почему я прошу о помощи. Я ответил на большинство других вопросов и хочу ответить на этот вопрос сейчас!

Редактировать 2: я думаю, что рассмотрение областей по простым факторам — это тупик. Это потому, что существует 10²⁴ различных областей. Прямоугольники с простыми площадями имеют 0 решений, прямоугольники с полупростыми площадями имеют 1 решение, если одно из простых чисел равно 2, в противном случае они имеют 0 решений. Но подсчет только всех полупростых решений занял бы слишком много времени, так как нам нужно было бы подсчитать все простые числа p, такие что 2*p ‹ 10²⁴, что невозможно.

Редактировать 3: я убрал код:

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

Я не думаю, что взлом кода грубой силы сработает. Помните, что достаточно просто подсчитать решения (x, w, h) каждой из этих трех подзадач. Последнее такое суммирование будет иметь ограничения 0 ‹ x ‹ N, 0 ‹ h ‹ N+1, x!=h, max(h, x²/h) ‹ w ‹ N+1, х|ш и хч|ч.

Я думаю, мы должны начать с предположения, что некоторое простое число p делится на x, w, h или даже x-h, а затем посмотреть, что мы можем сделать для других переменных. Если это сработает, можно рассмотреть p^k для произвольного k.