Код вопроса вычисляет 64-битный получить<2>(I)=$h:=f\,g\bmod n$ из входов:
- 64-битный
модуль_значение=$n$ с $n\in[2,\,2^{63}]$
- 64-битный
получить<0>(I)=$f$ с $f\in[0,\,2^{64}-1]$
- 64-битный
получить<1>(I)=$г$ с $г\в[0,\,2^{64}-1]$ и $f\,g<2^{64}\,n$, условие, которое выполняется, если $f,g\in[0,\,n-1]$ (что, я думаю, всегда имеет место в приложении).
Результат $ч$ остаток от евклидова деления $ф\,г$ к $n$. Он математически определяется $0\le h<m$ и $\существует q\in\mathbb Z,\ f\,g=q\cdot n+h$.
Код сначала вычисляет 128-битное г$=z:=f\,g$, то результат получить<2>(I)$=h:=z\bmod n$ к редукция Барретта:
- Он был предварительно вычислен (вне кода вопроса)
const_ratio$=r=\left\lfloor2^{128}/n\right\rfloor$
- $\hat q:=\left\lfloor z\,r/2^{128}\right\rfloor$, что правильно $q$ в пределах одного по умолчанию (примечание: $\шляпа q$ это конечное значение переменной
tmp1).
- $\шляпа ч:=z-q\cdot n$, что правильно $ч$ в пределах возможного превышения точно $n$ (Примечание: $\шляпа ч$ это конечное значение переменной
tmp3).
- $h:=\шляпа h-n$ когда $\шляпа ч\ге п$, или же $ч:=\шляпа ч$ в противном случае.
Код использует алгоритмы начальной школы для выполнения многозначной арифметики, перенесенной из базы $10$ основать $2^{64}$ (с оптимизацией и вариантом, подробно описанным в следующем разделе). Эквивалентом цифр являются так называемые конечности, здесь 64-битная.
Продукт г$=г$ выражается двумя конечностями г [0]$=z_0$ (низкий порядок), г [1]$=z_1$ (высокий порядок), таким образом, с $z=z_0+2^{64}\,z_1$, и $z_0, z_1\in[0,\,2^{64}-1]$.
Множитель Барретта const_ratio$=r=\left\lfloor2^{128}/n\right\rfloor$ аналогично выражается как две конечности const_ratio_0$=r_0$ и const_ratio_1$=r_1$, благодаря нижнему концу интервала в предусловии $n\in[2,\,2^{63}]$.
Промежуточный продукт $г\,г$ всегда соответствует трем ветвям (хотя в целом произведение двух количеств, выраженных в виде двух ветвей, потребовало бы четырех ветвей).
Предварительный коэффициент по умолчанию $\шляпа q$ подходит для одной конечности, так как $\шляпа д\ле д$ и $q<2^{64}$, с последним застрахованным входным условием $f\,g<2^{64}\,n$.
Предварительный остаток $\шляпа ч$ подходит для одной конечности, так как $\шляпа h<2n$ и $2n\le2^{64}$, с более поздним благодаря верхнему концу интервала в предусловии $n\in[2,\,2^{63}]$.
Детализация алгоритма кода (как задано в комментарий и награда):
Я буду использовать иллюстрацию в десятичном формате. В этой базе, поскольку $2\le n\le10/2$, const_ratio$=г$ может быть только $\влево\lfloor100/2\вправо\rfloor=50$, $\лево\lfloor100/3\право\rfloor=33$, $\лево\lfloor100/4\право\rfloor=25$, $\лево\lfloor100/5\право\rfloor=20$, но я притворюсь const_ratio$=г=29$ потому что это делает более интересный пример. По той же причине я буду использовать г$=г=34$, даже несмотря на то, что это нельзя получить как произведение двух цифр.
Продукт г получается в коде умножить_uint64 (получить<0>(I), получить<1>(I), z) как две конечности г[0] и г[1].
Мясо вычислений $\hat q:=\left\lfloor z\,r/2^{128}\right\rfloor$. Это аналог в базе $2^{64}$ из $9:=\влево\lfloor29\cdot34/100\вправо\rfloor$ в базе 10. Оба аргумента $29$ и $34$ к умножению двузначны, но достаточно малы, чтобы их произведение $986$ является трехзначным (а не четырехзначным), и нас интересует только третья цифра справа. Алгоритм начальной школы для вычисления $986:=29\cdot34$ будет представлен как
2 9 const_ratio
х 3 4 з
-----
1 1 6
+ 8 7
-------
9 8 6
В алгоритме начальной школы есть четыре однозначных умножения (которые выполняет код) и несколько дополнительных операций (которые код немного реорганизует):
4 раз 9, 36; записывать 6, держать 3;4 раз 2, 8; плюс 3 (сохранено), 11; напиши это.3 раз 9, 27; записывать 7, держать 2;3 раз 2, 6; плюс 2 (сохранено), 8; напиши это.
Первое из этих четырех умножений происходит во фрагменте кода умножить_uint64_hw64 (z [0], const_ratio_0, и перенос), который умножает младший член $г$ с младшим коленом $z$, как мы умножаем младшую цифру 4 из 34 с младшей цифрой 9 из 29. Обратите внимание, что «напишите 6" в данных обстоятельствах бессмысленно, так как любая цифра, которую он записывает, останется отделенной в правом столбце вычислений без какой-либо возможности повлиять на крайнюю левую цифру и будет игнорироваться, когда мы делим на 100 и округляем в меньшую сторону (эквивалентно, оставляем только третью цифру из справа). Вот почему младший 64-разрядный 128-разрядный продукт даже не вычисляется, как указано в вопросе. Эквивалент 3 в 36 хранится в нести.
Второе умножение происходит в умножить_uint64(z[0], const_ratio_1, tmp2), который умножает старший член $г$ с младшим коленом $z$, в результате чего две конечности tmp2; 64-битный температура [0] получает эквивалент 8 в 8, и тмп[1] получает эквивалент 0 за (обратите внимание на ведущий 0 подавляется в обычном написании десятичных целых чисел). Эквивалент 8 плюс 3 происходит в add_uint64(tmp2[0], перенос, &tmp1), с младшей цифрой 1 результата 11 в tmp1, и новый перенос 1 в выводе этой функции. Это используется как правый операнд в tmp3 = tmp2[1] + … (который пропускается в алгоритме начальной школы с конкретным примером, который я взял с тех пор, как 0 был подавлен), давая эквивалент левой 1 в 116. [Примечание о выводе add_uint64: это сгенерировано static_cast<беззнаковый символ>(*результат <операнд1), что сравнивает *результат и операнд, затем поворачивает истинный к 1, ЛОЖЬ к 0. Сделано после *результат = операнд1 + операнд2, который сообщает, вызвало ли это добавление перенос. Некоторые компиляторы распознают эту идиому, используют бит C слова состояния и повторно используют C в следующем добавлении].
Третье умножение происходит в умножить_uint64 (z [1], const_ratio_0, tmp2), который умножает младший член $г$ с ветвью высокого порядка $z$, с результатом на конечности в tmp2, как и мы 3 х 9 = 27. На этот раз нам нужны обе конечности/цифры: эквивалент 7 идет к тмп2[0] и эквивалент 2 идет к тмп2[1]. Вот сделал вариант алгоритма начальной школы: сразу добавил tmp1 (аналог среднего 1 в 116) на младшую конечность с add_uint64(tmp1, tmp2[0], &tmp1), выполняя эквивалент 1 + 7 = 8, без переноски. Результат 8 хранится в tmp1 потому что семантика add_uint64 нуждается в пункте назначения, но на самом деле он игнорируется, потому что нас не волнует средняя цифра в 986. Выход переноса на add_uint64 используется как правый операнд в перенос = tmp2[1] + …, выполняя эквивалент 1 + 0 = 1 в нашем примере. Несмотря на название нести, который содержит полноценную 64-битную конечность/цифру.
Четвертое умножение происходит в г[1] * const_ratio_1, который умножает старший член $г$ с ветвью высокого порядка $z$, как и мы 3 х 2 = 6. Здесь контекст гарантирует, что результат соответствует одной конечности, поэтому можно использовать собственный оператор C для умножения. Затем результат используется как левый оператор ¦ + tmp3 + перенос, выполняя эквивалент 6 + 1 + 1 = 8. Опять же контекст гарантирует эти значения $\шляпа q$, Хранится в tmp1, подходит для одной конечности/цифры.
затем tmp3 = z[0] - tmp1 * значение_модуля выполняет $\шляпа ч:=z-q\cdot n$. Контекст гарантирует, что математически точный результат соответствует одной конечности/цифре (хранящейся в tmp3) даже не смотря на $q\cdot n$ не. Это позволяет использовать родные операторы C, которые полностью пропускают вычисление ветви высокого порядка.
затем SEAL_COND_SELECT(tmp3 >= значение_модуля, tmp3 - значение_модуля, tmp3) вычисляет $ч$ от $\шляпа ч$ условно вычитая $n$ когда $\шляпа ч\ге п$. Оператор выбора скрыт в макросе.
Два примера для базы $2^{64}$ (значения в шестнадцатеричном формате с прямым порядком байтов с пробелом между конечностями):
модуль_значение 000076513ae0b1cd
const_ratio 00000000000229e6 7f4ca82ba3a115f1
получить<0>(I) 00005f0fd669f2c7
получить<1>(I) 000041a1f91ef16f
г 00000000185f2ae8 a455846cb7cf9b49
tmp1 000034bb854f9a8d
tmp3 00000fcebfd55b60
получить<2>(I) 00000fcebfd55b60
модуль_значение 686f4b7702a9c775
const_ratio 0000000000000002 7387d66ffd685b82
получить<0>(I) 536094611fa2b19b
получить<1>(I) 675ef5187093ff63
z 21aac8fcf31d6421 62e675ba16d513f1
тмп1 5287278703394bb1
tmp3 72b1d3d2b9f5e50c
получить<2>(I) 0a42885bb74c1d97
Примечание: для $n\in[2^{63}+1,\,2^{64}-1]$, количество $\шляпа ч$ может переполнить одну конечность, и код в его нынешнем виде перестанет работать. Например. для ввода $f=g=2^{32}$, мы получили $z=2^{64}$ таким образом $\шляпа q=0$ (для любой $n>2^{63}$), таким образом $\шляпа ч=г=2^{64}$ и вывод $0$ а не истинный $ч=2^{64}-n$. В полном источнике есть комментарий: «Класс Modulus представляет неотрицательный целочисленный модуль до 61 бита», поэтому такие проблемы для больших $n$ происходит только тогда, когда вызывающий код ошибается. Плюс, если я правильно понял, $n=2^{60}-2^{14}+1$ является главной целью.
Альтернатива: для чего-то, выполняющего ту же функцию, что и код вопроса, в 4 коротких строках кода вместо 11, для всех $n\in[1,\,2^{64}-1]$, не требующий предварительных вычислений, возможно, быстрее, но совместимый только с последними компиляторами x64+CPU, см. первый из эти фрагменты кода (второй - маленький вариант без ограничения $f\,g<2^{64}\,n$ ). Я не делаю заявлений о постоянном времени любого кода.
