QFACT - Quick FACtorization

алгоритм быстрого поиска целых делителей целого числа

 

Автор категорически запрещает использовать приведенный алгоритм для нанесения ущерба кому бы то ни было.

 

Постановка задачи:

Пусть имеется некоторое целое положительное число N такое, что оно является произведением двух других целых положительных чисел P и Q. Зная N найти P и Q.

 

Алгоритм с полиноминальным законом затрат времени от N искать нет смысла, поскольку для очень больших чисел время работы алгоритма будет физически нереализуемым. В качестве возможной основы рассмотрим алгоритм «ловли льва в пустыне» с логарифмическим законом затрат времени от N. Суть этого алгоритма состоит в следующем. Пусть в пустыне есть лев. Огородим пустыню по периметру забором, чтобы лев не убежал. Лев находится в огороженном участке. Перегородим участок забором на две примерно равные части. Если льва нет в одной части, то он находится в другой. Перегородим ту часть, в которой находится лев, на две примерно равные части забором и снова посмотрим, в какой из частей находится лев. Так будем делать до тех пор, пока  лев не окажется в маленьком загоне. Для алгоритма быстрой факторизации необходимо не только научиться огораживать участки, но еще и узнавать в каком из участков находится лев. Научимся это делать.

Задача факторизации или разложения целого числа на множители является НП-полной задачей. Из теоремы Гёделя следует, что НП-полнота проявляется как следствие невозможности разрешения логических противоречий в рассуждениях. Эти противоречия влекут за собой теоретическую неразрешимость задачи в рамках заданных аксиом и ограничений, определяющих некое пространство или его подпространство. Из теоремы Гёделя также следует, что полное описание исследуемой системы высказываний, а значит и разрешение всех противоречий, возможно только в рамках более мощной системы высказываний, включающей в себя (объемлющей) исследуемую систему высказываний. Cледовательно, для решения задачи быстрой факторизации необходимо так расширить систему высказываний, описывающей алгоритм факторизации, т.е. так расширить пространство, отменив некоторые аксиомы или ограничения, чтобы НП-полная задача стала разрешимой или даже тривиальной. Правильность этого вывода подтверждалась неоднократно. Классический пример с пятым постулатом геометрии Эвклида и геометриями Лобачевского и Римана.

В качестве иллюстрации необходимости снятия ограничений и понимания механизма возникновения НП-полноты задачи рассмотрим известный логический парадокс про приказ командира полка, в котором имеется один единственный брадобрей. Командир издал приказ своим подчиненным о запрете ходить небритыми и одновременно о запрете бриться самостоятельно. Неразрешимая проблема возникла у брадобрея: он не может ходить небритым и одновременно не может сам себя брить, не нарушив приказ командира. Неразрешимая, казалось бы, проблема брадобрея сразу разрешается, как только будет снято ограничение единственности брадобрея в полку. Пусть в полку будет хотя бы еще один брадобрей. Тогда брадобреи будут брить друг друга, все подчиненные будут побриты и приказ командира не будет нарушен. В расширенном полку противоречий у брадобреев уже не возникает.

Необходимо снять еще одно очень важное ограничение – ограничение целочисленности. Рассмотрим постановку задачи вычисления квадратного корня: пусть имеется некоторое число N такое, что оно является произведением двух других чисел P и P. Как видно, постановка этой задачи практически полностью совпадает с постановкой задачи факторизации числа N. Но, как исходное число, так и искомое при вычислении квадратного корня далеко не всегда являются целыми числами. На основе совпадения постановок задач можно предположить, что алгоритм вычисления квадратного корня является частным вырожденным случаем общего алгоритма разложения чисел на множители. Следовательно:

1. если существует алгоритм для общего случая, то существует и алгоритм для частного случая, т.к. алгоритм для частного случая может быть легко получен путем замены вырожденных участков графа полного алгоритма простыми дугами;
2. если алгоритм для частного случая не всегда дает точное решение, например, для случаев, когда целое число P не является точным квадратом целого числа N, то этим же свойством должен обладать и алгоритм общего вида.

 

Теперь можно приступить к описанию теоретической части алгоритма быстрой факторизации. Для краткости изложения будем опускать доказательства вполне очевидных свойств, которыми будем пользоваться при построении алгоритма. Вопросы теории полиномов, позиционных систем счисления, форматов представления кодов чисел и ошибок округления разрядной сетки для краткости изложения так же рассматривать не будем. Ежедневная практика доказывает их справедливость.

 

Свойство №1: если произведения двух пар положительных чисел равны и меньшее число первой пары меньше  меньшего числа другой пары, то большее число первой пары больше большего числа второй пары.

 

Это свойство устанавливает над числами A,B,C,D порядок A<C<D<B, если A*B=C*D и A<B, C<D, A<C. Из этого свойства следует, что и для среднего геометрического M=C=D=sqrt(A*B) двух чисел A и B справедливо A<M<B.

Введем операцию выравнивания кода или просто выравнивания. Под выравниванием кода числа X к числу N по основанию a будем понимать умножение или деление на a числа X до тех пор, пока выровненное число X' не будет удовлетворять условию aⁿ<N,X'<aⁿ⁺¹. Интервал [aⁿ:aⁿ⁺¹] назовем интервалом выравнивания, число a назовем основанием выравнивания, число n назовем порядком выравнивания, число X' назовем выровненным. Операция выравнивания имеет самостоятельный смысл приведения мантиссы числа к виду, необходимому для представления числа в формате с плавающей запятой в системе счисления по основания a. В дальнейшем будем применять термин мантисса для выровненных значений чисел. Всегда можно выполнить как прямое проецирование любого числа на интервал выравнивания, так и выполнить обратное проецирование.

Снимем ограничение единственности делителей P и Q числа N. Для этого в качестве делителей возьмем еще обычно игнорируемые из-за их тривиальности числа 1 и N.

Теперь все готово к доказательству очень важной теоремы. Эта теорема является ключом к построению алгоритма.

 

 Теорема «о загоне со львом»:

Выровненные значения одной из пар чисел P и Q или sqrt(a)*P и sqrt(a)*Q лежат на интервале [aⁿ:N], если для выровненного значения числа sqrt(N) справедливо aⁿ*N=sqrt(N)*sqrt(N) и aⁿ<sqrt(N)<N, или на интервале [N:aⁿ⁺¹], если для выровненного значения числа sqrt(N) справедливо aⁿ⁺¹*N=sqrt(N)*sqrt(N) и N<sqrt(N)<aⁿ⁺¹.

Доказательство:

Пусть имеется целое положительное число N, которое является произведением двух целых положительных чисел P и Q. Возьмем число a>1. Найдем порядок выравнивания n для числа N по основанию a. Найдем выровненное на интервал [aⁿ:aⁿ⁺¹] значение числа sqrt(N). Пусть a такое, что справедливо aⁿ*N=sqrt(N)*sqrt(N) и aⁿ<sqrt(N)<N для выровненного значения числа sqrt(N). Тогда, согласно свойству №1, или выровненные значения пары чисел P и Q должны удовлетворять aⁿ<P<sqrt(N)<Q<N для выровненных значений чисел, или выровненные значения пары чисел sqrt(a)*P и sqrt(a)*Q должны удовлетворять aⁿ<sqrt(a)*P<sqrt(N)<sqrt(a)*Q<N для выровненных значений чисел. Появление числа sqrt(a) не случайно, поскольку выровненные значения произведений sqrt(a)*sqrt(N)*sqrt(a)*sqrt(N)=a*N и sqrt(a)*P*sqrt(a)*Q=a*P*Q=a*N будут равны выровненному значению произведения sqrt(N)*sqrt(N)=P*Q=N. И одна из пар чисел P и Q или sqrt(a)*P и sqrt(a)*Q обязательно должна попасть в интервал, иначе не будет выполняться свойство №1. Что и требовалось доказать. Аналогично доказываем для случая aⁿ⁺¹*N=sqrt(N)*sqrt(N) и N<sqrt(N)<aⁿ⁺¹.

 

Важным следствием из теоремы «о загоне со львом» является существование интервалов мантисс чисел, не содержащих выровненные значения возможных нетривиальных делителей P и Q числа N, если существует такое основание системы счисления a, что выполняются условия теоремы. Однако, поскольку в теореме принимает участие число sqrt(a), то на выравнивающем интервале появляются два интервала, один из которых содержит выровненные значения пары чисел P и Q, а другой содержит выровненные значения пары чисел sqrt(a)*P и sqrt(a)*Q. Границы второго интервала вычисляются путем умножения границ основного интервала на sqrt(a). Если одна из границ интервала выравнивания попадает внутрь дополнительного интервала, то дополнительный интревал разбивается на два подинтервала по линии границы выравнивающего интервала и значения границ получившихся интервалов выравниваются. Например, если границы дополнительного интервала будут sqrt(a)*N<aⁿ⁺¹<sqrt(a)*aⁿ⁺¹, то появятся подинтервалы [aⁿ:sqrt(a)*aⁿ] и [sqrt(a)*N:aⁿ⁺¹]. Основной и дополнительные интервалы могут пересекаться и покрывать весь выравнивающий интервал. Если покрывается не весь выравнивающий интервал, то в оставшихся интервалах не могут находиться делители P и Q числа N. Поэтому, в алгоритме нас будут интересовать только те основания a, для которых можно выделить интервалы не содержащие делителей. Интервалы, содержащие делители или их произведения на sqrt(a) будем называть «загонами».

Убедиться в существовании хотя бы одного числа N, удовлетворяющего условиям теоремы можно посмотрев текст программы и протокол ее работы, приведенные в приложениях №1 и №2. Кроме того, как видно из приложения №3, выполнение условий теоремы «о загоне со львом» не является редким явлением и для других оснований, отличных от основания 2. В качестве исходного числа в примерах было взято недавно разложенное на множители 575-разрядное двоичное число. Аналогичная картина, но с обратным порядком кодов мантисс в двоичной системе счисления, наблюдается для числа из первого сообщения об алгоритме RSA. Во время работы над одной из версий алгоритма, работающего на основе поразрядного уравновешивания, было замечено, что коды мантисс чисел, которые брались для проверки работы алгоритма, всегда были либо больше, либо меньше кода числа N. Анализ такого поведения кодов чисел 1, N, P, Q и sqrt(N) привел к доказательству теоремы «о загоне со львом».

Порядок следования кодов мантисс выровненных чисел определяется только взаимными свойствами чисел a и N. Поскольку операция выравнивания – это всего-навсего операция нормализации мантисс чисел в формате с плавающей запятой в позиционной системе счисления с основанием a, то число 2 не может обладать особыми свойствами по сравнению с другими числами, которые так же могут быть взяты в качестве основания позиционной системы счисления и выравнивающего интервала.

Цель данной работы показать принципиальную возможность быстрого и малозатратного нахождения делителей числа N и как это делать, поэтому и не рассматриваются остальные варианты взаимного расположения мантисс чисел sqtr(N), N и N². Приведенное ниже описание алгоритма сделано для наглядности и простоты понимания в ущерб оптимальности алгоритма.

Для определения «загона», т.е. интервала со всеми делителями числа N, необходимо один раз вычислить sqrt(N), алгоритм вычисления которого является составной часть общего алгоритма и остается в качестве самостоятельного алгоритма вычисления квадратного корня при отбрасывании остальной части всего алгоритма.

 

Построим базовую процедуру «огораживания» целых делителей числа N. Учтем, что если число N можно представить n-разрядным числом в системе счисления с основанием a, то нужно найти не более n границ интервалов на основе теоремы «о загоне». Объясняется это тем, что с помощью операции выравнивания на интервал aⁿ:aⁿ⁺¹ проецируются интервалы 1:a, a:a², a²:a³, … a^n-1:aⁿ, aⁿ:aⁿ⁺¹ и, следовательно, обратная операция отображения «загонов» на интервал 1:N будет осуществляться делением на a границ найденных выровненных «загонов» соответствующее число раз. Для различных a нас будут интересовать только совместные интервалы, поскольку только они могут содержать внутри себя сами числа P и Q, а не проекции этих чисел. Функция «огораживания» encl имеет следующие параметры:

1.      N – число, для которого ищутся делители,

2.      a – основание выравнивания,

3.      encs – массив границ интервалов.

 

Функция однократного «огораживания»:

^{1.        Находим числа M=sqrt(N) и N2=N2}

^2.        Находим число n такое, что aⁿ<N2<aⁿ⁺¹

^3.        Выравниваем числа M и N на найденный интервал [aⁿ:aⁿ⁺¹]

^4.        Если нельзя однозначно определить интервал с делителями P и Q согласно теоремы «о загоне со львом», то прерываем работу функции encl и возвращаем массив encs

^5.        Находим границы «загона» по теореме «о загоне со львом»

^6.        Находим обратные проекции «загона» на интервалы 1:a,…,a^n-1:aⁿ

^7.        Находим совместные интервалы из encs и 1:a,…,aⁿ:aⁿ⁺¹

8.      Возвращаем уточненный список границ совместных «загонов»

 

Прежде чем приступить к описанию алгоритма, необходимо сделать несколько замечаний. Если в качестве нормы выравнивания брать степени одного числа, то вполне возможно, что не произойдет отсева лишних «загонов», либо отсев будет незначительным. Поэтому, в качестве норм выравнивания следует брать различные числа и/или их степени. Например, простые числа 2, 3, 5, 7, 11 и т.д. Однако, как представляется, может оказаться вполне достаточным взять лишь числа вида 2ⁱ+1. Размеры полученных целочисленных «загонов» могут быть значительными и простой перебор чисел из этих «загонов» может занять очень большое время. Для значительного уменьшения затрат времени на перебор чисел и физической нереализуемостью вычислительных затрат рассмотрим число (2ⁱ+1)/2ⁱ, где i=1,…int(n/2). При максимальном значении i=int(n/2) дробная часть (2ⁱ+1)/2ⁱ примерно равна 1/2ⁱ, т.е. ширине младшего двоичного разряда числа 2ⁱ. Поэтому алгоритм должен состоять из двух частей: целочисленной предварительной и рациональной окончательной. Теперь все готово для построения самого алгоритма быстрой факторизации или быстрого поиска целых делителей целого числа.

 

QFACT или алгоритм  быстрого поиска целых делителей целого числа:

1.      Найдем число k такое, что  1<N<2^k-1

2.      Cоздадим начальное множество «загонов» encs=[1:N]

3.      Для всех i=k-1,…,0 выполним отсев «загонов» encs=encl(N,2ⁱ+1,encs)

4.      Для всех i=1,…,k/2 выполним отсев «загонов» encs=encl(N,((2ⁱ+1)/2ⁱ),encs).

5.      Проверим на делимость числа N все целые числа, оставшиеся в «загонах», границы которых меньше sqrt(N).

 

Шаг №5 алгоритма является реализацией признака целочисленности делителей числа N.  Не обязательно до конца выполнять шаг №4. Скорее всего, нужно искать компромисс между ростом вычислительных затрат и простым перебором чисел на шаге №5. При вычислениях степеней числа a необходимо пользоваться алгоритмом поразрядного уравновешивания, сложность вычислений которого имеет логарифмический закон. В итоге, вычислительные затраты приведенного алгоритма растут по закону квадрата логарифма числа N, а при использовании предварительно вычисленных значений степеней чисел a будут расти медленнее квадрата логарифма и в идеале должны приближаться к логарифму в первой степени. Затраты на разложение 1000-разрядного числа на множители такие: нужно выполнить 1000 шагов №3 и 500 шагов №4. При выполнении шага №3, когда i=1, получится не более 1000 интервалов. Если начинать работу алгоритма в обратном порядке числа i, т.е. с широких «загонов», то число «загонов» будет всегда небольшим, если только число N не состоит из более чем двух простых нетривиальных сомножителей. Та же картина с интервалами будет наблюдаться и на шаге №4 при возрастании i. При выполнении шага №4 «загоны», границы которых имеют малые значения, будут постоянно отсеиваться. Нет смысла строить все «загоны». Нужно строить только те «загоны», которые являются совместными с оставшимися «загонами» из encs. Нет необходимости каждый раз производить все вычисления всех степеней числа a, поскольку они могут быть вычислены ранее. Для реализации алгоритма необходимо переписать библиотеку работы с числами в формате с плавающей запятой с управляемой разрядностью мантиссы и порядка. Или написать более оптимально сам алгоритм.

При работе рациональной части алгоритма появляется ошибка округления результата из-за ограниченности разрядной сетки. Если не ограничивать разрядную сетку, то сложность вычислений с полной разрядной сеткой может стать физически нереализуемой. Преодолевается это следующим образом. В случае, когда младшие разряды кода мантиссы очередного произведения вытесняются за пределы разрядной сетки, производится расщепление результатов на две последовательности чисел, одна из которых используется для нахождения границы aⁿ, а другая для нахождения границы aⁿ⁺¹, для которых всегда производится прибавление единицы в младший разряд мантиссы. Таким образом, числа в первой последовательности всегда будут меньше aⁿ, а числа во второй последовательности будут всегда больше aⁿ⁺¹. После нахождения границ «загонов» производится деление чисел на aⁱ из той же последовательности. Благодаря такому расщеплению степеней числа a, обратные проекции точных границ «загонов» всегда будут находиться внутри границ неточных «загонов».

По предварительным и достаточно грубым оценкам, основанным на времени работы модели целочисленной части алгоритма на языке bc, для оптимальной реализации алгоритма, в которой не выполняются повторяющиеся вычисления, как и в алгоритмах Ли, Винограда и быстрых сверток, на компьютере с процессором P6 с тактовой частотой 2.6 ГГц для 1000-разрядного двоичного числа время поиска целочисленных делителей должно составлять порядка 9-15 секунд.

Умозрительно приведенный алгоритм можно представить следующим образом. Возьмем стеклянную пластинку и нанесем на нее черной краской числовую ось. Отметим на ней другой краской числа log(1), log(a), log(N). Разобьем интервал log(1):log(N) на интервалы шириной log(a), начиная с точки log(1). Число log(N) должно находится внутри последнего интервала. Найдем логарифмы границ интервалов, содержащих все проекции делителей числа N на последний интервал по теореме «о загоне». Сотрем краску внутри «загонов». Скопируем на остальные интервалы границы «загонов» из последнего интервала и сотрем краску внутри границ этих «загонов». Возьмем еще одну пластинку и проделаем с ней ту же процедуру, но уже для другого числа a. Наложим одну пластинку на другую и посмотрим через них на свет. Останутся узкие пропускающие свет полоски, а некоторые интервалы совсем исчезнут. Если взять много пластинок, то останется малое число очень узких прозрачных полосок. Получилась своеобразная муаровая решетка. Эта решетка и явилась той идеей, которая была положена в основу алгоритма. Голографический метод нахождения делителей, основанный на муаровых решетках, оказывается гораздо более эффективным, чем решето Эратосфена.

 

Южно-Сахалинск-Москва

1999 - 2002 - 14 октября 2005

Юрий Аверьянов (Urix)

 

Copyright © Yuri I.Averyanov (Urix), 2005. All rights reserved.

Приложение №1

#!/usr/bin/bc -q

 

define f50 (a)

{

    auto r;

    r = 2^50;

    while (a/2 >= r) { a /= 2 }

    return (a);

}

 

define round (a, e)

{

    while (a/2 > e) { a /= 2 }

    return (a);

}

 

p  = 398075086424064937397125500550386491199064362342526708406385189575946388957261768583317

q  = 472772146107435302536223071973048224632914695302097116459852171130520711256363590397527

 

print "P = ",p,"\n"

print "Q = ",q,"\n"

print "N = ",p*q,"\n"

print "\n#############################################\n\n"

obase=2

print "P = ",p,"\n"

print "Q = ",q,"\n"

print "N = ",p*q,"\n"

print "\n#############################################\n\n"

obase=10

 

n1 = p*q

e  = 1;

n  = 0;

while (e < n1) { e *= 2; n += 1 };

e /= 2; n -= 1;

r  = round(sqrt(2*e*e),e)

m1 = round(sqrt(p*q*e*e),e)

 

print "n = ",n,"\n\n"

obase = 2

 

print "\n#############################################\n\n"

print "E1  = ",f50(e+1),"\n"

print "N^2 = ",f50(n1^2),"\n"

print "R   = ",f50(r),"\n"

print "N   = ",f50(n1),"\n"

print "P   = ",f50(p),"\n"

print "M   = ",f50(m1),"\n"

print "Q   = ",f50(q),"\n"

print "E2  = ",f50(e*2-1),"\n"

 

print "\nWhere M=sqrt(N), R=sqrt(2)*2^n\n"

 

quit

 

Приложение №2

$ ./qfact7.bc

P = 3980750864240649373971255005503864911990643623425267084063851895\

75946388957261768583317

Q = 4727721461074353025362230719730482246329146953020971164598521711\

30520711256363590397527

N = 1881988129206079638386972394616504398071635633794173827007633564\

22988859715234665485319060606504743045317388011303396716199692321205\

734031879550656996221305168759307650257059

 

#############################################

 

P = 1100110011101001010101000101011111110001001001111100010010010101\

01000110111001010111100001000001000000101001101011000110101011110111\

00000110000001001100011001001110100000101000000100000001100010101100\

10111011010100111000001000110011101110110101011111110011011110100110\

11100110100010010101

Q = 1111001101011100101110000010010111110000000000110011100101101111\

00010001111100010100100110100101100110100000111111100111001111011011\

01000011000101100010110010110011000001011011100011101011111111101101\

10100101110010101000010000001011111000100100110001010011011010110110\

10101010111001010111

N = 1100001011001011101100100100111111011011111110010010001110110110\

00010010011010001110001111110001000110100011100010010110110111100100\

01010111010010110011101110100101100001110011000011001011110101100101\

00101001001110001000011001001110001000100010001111101110111010110111\

00000100101000010111110011111101000010001101000101101011010001101000\

10010001101001100001010001110100011101011001100100111001110001101110\

01001001101010101111111001111111001001011001010101010100100011000111\

01001100000111010111111110111000110100100100110011010001010111001011\

001000111011010011001101000010100011

 

#############################################

 

n = 575

 

#############################################

 

E1  = 100000000000000000000000000000000000000000000000000

N^2 = 100101000011100101011100010101010101000000111111100

R   = 101101010000010011110011001100111111100111011110011

N   = 110000101100101110110010010011111101101111111001001

P   = 110011001110100101010100010101111111000100100111110

M   = 110111110100111110001110000111010010111010001000101

Q   = 111100110101110010111000001001011111000000000011001

E2  = 111111111111111111111111111111111111111111111111111

 

Where M=sqrt(N), R=sqrt(2)*2^n