Построение быстрых алгоритмов вычисления наибольшего общего делителя пар натуральных чисел и связанных алгоритмов тема диссертации и автореферата по ВАК РФ 05.13.11, кандидат наук Аль-Халиди Аркан Мохаммед Али

• Аль-Халиди Аркан Мохаммед Али
• кандидат наук
• 2020

Оглавление диссертации кандидат наук Аль-Халиди Аркан Мохаммед Али

ВВЕДЕНИЕ

ГЛАВА 1. ОБЗОР АЛГОРИТМОВ ВЫЧИСЛЕНИЯ НАИБОЛЬШЕГО ОБЩЕГО ДЕЛИТЕЛЯ

1.1. Классический алгоритм Евклида

1.2. Расширенный алгоритм Евклида

1.3. Бинарный алгоритм вычисления НОД

1.4. k-алгоритм вычисления НОД

1.5. Аппроксимирующий Л-арный алгоритм

1.6. Новая стратегия вычисления параметров х и у при малых значениях |а|

1.7. Обзор результатов Амера по программированию Л-арного алгоритма

1.8. Выводы по главе

ГЛАВА 2. РАЗРАБОТКА МОДЕЛИ ТЕСТИРОВАНИЯ РАЗЛИЧНЫХ АЛГОРИТМОВ ВЫЧИСЛЕНИЯ НОД

2.1. Анализ вычислительных возможностей языков C и Python

2.2. Программирование аппроксимирующего Л-арного алгоритма

2.3. Реализация аппроксимирующего Л-арного алгоритма на языке Python

2.4. Экспериментальные результаты вычислений на языке Python

2.5. Анализ полученных результатов

2.6. Выводы по главе

ГЛАВА 3. РЕАЛИЗАЦИЯ АЛГОРИТМОВ ВЫЧИСЛЕНИЯ НОД НА ЯЗЫКЕ C++

3.1. Методика экспериментальных вычислений

3.2. Экспериментальные результаты вычислений для различных Л

3.2.1. Случай для к =

3.2.2. Случай L = 500 при увеличивающих значениях к

3.3. Сравнение результатов для реализации на Python и C++

3.4. Выводы по главе

ГЛАВА 4. РАЗРАБОТКА АЛГОРИТМА ВЫЧИСЛЕНИЯ ОБРАТНЫХ ПО МОДУЛЮ ЭЛЕМЕНТОВ

4.1. Решение уравнения Безу с использование расширенного алгоритма Евклида

4.2. Вычисление обратных по модулю элементов с использованием k-арного алгоритма

4.3. Примеры вычисления обратных элементов по схеме k-арного алгоритма

4.4. Оценка производительности расширенного Л-арного алгоритма

4.4.1. Оценка частоты вариантов вычисления обратного элемента

4.4.2. Выводы по вариантам выпадения для пар различной длины

4.4.3. Оценка времени вычисления обратного элемента

ЗАКЛЮЧЕНИЕ

СПИСОК ИСПОЛЬЗОВАННЫХ ИСТОЧНИКОВ

СПИСОК СОКРАЩЕНИЙ И УСЛОВНЫХ ОБОЗНАЧЕНИЙ

ПРИЛОЖЕНИЯ

ВВЕДЕНИЕ

Современные криптографические системы работают с длинными числами и выполняют многочисленные рутинные операции арифметических вычислений и встроенные процедуры обработки данных. Одной из наиболее распространенных операций, входящей в различные числовые алгоритмы, является процедура вычисления их наибольшего общего делителя (далее -НОД) целых чисел [16]. Для вычисления наибольшего общего делителя чаще всего используется классический алгоритм Евклида. Однако известны и другие алгоритмы, общей целью создания которых было стремление уменьшить время вычисления НОД, сократить время его поиска.

Одним из самых перспективных таких алгоритмов является парный алгоритм поиска НОД, опубликованный в 1990-х годах Джонатаном Соренсоном [104]-[106]. По сравнению с классическим алгоритмом Евклида, данный алгоритм позволяет заметно ускорить время вычислений при аккуратном выборе параметра к. В дальнейшем самим Д. Соренсоном и другими учёными было предложено несколько вариантов улучшения и модификации к-арного алгоритма ([116]-[119]).

Одним из них является аппроксимирующий к-арный алгоритм, предложенный Ш.Т. Ишмухаметовым [78]- [82] Этот алгоритм позволяет добиться значительного ускорения процедуры вычисления НОД для использования в алгоритмах, работающих с длинными целыми числами. Достоинством нового алгоритма является значительное сокращение общего числа итераций при вычислении наибольшего общего делителя двух чисел, что приводит к значительному сокращению времен работы алгоритма.

Одним из недостатков работы к-арного алгоритма является появление так называемых «ложных факторов». Ложными факторами называются дополнительные множители, входящие в состав НОДа, вычисленного по схеме парного алгоритма, но не входящие в НОД исходных чисел. Для

удаления ложных факторов необходимо выполнить дополнительное вычисление с помощью алгоритма Евклида. Эта операция, хотя и занимает дополнительное время, но значительно меньше, чем полное время работы к-арного алгоритма. В работах Д. Долгова [61], [62], Еникеева Р.Р. [11]-[15], Sedjelmaci [109] и других авторов исследуются алгоритмы ускорения k-арного алгоритма и проблемы поиска и удаления ложных факторов.

В статье [60] Диксоном получены теоретические оценки среднего числа итераций в алгоритме Евклида. Эти оценки были экспериментально проверены и подтверждены в кандидатской диссертации Амера [4]. В работе исследовательницы из Франции Бриджит Вале [114] результаты Диксона были усилены, убрана ссылка на недоказанную обобщенную проблему Гильберта и получены средние оценки числа итераций для нескольких вариантов алгоритма Евклида, бинарного алгоритма и разных вариантов к-арного алгоритма.

В работах Ишмухаметова Ш.Т., Аль Халиди А., Амера И. [16], [1]-[4] изложены разные аспекты программирования аппроксимирующего к-арного алгоритма.

Учебники и учебные пособия [5],[6],[9],[19]-[30] содержат основные понятия и определения, связанные с рассматриваемой проблемой. В работах иностранных исследователей [38]-[50] содержатся результаты, относящиеся к теории вычислений НОД для натуральных чисел и полиномов.

Следующий цикл статей [51]-[123] посвящен известным работам русских и иностранных исследователей по различным аспектам теории алгоритмов, в которых исследуются проблемы, для которых проблема вычисления НОД входит как часть решения. В частности, вычисления точек эллиптических кривых, тесты простоты, вычисления в конечных полях, проблема Гильберта и другие важные математические проблемы основаны на алгоритмах, включающих в себя вычисление НОД как важный элемент.

Поэтому эти статьи также относятся к проблемной области нашего исследования.

Одной из важных проблем при выполнении вычислений с длинными числами является скорость выполнения операций вычисления НОД для целых чисел и элементов конечных полей. Для эффективной работы необходимо быстрое выполнение всех элементарных составляющих алгоритмов вычисления НОД - это подбор параметра к, выбор области изменения всех параметров, способы представления данных в памяти компьютера, эффективное программирование всех делателей, составляющих алгоритм.

В диссертации Исмаила Амера [4] разработана тестовая компьютерная модель вычисления наибольшего общего делителя и использованием ^ арного алгоритма и аппроксимирующего к-арного алгоритма. Им предложено целая серия улучшений и приемов ускорения этих алгоритмов, позволившая добиться ускорения вычисления НОД целых чисел размерности до 3000 бит (такие числа используются в современной криптографии) до 4-5 раз по сравнению с исходными реализациями, предложенными Соренсоном, Вебером и Ишмухаметовым.

Также им был построен и реализован алгоритм поиска строго псевдопростых чисел, основанный на аппроксимирующем алгоритме вычисления НОД [77]. Этот алгоритм работает примерно в пять раз быстрее своего аналога, основанного на схеме Евклида. Этим была доказана практическая ценность использования аппроксимирующего алгоритма в вычислениях с длинными числами.

Одной из проблем реализации к-арного алгоритма Вебера и Соренсона является ухудшение его работы при больших значения коэффициента р=А/В, где А и В - числа, поданные на вход очередной итерации алгоритма. Если его значение превышает значение параметра к, то следующая итерация выполняется намного хуже, чем в случае вычисления

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

Аппроксимирующий алгоритм обеспечивает значительное ускорение по сравнению с Л-арным алгоритмом Соренсона именно из-за того, что при больших значениях коэффициента следующая пара будет

опять иметь большое значение коэффициента , что позволяет аппроксимирующему алгоритму работать быстрее.

Целью нашей работы явилась дальнейшее ускорение аппроксимирующего алгоритма Ишмухаметова путем анализа теоретических и экспериментальных составляющих, а также разработка расширенной версии алгоритма, аналогичной расширенному алгоритму Евклида, для ускорения вычисления обратных по модулю элементов для нужд теории чисел и криптографии.

Одним из недостатков аппроксимирующего алгоритма является некоторое его замедление при значении ключевого параметра , высчитываемого на каждой итерации алгоритма, близком к Мы сумели преодолеть этот недостаток путем использования двух конкурирующих стратегий при реализации аппроксимирующего алгоритма. Это позволило увеличить значение коэффициента редукции для числа

вычисляемого на очередной итерации алгоритма с теоретической оценки р >3к/4 до оценки р > 1 , что позволило добиться 25-процентного ускорения работы алгоритма.

Теоретическая идея использования аппроксимирующего алгоритма для вычисления коэффициентов Безу в уравнении А и + ВV = сС принадлежит авторам статьи «Вычисление коэффициентов Безу для Л-арного алгоритма нахождения НОД», опубликованной в журнале Известия ВУЗов. Математика в 2017 году [79]. В нашей работе приведено полное исследование этого

7

алгоритма, разработаны шесть вариантов его реализации при различных сочетаниях свойств подаваемых на вход пар чисел, выполнена практическая реализация алгоритма и выполнены практические вычисления, доказывающие значительное ускорение общей работы процедуры поиска коэффициентов Безу и обратных по модулю элементов.

Для всех рассмотренных в диссертации алгоритмов были приведены подробные теоретические описания и обоснования, выведены оценки сложности и описана реализация на языках программирования Python и С. Для реализации арифметики длинных чисел использована библиотека длинных чисел MPIR [9]. Выполнено тестирование каждого из алгоритмов по времени выполнения на числах различной длины. Выполнены экспериментальные вычисления и разработана методика выбора оптимальных параметров, обеспечивающих быструю сходимость аппроксимирующего Л-арного алгоритма. Результаты экспериментов собраны в таблицы. Построены графики оценки скорости выполнения алгоритмов и числа итераций вычисления НОД для пар чисел различной длины.

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

Соответствие паспорту специальности.

Работа выполнена в рамках направления области исследований «Модели, методы и алгоритмы проектирования и анализа программ и программных систем, их эквивалентных преобразований, верификации и тестирования» паспорта специальности 05.13.11 - Математическое и

программное обеспечение вычислительных машин, комплексов и компьютерных сетей.

Введение диссертации (часть автореферата) на тему «Построение быстрых алгоритмов вычисления наибольшего общего делителя пар натуральных чисел и связанных алгоритмов»

Апробация работы.

Результаты, полученные автором, неоднократно представлялись на семинарах кафедры системного анализа и информационных технологий КФУ, на итоговых научных конференция КФУ, на 3-й Международной конференции - Брикс 2019 (Иннополис, июль 2019).

Публикации автора исследования.

По теме НКР автором подготовлены 9 публикаций, из которых 3 входят в список Scopus и Web of Science, и 6 - в список РИНЦ.

Научная новизна.

В работе предложены новые методы ускорения аппроксимирующего Агарного алгоритма: разработана модификация аппроксимирующего алгоритма с использованием двух конкурирующих стратегий, выполнена реализация модифицированного алгоритма на Python и C++, выполнены теоретические и экспериментальные оценки скорости сходимости алгоритма в зависимости от длины исходных чисел и значения параметра , исследовано влияние на скорость алгоритма различных предвычислений, в том числе таблиц обратных элементов по модулю , таблиц рядов Фарея.

Разработана эффективная реализация алгоритма Ишмухаметова-Мубаракова - Аль Анни вычисления коэффициентов Безу, построен быстрый алгоритм вычисления обратных по модулю элементов, превосходящий аналоги на числах до 3000 бит примерно в 5 раз.

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

Методы исследования.

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

Практическая ценность работы.

Разработаны комплексы программ реализации различных алгоритмов вычисления НОД на языках программирования Python и C++. Выполнены экспериментальные вычисления и разработана методика выбора оптимальных параметров, обеспечивающих быструю сходимость аппроксимирующего к-арного алгоритма. Результаты экспериментов собраны в таблицы. Построены графики оценки скорости выполнения алгоритмов и числа итераций вычисления НОД для пар чисел различной длины. Доказана важность и практическая значимость использования предтаблиц в процедуре вычисления НОД с помощью аппроксимирующего к-арного алгоритма.

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

ГЛАВА 1. ОБЗОР АЛГОРИТМОВ ВЫЧИСЛЕНИЯ НАИБОЛЬШЕГО ОБЩЕГО ДЕЛИТЕЛЯ

В этой главе мы выполним обзор основных алгоритмов вычисления наибольшего общего делителя НОД для пар натуральных чисел. Нами будут подробно изучены три основных алгоритма нахождения НОД - алгоритм Евклида, А-арный алгоритм Соренсона и аппроксимирующий алгоритм Ишмухаметова.

1.1. Классический алгоритм Евклида

Алгоритм вычисления наибольшего общего делителя (НОД) был открыт древнегреческими математиками и известен как алгоритм "взаимного вычитания". И хотя еще Аристотель упоминал об этом алгоритме, его приписали Евклиду.

На вход алгоритма подаются натуральные числа A и B, где А > В > 0 . Классический алгоритм Евклида поиска НОД (кратко, КАЕ) основывается на следующем рекуррентном равенстве:

GCD{A,B) = GCD{B,Amod В). (аббревиатура GCD для обозначения НОД является сокращение от английского Greatest Common Divisor и будет далее использоваться в тексте диссертации).

Алгоритм выполняется по шагам, называемым итерациями, в ходе каждой из которых вычисляется остаток , происходит переход к

новой паре и вычисление продолжается с новой парой. Вычисление

заканчивается, когда второй аргумент пары становится равным 0, тогда первый аргумент и есть искомый НОД:

ипл г ил Í а> если Ь = 0 Í Л 14

На языке программирования С++ процедура вычисления НОД имеет

вид:

int gcd(int A, int B) {

if (B == 0) return A; return gcd(B, A % B);

}

Оценивая сложность данного алгоритма, заметим, что за каждые две итерации цикла делимое уменьшается, по меньшей мере, в два раза. Это означает, что количество итераций есть О (п) , где п - длина входных чисел. Наихудшая оценка достигается на соседних числах ряда Фибоначчи, а среднее число итераций оценивается величиной (Dixon [60], Vallee [114]):

12

Number of iteration---r log2 B,

TC

где B - меньшее из заданных чисел (А , В ) . Длина А не влияет на асимптотическую оценку числа итераций, поскольку за одну итерацию происходит переход к новой паре ( где

Операцию вычисления остатка А тоd В по схеме (1.1.) можно выполнить за время О (п1 о g 2п) , так что общая производительность алгоритма оценивается величиной . Обозначим через

исходные числа, тогда процесс вычисления НОД по схеме Евклида представляет собой вычисление последовательности {А ¿}f= 0 , где очередной член последовательности вычисляется по схеме

ft 1 К 7')^)

где - целая часть от деления первого множества пары

(А ', А ' _ -,_) на второе. Обозначим через Q ' матрицу ( ^ Ц '), тогда переход к

паре от исходной пары (А выполняется по формуле

ft 1 ) = Q 'Q '_ - ■ -Qitt)

через цепочку матричных вычислений. Значения целых чисел q ' являются

небольшими, и по замечанию Дональда Кнута [18] более, чем в 70% случаев

12

меньше или равны 3. Поэтому значение ц^ на ьй итерации может быть вычислено путем анализа нескольких старших цифр чисел А ^ и АI_г. Поэтому почти во всех случаях можно выполнить такой анализ на несколько итераций вперед, не выполняя полного деления чисел А ^ на А ^ _ ±, и перемножив соответствующие матрицы @ ¿, выполнить переход на несколько итераций за одно вычисление

Эта идея лежит в основе алгоритма Лемера (D. Lehmer), позволившего ускорить алгоритм Евклида до 2-х и более раз. Сам Лемер не публиковал свой алгоритм, но его описание можно найти в разделе 4.5 тома 2 монографии Д. Кнута [18].

Еще одно приложение алгоритма Евклида - это построение цепной дроби, приближающей произвольное действительное число а (см. [16], [23])

Процесс построения такой дроби для рационального числа представляет собой ту же процедуру Евклида, и процесс построения цепной дроби - это процесс вычисления последовательности частных ц ¿, заканчивающийся при достижении значения цп, равному НОД исходных чисел. Это произойдем, когда соответствующее значение Ап окажется равным 0. Если ограничиться частичным вычислением, то на ьй итерации мы может получить приближение числа а подходящей дробью а^:

1

а = qt-\--

1

а ~ aL = qt -\--

Ч2

Точность аппроксимации будет выше при увеличении номера ¿. Отметим, что идея аппроксимации действительного или рационального

13

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

1.2. Расширенный алгоритм Евклида

Алгоритм Евклида можно расширить для нахождения по заданным a и b таких целых x и у, что

(1.2)

где d - наибольший общий делитель a и b. Числа x и у называются коэффициентами Безу, а само уравнение (1.2) уравнением Безу.

Лемма. Пусть для положительных целых чисел a и b (a > b) известны d = Н О Д( a, b) = Н О Д ( b, а то d b) , а также числа х ' и у ' , для которых:

d = х' • b + у' • (a mod b)

Тогда значения x и у, являющиеся решениями уравнения ах + b у = d, находятся из соотношений (2):

х = у ',у = х' - у' * | (1.3)

Через |_x J здесь обозначена целая часть числа x.

Функция , приведенная

ниже, по входным числам a и b находит d = НОД^, b) и коэффициенты Безу x, у. Для поиска неизвестных x и у необходимо рекурсивно запустить функцию GCDext(b, a mod b, d, x, у) и пересчитать значения x и у по выше приведенной формуле. Рекурсия заканчивается, когда b = 0. При b = 0 НОД^, 0) = a и a = a * 1 + 0 * 0, поэтому полагаем x = 1, у = 0.

Этот алгоритм может быть записан на языке программирования С++ следующим образом:

void gcdext(int a,int b, int *d, int *x, int *y){ int s; if (b == 0) {

*d = a; *x = 1; *y = 0; return;

}

gcdext(b,a % Ь^,х,у);

б =

у;

*х = б;

*у = *

х

:а / Ь)

у);

}

Для ручного выполнения расширенного алгоритма Евклида удобно воспользоваться таблицей с четырьмя столбцами (как показано ниже в примере), соответствующих значениям а, Ь, х, у. Процесс вычисления разделим на три этапа:

а) Сначала вычислим НОД(а, Ь), используя первые два столбца таблицы и соотношение (1). В последней строке в колонке а будет находиться значение й = НОД(а, Ь).

б) Занесем значения 1 и 0 соответственно в колонки х и у последней строки.

в) Будем заполнять последовательно строки для колонок х и у снизу вверх. Значения х и у в последнюю строку уже занесены на этапе б). Считая, что в текущей заполненной строке уже находятся значения (х', у'), мы при помощи формул (2) будем пересчитывать и записывать значения (х, у) в соответствующие ячейки выше стоящей строки.

Пример вычисления. Найдем решение уравнения 5 х + 3 у = 1 . Вычисление НОД(5, 3) и нахождение соответствующих х, у произведем в таблице 1.1:

Таблица 1.1. Нахождение коэффициентов х, у

а)

а В х у

5 3 -1 2

3 2 1 -1

2 1 0 1

1 0 1 0

в)

Порядок и направление вычислений указаны стрелками и буквами. Из таблицы 1 находим, что , то есть:

х = -1, у = 2. 15

*

Рассмотрим, например, как получены значения (х, у) для первой строки. Со второй строки берем значения (х' ,у' ) = ( 1 , — 1 ) . По формулам (1.3) получим:

х = у' = -1,

у = х - у

= 1 - ( " 1 )

= 1 + 1 = 2 .

Мы вернемся к этим формулам повторно в главе 3, где будем строить расширенную схему для ^-арного алгоритма.

1.3. Бинарный алгоритм вычисления НОД

Этот алгоритм был разработан в 1967 израильским программистом Р.Стейном [110]. Входными данными алгоритма являются нечётные натуральные числа. Обозначим их и, V. Вычисление по схеме бинарного алгоритма (рисунок 1.1) можно представить следующим образом:

while (u!= 0 && vl=0} f if (u %2==0) u-u/2; else if {v%2—0}

else t- abs{u - v)/2; if (u > v) then u=t else v=t; if {u = = 0) then t=v else t=u return t; }

Рисунок 1.1. Бинарный алгоритм

За каждую итерацию одно из исходных чисел уменьшается по меньшей

мере в два раза, поэтому число итераций этого алгоритма ограничено сверху

бинарной длиной большего из чисел исходной пары. По оценке Диксона [60]

оно сравнимо с числом итераций в классическом алгоритме Евклида, но

немного превышает его. Однако, из-за того, что операции, выполняемые на

шаге алгоритма, имеют оценку (вычисление разности и циклический

сдвиг), то алгоритм имеет асимптотическую оценку В целом,

16

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

Расширенная версия бинарного алгоритма для вычисления коэффициентов Безу была разработана М. Пенком и приведена в фундаментальной монографии Р.Крэнделла и П.Померанса ([19], с.522).

1.4. к-алгоритм вычисления НОД

Обобщая бинарный алгоритм, получим парный алгоритм поиска НОД. к-арный алгоритм был разработан в 1990 году Д.Соренсоном [104]-[106] и улучшен независимо друг от друга Вебером [117],[118] и Джебелеаном [84]-[86] в 90- годах 20 века. В основе этого алгоритма лежит теорема, доказанная Соренсоном.

Теорема (Соренсон [104]). Пусть заданы целые числа А > В > 0 , и к-

положительное натуральное число, взаимно простое с . Существуют целые числа х и у, удовлетворяющие соотношению |х| + |у| < 2 V к такие, что выполняется соотношение:

Ах + Ву= 0 то( к. (1.4)

Доказательство этой теоремы будет приведено ниже. Сам к-арный алгоритм работает следующим образом. Выберем значение к, равным степени 2, к = 2 5. Это позволит выполнять быстрее текущие операции с длинными числами типа сдвигов, сравнений, вычисления обратных элементов по модулю к и других быстрее, чем в случае других значений .

Получив на входе шага итерации два нечетных числа А и Б, находим целые числа и , удовлетворяющие теореме Соренсона. Полагая

, получим новую пару , имеющую такой же НОД, что и

НОД либо НОД, кратный НОД . Посторонние множители,

появляющиеся в НОД(В/С) (ложные факторы) должны быть изъяты на

последней стадии алгоритма путем дополнительного вычисления НОД исходных А и В с НОД, полученным как выход к-арного алгоритма.

Определяя , мы также должны проверить условие

(С = Н О Д ( С; к) ф 1 , и если оно не выполняется, то требуется поделить С на это .

Отметим, что уравнение (1.4) эквивалентно уравнению:

у = -цх тос( к, где ц = АВ~ 1 тоСС к (1.5)

Искомую пару х, у можно найти, перебирая значения х от 1 до V к и вычисляя два значения у из интервала [—V к -.V к] , равные у = -цх то с( к и у = к -(цх тосС к). Вся процедура вычисления очередной пары (А,В) , сводится к следующим вычислениям:

1. Сначала вычислим . Отметим, что в силу выбора к вычисление остатков по модулю к представляет собой просто операцию отделения младших разрядов чисел А и В длины .

2. Далее вычисляем

3. В цикле по переменной от 1 до вычисляем два значения

у = - цх то с( к и у = к - (цх то СС к ) до тех пор, пока не будет найдено значение пары , удовлетворяющей условиям теоремы Соренсона.

4. После этого вычисляет | | . Если полученное число является четным, то дополнительно сокращаем на двойку до тех пор, пока

оно не станет нечетным.

5. После этого можно перейти от пары к новой паре

В статье [100] Д.Соренсон предложил выполнять к-арный алгоритм в четыре этапа: этап предварительных вычислений, первое «пробное» деление, основная часть, второе «пробное» деление.

Этап предварительных вычислений: если основание алгоритма к уже

выбрано, то для всех натуральных чисел, меньших к, алгоритм

предусматривает вычисление их наибольшего общего с числом к делителя и

18

поиск обратного элемента по модулю к (соответствующие таблицы обозначим буквами G, I).

Помимо этого составляется специальная таблица A для быстрого вычисления чисел x и y и таблица для хранения некоторых малых общих делителей чисел u и v, которые отсекаются на стадии пробного деления. Все эти вычисления не зависят от переданных алгоритму чисел u и v, а потому могут быть выполнены заранее и только один раз.

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

Далее выполняется основной цикл, в ходе которого вычисляется число М . Это число не обязательно равно искомому НОД, но равно их кратному, то есть может содержать посторонние множители. Чтобы исключить их, выполняется второе пробное деление.

Второе пробное деление: эта стадия алгоритма является завершающей и использует результат, полученный после выполнения основной части алгоритма. Чтобы избавиться от посторонних множителей, делим на все малые делители, находящиеся в заранее составленной таблице, а затем умножается на результат первого пробного деления. В итоге получим НОД заданных чисел u и v.

Наиболее значимая часть алгоритма заключается в главном цикле, который представляет основную стадию алгоритма. На языке программирования C++ он может быть записан следующим образом:

int KARY(int u, int v, int k){

int u1, v1; int a, b, x;

while ((u != 0) && (v != 0)){

u1 = (int)(u % k); v1 = (int)(v % k); q=u1*Inv(v1,k);

Find(u,v,x,y,q,k);// Searh for (x,y) t = abs(a * u + b * v) / k; if (u > v) u = t; else v = t; }}}

}

Некоторые этапы вычисления дублируются на многих итерациях, например, вычисление обратных элементов, поэтому для ускорения общего вычисления возможно использование предварительно вычисленных таблиц. Следует, однако, отметить, что варианты, предложенные Д.Соренсоном, не всегда являются возможными в силу того, что размер таблиц становится таким большим, что превышает реальный выигрыш от их использования. Например, по расчетам Амера [4] хранение таблиц пар чисел (х , у) при занимает сотни террабайт и становится практически нереализуемым.

При реализации к-арного алгоритма первой проблемой встает способ выбора целых чисел х и у,удовлетворяющих уравнению (1.4). Согласно теореме Д.Соренсона, коэффициенты х и у можно выбрать 0 < |х| + |у| <

24k.

Для максимального быстродействия алгоритма выгодно также выбрать число x малым и положительным, число y достаточно большим отрицательным.

Сложность к-арного алгоритма для чисел и и v длиной п бит в оптимальном случае выбора параметра к равна, согласно оценке Д.Соренсона, О (п Л 2 /1 о g/) . Однако, эта оценка приведена для плавающего значения параметра /, которое выбирается по методологии Д.Соренсона примерно равным , где длина входных чисел. По мере уменьшения

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

1.5. Аппроксимирующий Л-арный алгоритм

Описание аппроксимирующего к-арного алгоритма (кратко, АКА) дано в статье Ш.Т.Ишмухаметова [78]. Этот алгоритм является естественным развитием -арного алгоритма Соренсона и отличается в способе выбора параметров х,у, удовлетворяющих соотношению (1.4). Рассмотрим кратко схему этого алгоритма.

Пусть даны натуральные числа Для простоты предположим,

что числа А,В - нечётные. Выберем параметр к равным степени двойки к = 2Вычислим а = А тоСС к, Ь = В тоСС к. Теперь по аналогии с к-арным алгоритмом Соренсона требуется найти такие числа х и у такие, что выполняется

Ах + Ву = 0 тос( к, а потом вычислить и заменить пару новой парой

(В ,С) .

Искомые значения коэффициентов выбираются так, чтобы

минимизировать сумму | | Очевидно, что для этого необходимо,

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

новый способ вычисления чисел используя приближённое значение величины г, вычисляемой по формуле:

г г — q^

а =

где г = А/В, ц=А В~ 1то сС к, скобки \.. .] обозначают целую часть числа. Рассмотрим процедуру поиска коэффициентов более подробно.

Обозначим через положительное значение . Согласно

(1.5) у = — цх (тосС к) или у = — цх + кБ для некоторого б Е X. Тогда:

С =

Ах+Ву

= В

rx-qx+ks

x + s

Обозначим через (3 дробь (г — q)/ к, получим

С = B\(x + s\,s EZ.

Пусть (( = а + ss, где s s - целое число, выбранное так, чтобы выполнилось |а\ < 0 ,5 . Тогда,

С = B\ax + x • ss + s\,s Е Z (1.6)

Наша задача состоит в том, чтобы подобрать параметры x, 1 <x <к, и целое s так, чтобы \ ax + x • ss + s \ приняло наименьшее возможное значение.

Для решения этой задачи будем использовать две возможные стратегии, первая из которых была базовой и разработана в статье [78] и реализована в кандидатской диссертации Амера [4], а вторая - разработана в этой работе. Рассмотрим сначала первую стратегию:

1. Сначала вычисляется приближенное значение г = А /B с точностью до 1 /к. Отметим, что в языках высоко уровня типа Java, Python, нет операторов приближенного вычисления отношений, поэтому здесь можно было бы выполнить только полное деление, которое занимает значительную часть общего времени вычисления на отдельной итерации. Нами в совместной статье [2] был разработан алгоритм приближенного вычисления этого параметра, позволивший добиться ускорения выполнения

аппроксимирующего алгоритма до полутора раз ([4], глава 3) только за счет ускорения вычисления этого параметра.

2. Потом вычисляется д, удовлетворяющее (1.5).

3. Далее вычисляем // =

r—q

, находим а и ищем рациональную

дробь т/п, со знаменателем ограниченным сверху значение /, приближающую а.

4. Для аппроксимации а используем дроби Фарея [23]. Дробь Фарея т /п представляет собой правильную дробь с числителем и знаменателем, меньшими величины к. Предположим, что удалось найти

т

аппроксимирующую дробь — « а с точностью г:

£ =

т

--а

п

тогда полагая x равным n и s равным - т — п • s s получим:

|ах + х • ss + s| = + г + ss) n + ss = |т + гп + п • ss + s| = |гп|. (1.7)

Лемма. Если аппроксимировать действительное число а из интервала:

4 к -2"

\к' k-l)

с точностью г < 1 / к2, поэтому выражение (1.7) не превысит значения 1 / к.

Рациональное приближение для параметра а можно получать разными способами, например, используя цепные дроби или ряды Фарея. Проблема построения наилучшего приближения изучена Р.Еникеевым [15]. В нашем случае наилучшее приближение обеспечивается использованием дробей Фарея. Рассмотрим процедуру построения аппроксимирующей дроби более подробно.

Определение 1. Пусть задана константа С>0. Дробями Фарея называется совокупность правильных рациональных дробей со знаменателем, ограниченным С.

Определение 2. Пусть даны две рациональные дроби г1=т1/п1 и г2=т2/п2. Медианой этих дробей называется дробь г, равная (т1+т2)/(п1+п2). Медиана г удовлетворяет условию г1<г<г2.

Пусть выбрана константа к> 1, и задано действительное число аЕ(0; 1/2). Будем строить правильную дробь со знаменателем, меньшим к, обеспечивающую наилучшее приближение а. Предположим, что

В качестве первого приближения возьмем пару соседних дробей вида (1/(п+1); 1/п), содержащую а. Сделаем индукционное предположение о том,

что на шаге I найдены две правильные дроби:

т1 т2

гг = — < а < гп — Щ п2

Проверим условие п1+п2к. Если оно выполнено, то завершаем вычисление, выбрав в качестве решения одну из дробей т1/п1 и т2/п2.

Иначе, построим медиану г=(т1+т2)/(п1+п2). Проверим условие г а. Если оно выполнено, то переходим к новому интервалу (г1; г), содержащему а. Иначе, перейдем к интервалу (г; г2). Переходим к шагу /+1.

Процесс конечен, так как на каждой итерации происходит увеличение знаменателя одной из ограничивающих дробей.

Приведем пример построения дроби Фарея для а = 0, 278:

ае(?4И?;7И^)=(0; 272:0 286):

Процедура заканчивается на третьем шаге, поскольку следующая медиана имеет знаменатель п=11+7=18>к. Наилучшим приближением для а=0, 278 является дробь 3/11, поскольку расстояние от левой границы до а равно примерно 0,006, а от а до правой границы 0,08. Однако с точки зрения производительности алгоритма такой выбор будет неверным. Действительно,

коэффициент р1 определен, как величина обратная произведению расстояния на значение знаменателя. Для левой дроби 3/11 это значение равно 1 / (0,00 6 • 1 1 ) = 1 5 ,2 . Для правой дроби 2/7 значение рх = 1/ (0,008 • 7 ) = , значит, выбор правой границы является более эффективным.

Найдем теоретическую границу между двумя дробями, для которой значение коэффициента р1 одинаково для левой и правой границ:

Значит, точка между равновесия между двумя дробями, на которой значения коэффициента редукции совпадают и равны медиане граничных точек. Оценим наименьшее (наихудшее значение) коэффициента редукции р^а) на интервале (1/(к-1);1/2).

Предположим, что значение а находится между двумя соседними дробями т1/п1 и т2/п2 такими, что п1 + п 2 < /.

Лемма 1. Расстояние между соседними дробями Фарея т1=п1 и т2=п2 равно 1/(п1+п 2).

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

Предположим, что на шаге I соседние дроби т -/п 1 и т 2 /п2, содержавшие а, удовлетворяют условию леммы, т.е. выполняется - т 2п 1 = 1. Следующий интервал, содержащий а, имеет вид ( т1=п1; ( 7т + т2)=( п1 + п2)) или (( т-,^+т2 )=( п1+ п 2); (т2=п2)). Тогда, /т1+т2 т-^щ + т2п1—т1п1—т1п2 \ 1

т1+т2

илипЛа — тЛ = т7 — ап7 а =-

Щ+п2

Ущ + п2 ^ п^щ + п2)

Аналогично, для правого интервала. Лемма доказана.

п^щ+п^'

Лемма 2. Значение коэффициента редукции при приближении а, удовлетворяющего условию |а|>1/к дробями Фарея со знаменателем, ограниченными к—1, больше или равно к.

Доказательство. Точка равновесия между дробями т1/п1 и т2/п2 равна медиане г этих дробей. В наихудшем случае значение а совпадает с этой медианой и расстояние от нее до границ интервала равно согласно лемме 1 ^1=1/п1(п1+п2) и d2=1/n2(n1+n2) соответственно. Соответствующие значения коэффициента редукции равны между собой и равны:

р1(а) = 1/(п^) = 1/(п^2) = п1 + п2 > к.

Последнее неравенство вытекает из того факта, что дроби т1/п1 и т2/п2 является конечными в процедуре построения дробей Фарея, поэтому сумма их знаменателей больше или равно к. В частности, при а=1/к значение р1(а) равно к. Лемма доказана.

При значениях < 1 / к, ближайшей к |а| дробью Фарея будет дробь 1 / ( к — 1), тогда в наихудшем случае, когда а~0, расстояние до 1/(к-1) приближается к 1/(к-1), и значение р1(а) окажется близким к 1, что означает плохую сходимость алгоритма на этой итерации.

Чтобы сгладить негативные последствия подобных случаев, определим вторую стратегию.

1.6. Новая стратегия вычисления параметров х и у при малых значениях |а|

При малых значениях |а| первая стратегия дает сбой, ухудшающий эффективность аппроксимирующего алгоритма, поэтому будем здесь использовать вторую стратегию. Теперь мы не станем выполнять аппроксимацию а дробью т=п, а определим ,

тогда C= B|а|, и коэффициент редукции p2(а)= B/C примет значение p2(а) =1/|а|.

При уменьшении la | значение p2(a) увеличивается, а значение pi(a) уменьшается. Найдем точку равновесия, когда эти значения совпадают (для простоты предположим, что a>0).

Учтем, что при a< 1/(k— 1) s=1/(k—1)—a и m/n=1/(k—1): pi (A,B)=p2(A,B) ^ s n=a ^ (1/(k—1)—a) . (k— 1)=a; откуда, a=1/k. Значит, стратегия 2 дает лучшую сходимость при lal<1=k. При |a|=1/k обе стратегии дают одинаковое значение коэффициента редукции p(A,B)=k.

Ниже в главе 3 мы дадим описание процедуры поиска дроби т/ п и реализации этого алгоритма на языках C++ и Python.

По теореме об асимптотической сложности аппроксимирующего k-арного алгоритма количество итераций на каждом этапе оценивается 0 (п 1 о g 2 п), где п - длина исходных чисел в битах. Асимптотическая оценка времени работы всего алгоритма дается величиной:

Список литературы диссертационного исследования кандидат наук Аль-Халиди Аркан Мохаммед Али, 2020 год

СПИСОК ИСПОЛЬЗОВАННЫХ ИСТОЧНИКОВ

Работы автора по теме диссертации Публикации Scopus

1. Al Halidi Arkan M. A Study Of Euclidean K-Ary GCD Algorithm / N.A.Antonov, S.T.Ishmukhametov, Al Halidi Arkan M., M.E. Maiorova // Revista Publicando.-2017.-Vol.4, №13 (1).

2. Аркан Аль Халиди Мухаммед. Эффективное программирование процедуры вычисления НОД натуральных чисел // A.М.Аркан, Ш.Т.Ишмухаметов / Известия вузов. Математика, №7, 2020.

3. Al Halidi Arkan M. Effective Computations of Inverse by Module Elements with Approximating K-Ary GDD Algorithm by Ishmukhametov// Al Khalidi Arkan M. // IOP Journal of Physics, 2020, (в печати).

Публикации РИНЦ

4. Al Halidi Arkan M. Эффективные алгоритмы расчета НОД для длинных натуральных чисел/Альхалиди А.М.// Вестник современных исследований, № 8 - 2019,c.1-6.

5. А.М.Аль Халиди, И. Об одном обобщении k-арного алгоритма Евклида/ И.Амер, М.К.Аль Анни, А.М.Аль Халиди // Международная конференция по алгебре, анализу и геометрии, КФУ. — 2016, c. 89-90.

6. Al Halidi Arkan M. Разработка алгоритм k-арного для расчета НОД на высокой скорости / Альхалиди А.М.// МЦНС «Наука и просвещение» В75 -2019, c.11-14.

7. Аль Халиди Аркан Мoхаммед. Ускорение и разработка k-арного алгоритма для расчета НОД / Альхалиди А.М., Аль Тулайхи Мукдад Мухаммид Усман// Вестник современных исследований, № 10 - 2019, c. 1-5.

8. Аль Халиди Аркан Мoхаммед. Сравнение разных алгоритмов вычисления НОД / Альхалиди А.М.// МЦНС Наука и Просвещение. - 2019, c. 12-16.

9. Аль Халиди Аркан Мохаммед. Эффективная реализация аппроксимирующего k-арного алгоритма вычисления НОД в библиотеке MPIR / И.Амер, Альхалиди А.М., Ш.Т.Ишмухаметов // МЦНС Наука и Просвещение. — 2019, с. 10-14.

Литература на русском языке

1. Амер И., Альхалиди А.М., Ишмухаметов Ш.Т. Эффективная реализация аппроксимирующего k-арного алгоритма вычисления НОД в библиотеке MPIR // МЦНС Наука и Просвещение. 2019. — С. 10-14.

2. Амер И. Об ускорении k-арного алгоритма вычисления НОД натуральных чисел// И.Амер, Ш.Т.Ишмухаметов /Ученые записки Казан.унив. (сер.физ.-мат.науки), 2019, Том 161 Книга 1, с. 110-119.

3. Амер Исмаил Ф. О реализации аппроксимирующего k-арного алгоритма вычисления НОД // Вестник современных исследований. — 2018. — Т. 8.1(23). — С. 224-228.

4. Амер Исмаил. Разработка эффективных алгоритмов вычисления НОД натуральных чисел для криптографии и теории чисел, диссертация на соискание ученой степени кандидата технических наук, Казань, 2020, 135 с.

5. Акритас А. Основы компьютерной алгебры с приложениями. М.: Мир, 1994. С. 58-60.

6. Библиотека длинных чисел MPIR. [Электронный ресурс]. Режим доступа: http://www.mpir.org/ (Дата обращения: 12.01.2020).

7. Долгов Д.А. Comparative Analysis of Module Reduction Calculation Algorithms. Research Journal of Applied Sciences, 2015, 10, 442-446.

8. Долгов Д.А. GCD calculation in the search task of pseudoprime and strong pseudoprime numbers, Lobachevskii Journal of Mathematics, 2016, Volume 37, Issue 6.

9. Гельфонд А.О. Решение уравнений в целых числах (Наука, Москва, 1983), с 9-10.

10. Еникеев Р.Р. Анализ алгоритма Скуфа подсчета количества точек на эллиптической кривой // Компьютерные науки и информационные технологии: материалы международной научной конференции. Саратов: Изд-во ИЦ Наука, 2016. С. 153-157.

11. Еникеев Р.Р. Аппроксимация вещественного числа рациональным с ограниченным знаменателем // Ломоносов-2018: сборник тезисов XXV Международной научной конференции студентов, аспирантов и молодых ученых: секция «Вычислительная математика и кибернетика».- Москва: Издательский отдел факультета ВМК МГУ; МАКС Пресс, 2018. - С. 79-80.

12. Еникеев Р.Р. Доказательство корректности и сертификат НОД // Компьютерные науки и информационные технологии: материалы международной научной конференции. Саратов: ИЦ Наука, 2018. - С. 130133.

13. Еникеев Р.Р. Модулярное деление и вычисление обратного элемента по модулю степени двойки // Компьютерные науки и информационные технологии: материалы международной научной конференции. Саратов: ИЦ Наука, 2018. - С. 133-136.

14. Еникеев Р.Р. Нахождение коэффициентов Безу с помощью расширенного обобщенного бинарного алгоритма //Информационные технологии и математическое моделирование (ИТММ-2017). - 2017. - С. 284-291.

15. Еникеев Р.Р. Приближение вещественного числа рациональным в аппроксимирующем к -арном алгоритме // Учен. зап. Казан. ун-та. Сер. Физ.-матем. науки. 2019. - Т. 161, кн. 2. - С. 250-262. - ёо1: 10.26907/25417746.2019.2.250-262

16. Ишмухаметов Ш.Т. Методы факторизации натуральных чисел/ Ш.Т.Ишмухаметов. - Казань, 2011, с. 189.

17. Ишмухаметов Ш.Т. Вычисление коэффициентов Безу для k-арного алгоритма нахождения НОД// Ш.Т.Ишмухаметов, Б.Г.Мубараков, Маад Камаль Аль-Анни/ Изв. ВУЗов. Математика, 2017, №11, с.30-38.

18. Кнут Д.Э. Искусство программирования. Том 2 (3-е издание). — М.: Вильямс, 2001. — С. 367-408.

19. Крэнделл Р., Померанс П. Простые числа: криптографические и вычислительные аспекты// УРРС, М., 2011.

20. Левитин А.В. Алгоритмы. Введение в разработку и анализ. - М.: Вильямс, 2006. - С. 240-242.

21. Фаддеев Д.К. Лекции по алгебре. Учебное пособие для вузов. - М.: Наука, 1984.

22. Хассе Г. Лекции по теории чисел. — М.: Изд. иностранной литературы,1953. — С. 29-30.

23. Хинчин А.Я. Цепные дроби (3-е изд.). - М.: ГИФМЛ, 1961. - 112 с. Литература на английском языке.

24. Akhavi, Ali; Vallée, Brigitte (2000). AverageBit-Complexity of Euclidean Algorithms", Proceedings ICALP'00, Lecture Notes Computer Science 1853: 373387.

25. Amer I. On acceleration of the k-ary GCD Algorithm// I.Amer, S.T.Ishmukhametov/ International Conference on Innovations in Engineering, Technology and Sciences" (ICIETS), IEEE series N 45147, Sept.2018.

26. Amer I. On acceleration of the k-ary GCD Algorithm// I.Amer, S.T.Ishmukhametov./ International Conference on Innovations in Engineering, Technology and Sciences (ICIETS), IEEE series N45147, Sept.2018.

27. Amer I., Ishmukhametov S.T. Analysis of the k-ary Euclid for tuples of integers// Journal of Physics: Conf. Series, 1352 (2019), p.1-4 doi:10.1088/1742-6596/1352/1/012001 https://iopscience.iop.org/article/10.1088/1742-6596/1352/1/012001/pdf.

28. Amer I., Ishmukhametov S., Rubtsova R. Lenstra Factorization Method Convergence Investigation on Elliptic Curves, Research Journal of Applied Sciences 10 (8), (2015), pp. 365-370

29. Apostol Т. Introduction to Analytic Number Theory, 3rd printing. SpringerVerlag, 1986.

30. Arazi O., Qi H. On calculating multiplicative inverses modulo 2m // IEEE Transactions on Computers, 57(10):1435-1438, 2009.

31. Atkin A. and Morain F. Finding suitable curves for the elliptic curve method of factorization. Math. Сотр., 60: 399-405, 1993.

32. А. Atkin and F. Morain. Elliptic curves and primality proving. Math. Сотр., 61:29-68, 1993.

33. Е. Bach. Analytic Methods in the Analysis and Design of Number Theoretic Algorithms. А 1984 АСМ Distinguished Dissertation. The MIT Press, 1985.

34. Bach Е. The complexity of number-theoretic constants. Inform. Process. Letters, 62: 145-152, 1997.

35. Bach E, and Shallit J. Algorithmic Number Theory, volume 1. MIT Press, 1996.

36. Bailey D. and Crandall R. Random generators and normal numbers. Experiment. Math ., 11: 527-546, 2002.

37. Bailey D., Borwein, J., Crandall R., and Pomerance С. On the binary expansions of algebraic numbers. Journal de Theorie des Nombres, Bordeaux, 16: 487-518, 2004.

38. Barrett Р. Implementing the Rivest Shamir and Adleman public key encryption algorithm on а standard digital signal processor. In А. Odlyzko, editor , Advances in Cryptology, Proc. Crypto '86, volume 263 of Lecture Notes in Computer Science, pages 3 1 1-323. Springer-Verlag, 198.7

39. Bernstein D. Doubly focused enumeration of locally square polynomial values. In High primes and misdemeanours: lectures in honour of the 60th birthday

of Hugh Cowie WilliaTS, volume 41 of Fields Inst . Commun. , pages 69-76. Amer. Math . Soc, 2004.

40. Bernstein D. Fast multiplication and its applications. In J . Buhler and P. Steinhagenn, editors Cornerstones in Algorithmic Number theory (tentative title), a Mathematical Sciences Research Institute Publication, Cambridge University Press.

41. Bleichenbacher D. Efficiency and security of cryptosystems based on number theory. PhD thesis, Swiss Federal Institute of Technology Zurich, 1996.

42. Boneh D. Twenty years of attacks on the RSA cryptosystem. Notices Amer. Math. Soc., 46: 203-213, 1999.

43. Borwein J, Bradley D., and Crandall R. Computational strategies for the Riemann zeta function. J. Comp. App. Math., 121: 247-296, 2000.

44. Bojanczyk A. W., Brent R.P. A systolic algorithm for extended gcd computation // Comput. Math. Applic. vol. 14, no. 4, P. 233-238, 1987.

45. Bredihin B. Applications of the dispersion method in binary additive problems. Dokl. Akad. Nauk. SSSR, 149: 9-11, 1963.

46. Brent R. P. Analysis of the binary Euclidean algorithm, in New Directions and Recent Results in Algorithms and Complexity (edited by J. F. Traub), Academic Press, New York, 1976, 321-355.

47. Brent, Richard P. Twenty years' analysis of the Binary Euclidean Algorithm, Millenial Perspectives in Computer Science: Proceedings of the 1999 OxfordMicrosoft Symposium in Honour of Professor Sir Antony Hoare, Palgrave, NY: (2000), pp. 41-53.

48. Brent R. On the period of generalized Fibonacci recurrences. Math. Comp. 63: 389-401, 1994.

49. Brent R. Factorization of the tenth Fermat number. Math. Comp. 68: 429451, 1999.

50. Bressoud D. and Wagon S. A Course in Computational Number Theory. Key College Publishing, 2000.

51. Buchmann J, Jacobson, М. Jr and Teske Е. On some computational problems in finite groups. Math. Comp. 66: 1663-1687, 1997.

52. Buhler J., editor. Algorithmic Number Theory: Proc. ANTSIII, Portland, OR, volume 1423 of Lecture Notes in Computer Science. Springer-Verlag, 1998.

53. Buchberger B., Jebelean T. Parallel Rational Arithmetic for Computer Algebra systems: Motivating Experiments. 1992.

54. Computing multiplicative inverses [Электронный ресурс] - Режим доступа: goo.gl/Js5kYp. [Дата обращения: 11.04.2019].

55. D. Coppersmith. Small solutions to polynomial equations and low exponent RSA vulnerabilities, J. Cryptology, 10: 233-260, 1997.

56. Couveignes J., Dewaghe L, and Morain F. Isogeny cycles and the Schoof-Atkin-Elkies algorithm. J. Cryptology, 10: 233-260, 1997.

57. Davenport H. Multiplicative Number Theory (second edition). SpringerVerlag, 1980.

58. Davis M. Hilbert 's tenth problem is unsolvable. Amer. Math. Monthly, 80:233-269, 1973.

59. Diffie W. and Hellman M. New directions in cryptography. IEEE Trans. Inform. Theory, 22: 644-654, 1976.

60. Dixon J. The number of steps in the Euclidean algorithm, Journal of Number Theory (2), 414-422 (1970).

61. Dolgov D.A. An extended Jebelean-Weber-Sedjelmaci GCD algorithm // Chebyshevskii Sb., Kazan, KFU, 2018, Volume 19, Issue 2, pp. 421-431.

62. Dolgov D., GCD calculation in the search task of pseudoprime and strong pseudoprime numbers, Lobachevskii Journal of Mathematics 37 (6), 734{739 (2016).

63. Dubner D. and Gallot У. Distribution of generalized Fermat numbers. Math. Comp. 71: 825-832, 2002.

64. Dumas J. -G. — On Newton-Raphson iteration for multiplicative inverses modulo prime powers. arxiv (site).

65. Dusse S. R, Kaliski B. S. — A Cryptographic library for the Motorola DSP56000 // EUROCRYPT 90, pages 230-244, 1990.

66. Elliptic curve primality. [Электронный ресурс] - Режим доступа: https://en.wikipedia.org/wiki/Elliptic_curve_primality. - [Дата обращения: 05.03.2019].

67. Farey sequence. [Электронный ресурс] - Режим доступа: https://en.wikipedia.org/wiki/Farey_sequence, свободный. - [Дата обращения: 01.05.2019].

68. Friedlander J, Pomerance С, and Shparlinski L. Period of the power generator and small values of Carmichael's function. Math. Comp. 70: 1591-1605, 2001.

69. Galiev A, Prokopyev N, Ishmukhametov S, Latypov R.H., Stolov E.L Archain: A Novel Blockchain Based Archival System// Proceedings of the 2nd World Conference on Smart Trends in Systems, Security and Sustainability, WorldS4 2018. P.308-312.Gardner M. Mathematical games: а new kind of cipher that would take millions of years to break. Scientific American, August 1977.

70. Gladman B., Hart W., Moxham J., et al. MPIR: Multiple Precision Integers and Rationals, 2015. version 2.7.0, A fork of the GNU MP package (T. Granlund et al.) URL: http://mpir.org

71. Golomb S. Combinatorial proof of Fermat's 'little theorem'. Amer. Math. Monthly, 63, 1956.

72. Gordon D. and Pomerance C. The distribution of Lucas and elliptic pseudoprimes. Math. ^тр. 57: 825-838, 1991.

73. Grantham J. А probable prime test with high confidence. J. Number Theory, 72: 32-47, 1998.

74. Guy R. Unsolved Problems in Number Theory. Second edition, volume 1 of Problem Books in Mathematics. Springer-Verlag, 1994.

75. Hardy G. and Wright E. An Introduction to the Theory of Numbers. Fifth edition. Clarendon Press, Oxford, 1979.

76. J Hardy G.H., Wright. E.M. An Introduction to the Theory of Numbers (Fourth Edition). - Oxford: Oxford University Press, 1975. - 433 p.

77. Ishmukhametov S.T., Mubarakov B. On practical aspects of the Miller-Rabin Primality Test, Lobachevskii Journal of Mathematics, October 2013, Volume 34, Issue 4, pp 304-312

78. Ishmukhametov S.T. An Approximating k-ary GCD Algorithm// S.T.Ishmukhametov. / Lobachevskii Journal of Mathematics, 2016, Volume 37, Issue 6, pp. 723-729

79. Ishmukhametov S.T. Mubarakov B.G., and Al-Anni K. M. Calculation of Bezout's Coefficients for the k-Ary Algorithm of Finding GCD // Russian Mathematics, 2017, Vol. 61, No. 11, pp. 26-33. c Allerton Press, Inc., 2017.

80. Ishmukhametov S.T., Rubtsova R.G. A parallel computation of the GCD of natural numbers// Параллельные вычислительные технологии - XI международная конференция, ПаВТ'2017, г. Казань, 3-7 апреля 2017 г. Челябинск: Издательский центр ЮУрГУ. - 2017, p.120-129.

81. Ishmukhametov S.T., Rubtsova R.G. A fast algorithm for counting GCD of natural numbers, тезисы докладов международной конференции «Алгебра, анализ и геометрия«, Казань, КФУ, 2016, с.52-53.

82. Ishmukhametov S.T., Rubtsova R.G. Saveleyv N. The Error Probability of the Miller-Rabin Primality Test, ISSN 1995-0802, Lobachevskii Journal of Mathematics, 2018, Vol. 39, no.7, pp. 1010-1015., Pleiades Publishing, Ltd., 2018.

83. Janson S., Resultant and discriminant of polynomials (Lecture Notes, 2007).

84. Jebelean T. A Generalization of the Binary GCD Algorithm// Proc. of Intern.Symp.on Symb.and Algebr.Comp.(ISSAC'93), pp.111-116 (1993).

85. Jebelean T.^ Double-Digit Lehmer-Euclid Algorithm for Finding the gcd of Long Integers, J.Symb. Computations, (1995), 19, pp. 145-157.

86. Jebelean T.An Algorithm for Exact Division, J. Symbolic Computation, Vol. 15 (1993).

87. Jones G.A., Jones J.M. Elementary Number Theory. — Berlin: SpringerVerlag, 1998. — p. 7-11.

88. Johnson D, Menezes А., and Vanstone S. The elliptic curve digital signature algorithm (ECDSA). International Journal of Information Security, 1: 36-63, 2001.

89. Ko? Q.K. A New Algorithm for Inversion mod pk // IACR Cryptology ePrint Archive, 2017. — Report 2017/411. — 19 p.

90. Knuth, Donald (1998), Seminumerical Algorithms, The Art of Computer Programming, 2 (3rd ed.), Addison-Wesley, ISBN 978-0-201-89684-8.

91. Kurosh A.G., Higher Algebra (MIR Publishers, 1980).

92. Lehmer D. Euclid's algorithm for large numbers, Am. Math. Mon., 45: pp.227-233, 1938.

93. Li B, Liu L, Zhi L., A structured rank-revealing method for Sylvester matrix, Journal of Computational and Applied Mathematics 213 (1), 212-223 (2008).

94. List of C++ multiple precision arithmetic libraries. [Электронный ресурс] -Режим доступа:

https://en.wikipedia.org/wiki/List_of_C++_multiple_precision_arithmetic_libraries - (Дата обращения: 24.05.2019).

95. Matula D., Fit-Florea A. and Thornton M. — Table lookup structures for multiplicative inverses modulo 2/sup k/ — Proc.17th IEEE Symp. Computer Arithmetic, (2005), pp.156-163.

96. Menezes, Alfred; van Oorschot, Paul C.; Vanstone, Scott A. (October 1996). Handbook of Applied Cryptography. CRC Press. ISBN 0-8493-8523-7.

97. Moenck R. Fast computation of GCDs. In Proc. 5th Annual АСМ Symposium оп the Theory of Computing, pages 142-151, 1973.

98. Moller. On Schonhage's Algorithm and Subquadratic Integer gcd Computations, 2007.

99. Montgomery P. Modular multiplication without trial division. Math. Comp., 44: 519-521, 1985.

100. Montgomery, P.L. — Modular Multiplication Without Trial Division — Mathematics of Computataions vol 44 number 170 april 1985, p. 519-521.

101. Peterson L. Great computations. Science News, 157(10): 152-153, March 4, 2000.

102. Purdy G.B. A carry-free algorithm for finding the greatest common divisor of two integers. // Comput. Math. Applic. 9, p. 311-316, 1986.

103. Shoup V.A Computational Introduction to Number Theory and Algebra (Version 2). - Cambridge: Cambridge University Press, 2008. - p. 21.

104. Sorenson J. The k-ary GCD Algorithm//J.Sorenson/ Universitet of Wisconsin-Madison, Tecn.Report (1990), pp.1-20.

105. Sorrenson J. Two fast GCD Algorithms// J.Sorenson/ J.Alg. 16, №1, 110144 (1994).

106. Sorenson, J. An analysis of the generalized binary GCD algorithm, in: A. van der Poorten, A. Stein (Eds.), High Primes and Misdemeanors, Lectures in Honour of Hugh Cowie Williams, Banff, Alberta, Canada, AMS Math. Review 41, 254-258 (2004).

107. Sorenson, J. Webster, J. Strong pseudoprimes to twelves prime bases, Math. Comput. 2016, pp.1-17 DOI: 10.1090/mcom/3134

108. Schonhage A. Schnelle Berechnung von Kettenbruchentwicklungen, Acta Informatica, 1971, N1, 139-144.

109. Sedjelmaci Sidi Mohamed. Jebelean-Weber's algorithm without spurious factors // 2007, Information Processing Letters 102(6):247-252 Stein, J. Computational problems associated with Racah algebra", Journal of Computational Physics, 1 (3): 397-405, 1967.

110. Stein, J. (1967), "Computational problems associated with Racah algebra", Journal of Computational Physics, 1 (3): 397-405.

111. Stokes Jon (2007). Inside the Machine: An Illustrated Introduction to Microprocessors and Computer Architecture. No Starch Press. p. 117. Robert

Reese; J. W. Bruce; Bryan A. Jones (2009). Microcontrollers: From Assembly Language to C Using the PIC24 Family. Cengage Learning. p. 217.

112. Stehle, Damien; Zimmermann, Paul. A binary recursive gcd algorithm", Algorithmic number theory (PDF), Lecture Notes in Comput. Sci., 3076, Springer, Berlin, pp. 411-425 (2004).

113. Van der Waerden B.L, Emfuhrung in die AIgebrmsche Geometrze, 2nd ed. (Springer Verlag, New York, 1973).

114. Vallee B. Dynamical analysis of a class of Euclidean algorithms, Theor. Comp. Science, 297 (203), pp. 447-486.

115. Vardi I. Computational Recreations in Mathematica. - Boston: Addison-Wesley, 1991. - p. 155.

116. Vinberg E.V., A Course in Algebra (American Mathematical Society, 2003).

117. Weber K. The accelerated integer GCD algorithm// ACM Trans.Math.Software, 21, №1, 1-12 (1995).

118. Weber K. The Accelerated Integer GCD Algorithm, ACM Trans.Math. Software, 21, №1, р.111-122 (или 1-12) (1995).

119. Weber K., The accelerated integer GCD algorithm, ACM Transactions on Mathematical Software 21 (1), 111-122 (1995).

120. Multiple Precision Integers and Rationals Library [Электронный ресурс]. Режим доступа: www.mpir.org (Дата обращения: 07.01.2020).

Список сокращений и условных обозначений

АКА—аппроксимирующий -арный алгоритм. КАЕ—классический алгоритм Евклида. НОД—наибольший общий делитель.

DSS (Digital Signature Standard)—американский стандарт цифровой подписи.

GCD (Greatest Common Divisor)—наибольший общий делитель.

MPIR (Multiple Precision Integers and Rationals)—библиотеки длинных чисел.

RSA—алгоритм шифрования, который предложили Райвест (Rivest), Шамир

100

(Shamir), и Адлеман (Adleman).

GMP (Good Manufacturing Practice) — надлежащая производственная практика.

Список рисунков

1.1 Бинарный алгоритм....................................................................16

2.1. Процедура для массива a[0]...a[6] и опорного элементаp=a[3]............48

2.2. Диаграмма парного алгоритма.....................................................56

2.3. Диаграмма аппроксимирующим парным алгоритмом........................57

2.4. Алгоритм Соренсона при £=64......................................................58

2.5. парный алгоритм Д.Соренсона....................................................59

2.6. Аппроксимирующий парный алгоритм..........................................59

2.7. Диаграмма времен вычисления при различных к..............................60

2.8. Диаграмма времен вычисления при различных к...............................60

3.1. Диаграммы зависимостей времени вычисления для к = 1 6.................64

3.2. Диаграмма выигрыша времени вычисления для к = 16......................65

3.3. Сводная диаграмма числа итераций для к = 16...............................65

3.4. Диаграмма выигрыша числа итераций для к = 16..............................66

3.5. Сводная диаграмма времени вычисления для Ь = 5 0 0.......................67

3.6. Диаграмма выигрыша времени вычисления для Ь = 5 0 0....................68

3.7. Сводная диаграмма числа итераций для Ь = 5 00..............................68

3.8. Диаграмма выигрыша времени вычисления для Ь = 5 0 0....................69

4.1. Данные о числе вариантов при L=50..............................................85

4.2. Данные о числе вариантов при L=100............................................85

4.3. Данные о числе вариантов при L=200............................................86

4.4. Данные о числе вариантов при L=500............................................87

Список таблиц

1.1. Нахождение коэффициентов x, y.................................................16

1.2. Пример вычисления НОД с использованием алгоритма АКА..............28

2.1. (Kary) Оценка числа итерации в зависимости от k............................50

2.2. (Kary) Оценка скорости вычисления в зависимости от k...........................51

2.3. (Kary-Inverses with Mas) Число итераций в зависимости от k...............51

2.4. (Kary-Inverses with Mas)Время вычисления в зависимости от к............51

2.5. (KaryX,Ywith Mas) Число итераций в зависимости от k......................51

2.6. (Kary X, Y with Mas ) Время в зависимости от к.................................52

2.7. (Алгоритм АКА) Число итераций в зависимости от k.........................52

2.8. (Алгоритм АКА) Время вычисления в зависимости от k.....................52

2.9. (Алгоритм АКА + сортировка)Число итераций в зависимости от k......53

2.10. (Алгоритм АКА + сортировка) Скорость вычислений в зависимости от k................................................................................................53

2.11. (NEW-AKA) Оценка итерации вычисления НОД для 200 пар. Длина (1000-3000) десятичных разрядов в зависимости от k..............................53

2.12. (NEW-AKA) Оценка скорости вычисления НОД для 200 пар. Длина (1000-3000) десятичных разрядов в зависимости от к.............................54

2.13. Классический алгоритм Евклида.................................................55

2.14. k-арный алгоритм Соренсона.....................................................55

2.15. Аппроксимирующий k-арный алгоритм........................................55

3.1. Сводная таблица времени вычисления для к = 1 6............................64

3.2. Сводная таблица числа итераций для к = 16...................................65

3.3. Сводная таблица времени вычисления для L = 5 00...........................67

3.4. Сводная таблица числа итераций для L = 5 00.................................68

4.1. Пример решения уравнения Безу..................................................74

4.2. Прямой ход вычисления НОД по схеме аппроксимирующего алгоритма....................................................................................78

4.3. Обратный ход k-арного алгоритма..............................................79

103

4.4. Прямой ход аппроксимирующего алгоритма.................................80

4.5. Сводная таблица вариантов вычисления обратного элемента..............82

4.6. Данные о числе вариантов при L=50..............................................84

4.7. Данные о числе вариантов при L=100............................................85

4.8. Данные о числе вариантов при L=200............................................86

4.9. Данные о числе вариантов при L=500............................................86

4.10. Общее время вычисления для 100 пар чисел длины L=50..................88

4.11. Общее время вычисления для 100 пар чисел длины L=100................88

4.12. Общее время вычисления для 100 пар чисел длины L=200................88

4.13. Общее время вычисления для 100 пар чисел длины L=500................89

4.14. Общее время вычисления для 100 пар чисел длины L=1000...............89

ПРИЛОЖЕНИЯ

Процедуры на языке C

Приложение 1. Процедура генерации п пар взаимно-простых чисел размерности L.

void Generate(int n, int L) {

int i = 0, l = L / 5, r, r1, r2; uint64_t d; for (i = 0; i < n; i + +) {

r = 50000 + rand() % 50000;

mpz_set_ui(A, r); r = 10000 + rand() % 4 0000; mpz_set_ui(B, r); for (int j = 2; j < l; j++) { r = rand() % 100000; mpz_mul_ui(A, A, 100000); mpz_add_ui(A, A, r); r = rand() % 100000; mpz_mul_ui(B, B, 100000); mpz_add_ui(B, B, r);

}

r = rand() % 100000; if (r % 2 == 0)r--;

mpz_mul_ui(A, A, 100000); mpz_add_ui(A, A, r); r = rand() % 100000; if (r % 2 == 0)r--;

mpz_mul_ui(B, B, 100000); mpz_add_ui(B, B, r); mpz_set(A1, A); mpz_set(B1, B); d = EuclidM(A1, B1, C); d = mpz_cmp_ui(C,1); if (d > 0) {mpz_fdiv_q(A, A, C); mpz_fdiv_q(B, B, C);} mpz_init_set(Am[i], A); mpz_init_set(Bm[i], B);

}

}

Приложение 2. Процедура вычисления наибольшего общего делителя по схеме Евклида для библиотеки MPIR

uint64_t EuclidM(mpz_t A , mpz_t B, mpz_t & ) { mpz_t A1, B1;

mpz_init(A1);mpz_init(B1); mpz_set(A1, A); mpz_set(B1, B); long i = 0; size_t l = mpz_size(B1); while (l>0) {i++;

mpz_fdiv_r( , A1, B1); l = mpz_size(C); if (l ==_0) { mpz_set (C, B1); return i; } mpz_set(A1, B1); mpz_set(B1, C);

}

mpz_set( ,B1); return i;

mpz_clear(A1); mpz_clear(B1);

Приложение 3. Процедура решения уравнения Безу (расширенный алгоритм Евклида) для 64-битовых чисел типа unsigned long long.

#define ull unsigned long long

// Ext_GCD solves the Bézout's identity and returns positive x,y and z=sign(y).

ull Ext_GCDT(ull a, ull b, ull & x, ull & y, int & zn) { ull d, x1, y1,t; int z; if (a%b == 0) {

d = a; x = 0; y = 1; z =1; return b;

}

d = Ext_GCDT(b, a%b, x1, y1, z); x = y1; t = (a / b) * y1; y = t - x1; z = -z; return d;

}

Приложение 4. Процедура решения уравнения Безу (расширенный алгоритм Евклида) для длинных чисел типа mpz_t.

int ExtEuclidT(mpz_t a, mpz_t t, mpz_t & ], mpz_t & y) {//у=ЬЛ-1 mod a, answer=1 iff a<0 & b>0

mpz_t t,a1,b1; mpz_init(t); mpz_init(a1); mpz_init(b1);

mpz_set(a1, a); mpz_set(b1, b);

mpz_t Div[5000];

int i = -1, s, z;

while (true) {

i++; mpz_fdiv_qr( , x, a1, b1); mpz_init_set(Div[i], y); s = mpz_cmp_ui( , 0);

if (s == 0) break; mpz_set(a1, b1); mpz_set(b1,

x);

}

mpz_set_ui(x, 0); mpz_set_ui(y, 1); if (i % 2 == 0)z = 1; else z = -1;

while (i > 0) { i--; mpz_set(t, y); mpz_mul( , y,

Div[i]);

mpz_clear(Div[i]); mpz_add( , y, x); mpz_set(x, t);

}

mpz_clear(t); mpz_clear(a1); mpz_clear(b1); return z;

Приложение 5. Процедура быстрого возведения в степень по модулю для чисел формата unsigned long long.

uint64_t PowMod( uint64_t a, uint64_t b, uint64_t n){ uint64_t k=1,c; while (b) { if (b%2==0) {

b /= 2; c=a*a; a = c%n; }

else { b--;

k = (k*a)%n;

}

}

return k;

}

Приложение 6. Процедура быстрого возведения в степень по модулю для чисел формата mpz_t.

void PowModT(mpz_t a, long b, mpz_t n, mpz_t & ) { mpz_t d; mpz_init(d); mpz_set_ui(c, 1); while (b) {

if (b % 2 == 0) { b /= 2;

mpz_mul(d,a,a); mpz_fdiv_r( ,d,n);

}

else {

b--;

mpz_mul( , c, a);

mpz_fdiv_r( , c, n); }

}

return; mpz_clear(d);

}

Приложение 7. Процедура Farey построения дроби Фарея для заданного числа формата float из интервала (0; 1/2).

void Farey(long k, double alp, long& m, long& n) { int m1, n1, m2, n2; double x; x = 1 / alp;

m1 = 1; m2 = 1; n2 = int(x); n1 = n2 + 1; while (n1 + n2 < k) { m = m1 + m2; n = n1 + n2; x = double(m) / n; if (x < alp) { m1 = m;

n1 = n ; }

else {

m2 = m; n2 = n; }

}

x = double(m1) / n1 + double(m2) / n2 - 2 * alp; if (x>0) {

m = m1; n = n1; }

else {

m = m2; n = n2; }

return;

}

Приложение 8. Процедура быстрого вычисления битовой длины числа формата unsigned long long.

int blen( uint64_t f) {

int l=1, h2,h = (f < 4294967296);

uint64_t b; if (h == 1) {h2 = (f < 65536); if (h2 == 1) {

l = 16; b = 32768; while (f < b) { b->1; l--;

}

}

else {l = 32; b = 65536; while (f < b) { b->1; l--;

}

}

}

else {

h2 = (f < 4295032832); if (h2 == 1) {

l = 48; b = 2147516416;

while (f < b) {

b->1; l--;

}

else {l = 64;

b = 9223372036854775808;// 2Л63 while (f < b) { b->1; l--;}}

}

return l;

}

Приложение 9. Процедура вычисления НОД по схеме аппроксимирующего алгоритма с использованием двух стратегий

int AKA(mpz_t A, mpz_t B, mpz_t & C, int s) { mpz_t A1,B1,T1,T2,Tmp;

mpz_init(A1); mpz_init(B1); mpz_init(T1);

mpz_init(T2); mpz_init(Tmp);

long k = s*s, b1, q, m, x, t;

int s0, salp, i = 0, a, b, l, d;

uint64_t y,f1, f2, f3, f6 = 65356;

float alp, beta, r;

float fmin = 1. / double(k);

size_t l1 = mpz_size( ); size_t l2 = mpz_size( ); mpz_set(A1, A); mpz_set(B1, E);

//************* Main Cycle

********************************

{

while (l1 > 1 i++;

f1 = A->_mp_d[l1 - 1]; f2 = (B->_mp_d[l2 - 1] //Sub 1. Computation r=A/B if (l1 == l2) { if (l < 48) {

f1 = (f1 << (64 - l) f2 = (f2 << (64 - l) r = float(f1) / f2;

l = blen(f1'

+(( ->_mp_d[l1 - 2]) >> l); + (( ->mpd[l1 - 2]) >> l)

}

}

64

l)) + (( -> mp d[l1

else

// case len(A)>len(B) {

f1 = (f1 << ( f2 = f2 >> l;

}

r = float(f1) / f2; // Sub 2. Computation q = ABA{-1}mod k a = ( ->_mp_d[0]) % k; b = ( ->_mp_d[0]) % k; q = (a*InvK[b]) % k; // Sub 3. Computation alpha=int((r-q)/k) beta = (r - q) / k;

2]

>> l)

if (beta > 0) {

s0 = int(beta); alp = beta - s0; salp = 1;// sign(alpha)positive if (alp > 0.5) { s0 + +;

alp = 1 - alp;

salp = -1; // sign(alpha)negative }

}

else {

s0 = int(beta) - 1; alp = beta - s0; salp = 1; if (alp > 0.5) { s0 + +;

alp = 1 - alp; salp = -1;

}

}

// fmin states a level below which we apply strategy 2, otherwise,

// strategy 1. By strategy 2 we replace standard Farey procedure by

// direct computation x=1, y=-q-k*s0;

if (alp < fmin) { // strategy 2 x = 1;

f3 = q + k*s0; //f3=-y

}

else {

Farey(k, alp, m, x); // strategy 1 m = m*salp;

y = q*x + (s0*x + m)*k;

}

mpz_mul_ui(T1, A, x); // Ax mpz_mul_ui(T2, B, y);// By b1 = mpz_cmp(T1, T2); if (b1 == 0) {

//----------- Case Ax+By=0 --------------

mpz_fdiv_q_ui(T1, x); m = (T1->_mp_d[0]) % 2; while (m == 0) {

mpz_div_ui(T1, T1, 2); m = (T1->_mp_d[0]) % 2;

}

t= EuclidM(A1, T1, C); mpz_set( , C); t = EuclidM(B1, B, C); return i + t;

}

if (b1 > 0)

mpz_sub( , T1, T2);

else

mpz_sub( , T2, T1); mpz_fdiv_q_ui(Tmp, C, k); f1 = Tmp->_mp_d[0]; d = 1; while (f1 % 2 == 0) {

d *= 2; f1 /= 2; }

if (d>1)

mpz_div_ui(Tmp, Tmp, d); d = mpz_cmp(B, Tmp); if (d == 0) {

mpz_fdiv_r( , A1, B); m = A->_mp_d[0]; while (m % 2 == 0) {

mpz_div_ui(A, A, 2); m = ( ->_mp_d[0]) % 2;

}

mpz_fdiv_r( , B1, A); m = ( ->_mp_d[0]) % 2; while (m == 0) {

mpz_div_ui(B, B, 2); m = ( ->_mp_d[0]) % 2;

}

t = EuclidM(A, B, C); return i + t;

}

if (d > 0) {

mpz_set( , B); mpz_set( , Tmp); l1 = l2;

l2 = mpz_size( );

}

else {

mpz_set( , Tmp); l1 = mpz_size( );

}

fflush(Ans);

}

f1 = A->_mp_d[0]; f2 = B->_mp_d[0]; if (f2 == 0) { //Case B=0. Remove extra factors t = EuclidM(A1, A, C); i += t;

mpz_set(B, C); t = EuclidM(B1, B, C); return i + t;

}

if (f1 < f2) { f3 = f1; f1 = f2; f2 = f3;

while (f2 > 0) { i++;

rr = float(f1) / f2; if (f1%f2 == 0) {

mpz_set_ui( , f2); t = EuclidM(A1, A, C); i += t;

mpz_set(B, ( ); t = EuclidM(B1, B, C); return i + t;

}

a = f1 % k; b = f2 % k; q = (a*InvK[b]) % k; beta = (rr - q) / k; if (beta > 0) {

s0 = int(beta); alp = beta - s0; salp = 1; if (alp > 0.5) {

s0++; alp = 1 - alp; salp = -1;

}

}

else {