Функция на pari/gp, которая перечисляет все подходящие прямоугольники для отрезка длины n, или возврашщает пустой вектор если подходящих нет. Авторство Qwen
rect(n) = {
\\ Проверка входных данных
if(n < 1 || n != floor(n), return("Ошибка: n должно быть натуральным числом"));
\\ 1. Вычисляем правую часть уравнения
my(P = 2^(n+1) + (n-1)^2);
\\ 2. Целевая четность для делителей x и y (должна совпадать с четностью n-1)
my(target_parity = (n - 1) % 2);
\\ 3. Инициализируем список для хранения решений
my(sols = List());
\\ 4. Получаем отсортированный вектор всех натуральных делителей числа P
my(divs = divisors(P));
\\ 5. Перебираем делители
for(i = 1, #divs,
my(x = divs[i]);
\\ Условие: x >= n + 1 (чтобы m >= n)
if(x < n + 1, next());
\\ Условие: x <= y (чтобы не получать симметричные дубликаты на этом этапе)
\\ Как только x^2 > P, мы прошли квадратный корень, можно прерывать цикл
if(x^2 > P, break());
my(y = P / x);
\\ Условие четности: и x, и y должны иметь ту же четность, что и n-1
if((x % 2) != target_parity || (y % 2) != target_parity, next());
\\ 6. Вычисляем m и k
my(m = (x + n - 1) / 2);
my(k = (y + n - 1) / 2);
\\ Добавляем пару (m, k)
listput(sols, [m, k]);
\\ Если m != k, добавляем также симметричную пару (k, m)
if(m != k,
listput(sols, [k, m]);
);
);
\\ 7. Преобразуем список в вектор и сортируем его для красивого вывода
my(res = Vec(sols));
return(vecsort(res, [1, 2]));
}
Запуск:
Код:
? rect(43)
time = 1 ms.
%43 = [[3754, 1178153386], [7438, 592968406], [18686, 235630694], [31790, 138438326], [37106, 118593698], [158866, 27687682], [27687682, 158866], [118593698, 37106], [138438326, 31790], [235630694, 18686], [592968406, 7438], [1178153386, 3754]]
?
Нормально работает до n около 200, потом тормозит если не повезло (факторизация чисел порядка
).
-- добавлено через 13 минут --И даже есть размеры в которых равенство.
Это, кажется, ничего не меняет. Ну есть и есть.
