Тогда, во-первых, любое a_k по модулю n является квадратом известного нам вычета (m+k).

Но, во-вторых (и это ключевой шаг), можно попробовать найти произведение a_k, которое будет полным квадратом.

А именно — предположим, что у какой-то их части оказались одни и те же простые делители. Тогда, чтобы произведение каких-то a_{k_j} было полным квадратом, достаточно их перемножить так, чтобы все степени простых оказались чётными.

Понятно, что если мы раскладываем на множители, скажем, 100-значное число n — то a_k будут порядка корня из n, то есть 50-значными. И раскладывать каждое из них на множители это опять та ещё задача.

Но — мы выберем некоторый набор Base "не слишком больших" простых (который называется базой ). Скажем — первую сотню тысяч, или миллион.

И если какое-то a_k не раскладывается в произведение простых из Base — то просто выбросим его. А если раскладывается — запомним (и назовём такое число " гладким ").

И когда мы, просто увеличивая k, наберём достаточно много гладких чисел (больше, чем простых в базе Base) — то вектора степеней разложения (координаты соответствуют простым из базы) начнут быть зависимыми по модулю 2. И произведение соответствующих a_{k_j} будет полным квадратом, и мы будем знать, чьим.

Вот мы и получили совпадение: построенное произведение a_{k_j} будет и настоящим полным квадратом в обычных целых числах (так что A — обычный корень из этого произведения), и будет по модулю n сравнимо с квадратом произведения (m+k_j) — соответственно, B это произведение (m+k_j).

Игрушечный пример: факторизуем число 64643.

Оказывается, что взяв три небольших k=2,3,7, мы уже видим (mod 2)-линейную зависимость между векторами степеней:

Находим A и B:

И находим делители —

Да, нам может не повезти. Никто не обещал, что первое же найденное соотношение окажется нетривиальным. Я бы (безответственно и навскидку) сказал, что шансы 50 на 50: можем найти нетривиальное соотношение, а можем тривиальное. Поэтому гладких чисел нужно находить не просто |Base|+1, а немного побольше (чтобы хоть одно соотношение да сработало)...

Второе — что в базу можно брать не все простые, а только примерно половину. Потому что если простое число p делит a_k=(m+k)^2-n, то
n=(m+k)^2 mod p,
поэтому n является квадратичным вычетом по модулю p.

Ответ на вопрос "является ли n квадратичным вычетом по модулю p" называется символом Лежандра.
И его вполне можно быстро найти. Вручную я, скорее всего, воспользовался бы мультипликативностью и квадратичным законом взаимности. Но для программирования, пожалуй, проще всего сказать, что символ Лежандра a над p равен a^{(p-1)/2} (ибо малая теорема Ферма). А k-ю степень мы вычисляем за log(k) умножений, ибо череда возведений в квадраты.

Да, в качестве рекламы — про квадратичный закон взаимности (связывающий символы Лежандра p по модулю q и q по модулю p) мне очень нравится рассуждение отсюда — https://users.mccme.ru/smirnoff/papers/qrecip01.pdf

Но так или иначе, символ Лежандра быстро ищется — так что можно быстро понять, какие простые у нас не могут появиться в принципе. Собственно, именно поэтому повторение 2, 19, 37, 47 в примере выше не столь удивительно (хотя нам и повезло — уже три вектора в четырёхмерном пространстве оказались линейно зависимыми — но куча других простых в первых десятках тут не появилась вообще).

Ну и последнее — это выбор размера базы. Возьмём слишком маленькую — и слишком мала будет доля гладких чисел, придётся перебирать много разных k. Скажем, много ли чисел вида
2^{p_1}*3^{p_2}*5^{p_3}
среди первого миллиона?

Возьмём слишком большую базу — и нужно будет много гладких чисел, опять-таки, придётся много перебирать.

Но на формулировке этой задачи оптимизации я, пожалуй, хочу в этом рассказе на сегодня остановиться.