HTTPS-запросы в OpenComputers ещё можно посылать без ухищрений, а вот другие протоколы поверх TLS уже просто так не получится. Это, например, FTP и IRC, но, в принципе, может быть и XMPP (?), NNTP (??), SMTP (????) или даже SIP для совсем безумных. Да и HTTPS-то тоже не без ограничений. Когда история только начиналась лет десять назад, в OC нельзя было указывать свои хедеры. Сейчас можно, но, как я понимаю, всякие хипстерские PATCH-запросики отправить не получится, и попытки это исправить провалились.
В общем, у реализации TLS прямо для OC на Lua польза есть. Так что где-то в прошлом тысячелетии я сделал libtls. Повезло, что в теме не сказать чтобы сильно разбирался и не представлял масштаба начинания. Качество реализации оказалось соответствующим. Поэтому, когда в 2018 году вышла новая версия протокола (1.3), было желание переписать всё с нуля и по-нормальному.
Прошло пять лет. На месте индустрия не стояла и как раз с середины 2010-х годов стала внедрять алгоритмы, обеспечивающие forward secrecy. Этот термин означает, что даже если кто-то набрутфорсит или сольёт ключи от одной TLS-сессии, то в расшифровке последующих сессий ему успех этот никак не поможет, потому что каждый раз ключи генерируются другие. И так уж вышло, что эти алгоритмы работают либо долго, либо на эллиптических кривых.
Как-то в IRC один человек пинганул меня, что очень хотелось бы в libtls эллиптическую криптографию. Разумеется, трогать libtls желания у меня никакого не было, поэтому вместо этого потратил 3 недели и запилил новую библиотеку libtls13. Это реализация TLS 1.3 для OpenComputers, но, кроме того, ещё и библиотека всяких криптографических алгоритмов на чистом Lua. Там и AES-GCM, и ChaCha20-Poly1305, и SHA2, и RSA, и, собственно, эллиптическая криптография: X25519, Ed25519. Когда я зарелизил, этот же человек попросил добавить ещё и ECDSA.
ECDSA — это алгоритм электронной подписи на эллиптических кривых. Отправитель генерит пару ключей: секретный и публичный. Секретный он, хм, держит в секрете, потому что им он подписывает сообщения. Публичный же пускает по ветру, чтобы любой, кто получает сообщение и подпись, мог последнюю проверить. Если кто-то сообщение подменит, то проверка эта провалится. (Либо же этот кто-то накрутил удачи столько, что нашёл коллизию хэшей SHA-256, чему можно только позавидовать: пока таких людей не обнаружилось.)
Тут нужно сказать, что ECDSA в libtls13 уже был, но только на дата-карточке. ECDSA работает с эллиптическими кривыми, и карточка умеет в две: secp256r1 (256-битные ключи) и secp384r1 (384-битные), — из которых в TLS юзабельна была только первая (потому что для 384 бит карточка не тот хэш юзает).
И, естественно, человеку нужна была именно secp384r1.
Вот про то, как я реализовывал на Lua алгоритм ECDSA для кривой secp384r1, и поговорим. Ниже:
определение эллиптических кривых конечные поля и группы реализация модульной арифметики для фиксированного модуля а также операций в группе эллиптической кривой махинации с битами реализация собственно ECDSA (что, на самом деле, не так и интересно) пара слов об ASN.1 — как JSON, только от большого комитета из literally 1984 сколько-то про тестирование Всё это приправлено математикой и сложными терминами, чтоб сложней было понять.
1. Эллиптические кривые
Криптография на эллиптических кривых — это такое очень конкретное приложение абстрактной алгебры. Поэтому, когда гуглится какая-нибудь статья про них, автор либо считает, что читатель уже её основы получил на первом курсе вуза, либо пытается (не сильно успешно) излагать первый курс вуза сам. Мне кажется, что средний игрок в Minecraft, даже с OC, не сильно соотносится с прилежным студентом математической специальности. Поэтому в посте дальше попытаюсь (не сильно успешно) изложить элементы абстрактной алгебры.
Итак. Взглянем-ка на вот этот набор символов.
y2 = x3 + ax + b.
Это уравнение. Тут есть буквы a и b — это какие-то константы. Они берутся от балды и не меняются. Ещё есть буквы x и y — это переменные. Они берутся от балды и вечно меняются.
Уравнений много всяких, но вот конкретно это называется уравнением эллиптической кривой. Кривая там получается, если взять все решения (x, y) уравнения и намалевать ими по плоскости. Вот как это выглядит (пикча с википедии).
Странно, конечно, что в первых трёх пикчах пузырь отчленён от усов, но в глазах математика кривая на каждой пикче ровно одна. Эллипсов нигде не видно. Впрочем, с ними эти кривые почти никак не связаны.
Есть нюанс, правда. Если так получилось, что 4a3 + 27b2 = 0, то математики говорят, что мы налажали. На графике это выражается тем, что кривая самопересекается, имеет острые шипы (как на рисунке выше при a = b = 0) или изолированные одноточечные выбросы. К счастью, здесь a и b выбрали за нас в американском институте, и там всё гладенько.
2. Кольца и поля
Выше все числа были вещественными. В электронных вычислительных машинках работать с ними как-то несподручно. Поэтому сделали два гениальных хода.
Во-первых, в уравнении ограничили a, b, x и y целыми числами. С этим уже попроще, но их как была бесконечность, так и стала. Так что, во-вторых, после каждой операции берётся остаток от деления на какое-то разумное число q. После каждой — это значит после каждого сложения, умножения и возведения в степень, то есть вот так:
(y % q)(y % q) % q = (((((x % q)(x % q)(x % q) % q) + ((a % q)(x % q) % q)) % q) + (b % q)) % q.
Каждый раз писать так влом, поэтому переопределим сложение, умножение и возведение в степень так, чтобы они сами работали по модулю q. И когда теперь я напишу y2 = x3 + ax + b, об этом вспомним и фактически делать будем то, что в большом уравнении.
Думаю, то, что мы переопределили сложение, программиста на Lua сильно смущать не должно. Во-первых, с целыми числами и так все эти операции совершаются по модулю 264 (потому что числа 64-битные), и поэтому плюс не такой, как в школе, где числа бесконечны, уже не впервые встречается. Во-вторых, норкоманы спокойно переопределяют метаметоды __add и __sub и не сильно парятся об этом. Но программисты на Lua не единственные люди, которые такими вещами занимаются. Как и следовало ожидать, первой такой эзотерикой озаботились математики.
Когда мы адекватно определяем сложение и умножение, вместе эти операции они называют кольцом. Само собою, как и эллиптические кривые на эллипсы не похожи, так и кольца эти с фигурой толком не связаны. Все эти имена — многовековое легаси, так что ничего удивительного. Кольцо — это такой абстрактный интерфейс, если по-человечески выражаться, с операциями + и *. Там надо, чтоб ещё кучка правил выполнялась, которые можно загуглить при желании (ну либо на всю жизнь запомнить, как я...). Польза от колец этих одна: они себя ведут приятным образом.
Вот, например, сложение/умножение целых чисел по модулю — это кольцо. И поэтому работать по модулю можно вполне успешно: можем вычитать, выполняются привычные утверждения вроде (x + y)2 = x2 + 2xy + y2 и т. д.
Кольца — это прекрасно, но возникают вопросы с делением. Деление здорового человека выглядит так: если a / b = c и b ≠ 0, то bc = a. То есть деление — обратная операция к умножению. В Lua что /, что // — деление курильщика. Первое нам бесполезно, потому что выдаёт числа с плавающей точкой, когда мы работали с целыми. Второе не обладает свойством выше.
К вопросу подойти можно с тыльной стороны. Поделить на b — то же самое, что умножить на обратное к b. Оно обозначается b −1. В школе нужно просто взять и поделить единицу о b, но вот прям здесь это выглядит подозрительным. Потому что сначала деление определили через умножение на обратное, а потом обратное считаем через деление.
Поэтому математики подсуетились и придумали поле. Это как кольцо, только ещё нужно для каждого ненулевого элемента определить ему обратный по умножению. То есть каждому b, кроме нуля, сопоставить b −1, чтобы b −1b = 1. И тогда можем спокойно делить, как ни в чём не бывало. Вся задача теперь состоит в определении обратного.
И для модульной арифметики возникают проблемы. Оказывается, что, например, по модулю 264, как в Lua, для чётных чисел обратных в принципе не может существовать. Хотя вот для 3, например, обратный элемент имеется и равен −6148914691236517205 (можете проверить). Это означает, что просто так модульная арифметика — кольцо, но не поле. Поэтому, как бы деление ни попытались определить, нормально на всех (ненулевых) числах оно не заработает.
Но как только модуль становится простым числом p, возникают чудеса. Потому что там вот как раз каждый элемент можно обратить. Например, если p = 13, то вот как это выглядит:
1−1 = 1 2−1 = 7 3−1 = 9 4−1 = 10 5−1 = 8 6−1 = 11 7−1 = 2 8−1 = 5 9−1 = 3 10−1 = 4 11−1 = 6 12−1 = 12 Никакой ракетной магии: для каждого числа просто перебираем 12 ненулевых остатков от деления на 13 и ищем тот, умножение с которым даёт остаток 1.
Модульная арифметика с простым модулем p является полем. Поэтому в нём верны все алгебраические преобразования, включая использующие деление. Например, a / b + c / d = (ad + bc) / bd, тоже со школы знакомое, справедливо, даже если все операции делать по модулю простого числа (а деление производить через умножение на обратное). Математики несколько веков без компьютера сидели и от скуки кучу таких преобразований напридумывали, и их все можно юзать в любом поле.
Короче говоря, поле — это круто, поэтому дальше везде все операции будем проводить по модулю простого числа.
3. Группа эллиптической кривой
В первой главе мы взяли уравнение эллиптической кривой и полюбовались на графики. Во второй сказали, что работать будем по модулю, чтоб компуктер не взорвался, и уточнили, что модуль простой, чтоб могли делить и властвовать. Количество текста, которое посвятили делению, несмотря на то что до этого оно нигде не фигурировало, должно наталкивать на мысль, что всё ещё впереди. Поэтому продолжаем углубляться в эллиптическую нору.
Точки эллиптической кривой — это пары чисел (x, y), которые решают уравнение кривой. Вот это уравнение ещё раз, чтоб не скроллили туда-сюда:
y2 = x3 + ax + b.
Теперь давайте придумывать, что с точками можно сделать.
Для начала можно увидеть, что y фигурирует в уравнении единожды и в квадрате. Это значит, что если (x, y) уравнение решает, то и (x, −y) тоже вполне сойдёт, потому что квадрат минус всё равно съест. Ровно поэтому на самой первой пикче все кривые были симметричны горизонтальной оси. Так что первая операция точке P сопоставит ей противоположную, у которой координата y заменена на −y. Эту операцию вполне логично мы обозначим минусом: −P.
Если есть минус, то должен быть и плюс. Это, конечно, не всегда так по жизни, но у нас абстрактные размышления, поэтому можем позволить. Плюс будет брать две точки и выплёвывать третью по какой-то схеме.
Самое простое — это нейтральная точка, или «ноль». Её обозначим буквой O. И скажем, что P + O = O + P = P. Примитивные манёвры: сложение с нулём отдаёт исходную точку. Заодно свяжем и плюс с минусом: P + (−P) = O. Всё, дело за малым — определить плюс для остальных точек и узнать, какую именно точку мы обозначили за O.
Мозг прежде всего тянется складывать точки покоординатно. Типа (a, b) + (c, d) = (a + c, b + d). Удобно, красиво, вот только в корне неверно. Слева от равенства были точки на нашей кривой. Справа от равенства точка очень вряд ли на ней окажется. Так как мы тысячу слов заморачивались с эллиптическими кривыми, такой плюс нас не устраивает.
Обратимся к геометрии. Для этого возьмём рандомно, например, a = −1, b = 1 и построим график.
Ой. По модулю 31 кривая — набор раскиданных точек. Как-то не очень информативно получается. Лучше тогда забыть про модуль и сообразить пока без него.
Вот у нас две рандомных точки. Ну, например, проведём через них прямую.
Вот она, третья точка-то! Можем, она и будет результатом сложения? То есть P + Q = R. Также получим P + R = Q.
Но теперь вспомним про «ноль». Пусть Q = O, так как буквы похожи. Второе равенство говорит, что P + R = O, то есть P = −R. А по первому P = R. Получается, P = −P: любая точка равна противоположной! Это как-то совсем не то, что мы бы хотели.
Решение просто — P + Q = −R. Иначе говоря, сумма трёх точек на прямой равна нулю. Тогда у нас всё будет работать шикарно.
Свяжем геометрию с алгеброй: компьютеру нужно координаты R высчитывать по формуле. Если вспомнить уравнение прямой, подставить координаты P и Q и разрешить относительно R, то получим вот такое:
s = (Q.y − P.y) / (Q.x − P.x),
xʹ = s2 − P.x − P.y,
yʹ = −P.y − s(xʹ − P.x).
Вот (xʹ, yʹ) — это и есть наша сумма.
Формулы всем прекрасны, кроме случая, когда у P и Q координата x одинакова. Потому что тогда в формуле для s будем делить на ноль, а это неприятно.
Формула считает угол наклона прямой между P и Q. Прямая эта называется секущей, а предел секущей, как известно, есть касательная. Поэтому при P.x = Q.x приравняем s углу наклона касательной — то есть производной:
s = (3P.x2 + a) / (2P.y).
Но деление-то тут всё ещё есть, и если y в точке P занулится, то, пока P + P считать будем, делить станем на ноль. Правда, если P.y = 0, то P = −P (мы ведь именно у этой координаты знак меняем). Значит, P + P = O. Тоже всё, в принципе, решается дополнительными проверками.
Кстати, пока мы здесь, в вещественных числах, было бы неплохо понять, а какие координаты у «нулевой» точки O. Ситуация получается очень весёлая, потому что их нет. Их украл товарищ Вейерштрасс. Вот как это вышло. Чтобы посчитать P + P, нам нужно взять касательную к P (см. справа).
Если будем двигать P к левому краю, пересечение с кривой будет всё выше и выше. В конце концов пересечение улетит на бесконечность, а касательная станет вертикальной. Как мы помним, результатом суммы будет O. Получается, нулевая точка O находится на бесконечности, как бы странно это ни звучало. В простых координатах (x, y) бесконечность не уместится. Мы потом придумаем, что с этим делать, когда будем реализовывать всё вышеописанное на Lua.
К слову, о Lua. Пока бы вернуться всё-таки к модульной арифметике. Как вся арканная геометрия нам поможет? А мы просто возьмём те же формулы. Посчитать по ним координаты мы можем, потому что они используют сложение, вычитание, умножение и деление, а все эти операции определены в поле. А результат будет решать уравнение потому, что оно тоже юзает те же операции.
Работа в поле (операции по модулю простого числа) нас спасает просто по всем пунктам. Получается, абстракции, которые математики наплодили, несут какую-то пользу? Удивительно...
Раз уж об абстракциях вспомнили, закинем ещё одну, чтоб объяснить название главы. Точки кривой с операцией их сложения математики назовут группой. Причём, так как порядок сложения точек не влияет на результат (по-умному это называется коммутативностью), группу назовут коммутативной (или абелевой). Поэтому, когда будем реализовывать операции на Lua, применим трюки из теории групп.
Итак, что мы получили.
Можем искать противоположные точки. Можем точки складывать между собой, оставаясь на кривой. Весь этот зоопарк называется группой. И под конец огребли дурное предзнаменование, что это нам ещё понадобится...
4. Арифметика в полях
Пора б уже выдвигаться на поле и пытаться сооружать на Lua всё, что до этого описывали. Для начала познакомимся с кривой, с которой будем работать. Она называется secp384r1, потому что:
определена она группой SEC над полем, образованным модульной арифметиков по модулю простого (p) числа длиной в 384 бит и параметры выбирались рандомно
Как уже говорил, кривая задаётся значениями a и b. У это кривой они такие:
a = −3,
b = 0xb3312fa7 e23ee7e4 988e056b e3f82d19 181d9c6e fe814112 0314088f 5013875a c656398d 8a2ed19d 2a85c8ed d3ec2aef.
Размер константы b внушает уважение, конечно. Осталось ещё p обозначить. Оно вот такое необычное:
p = 2384 − 2128 − 296 + 232 − 1.
Длина модуля — 384 бит, и в одно целое число в языке Lua (всего с 64 битиками) оно никак не влезет. Поэтому эти числа мы побьём на кусочки и будем хранить в таблице. Нужно только понять, какого размера брать куски. Это зависит от того, какие операции собираемся проделывать.
Нам нужны сложение и вычитание. В Lua с этим проблем не возникает, потому что + и - определены для 64-битных чисел.
Нужны умножение и деление. Деление, как уже говорил, — просто умножение на обратное. Обращение чисел, как я потом покажу, мы тоже реализуем с помощью умножения. Таким образом, нужно лишь обеспечить добротное умножение, и всё будет в шоколаде.
В Lua оператор * берёт два 64-битных числа и отдаёт 64-битное произведение. Если немного подумать, то 263 × 263 = 2126, что в 64 бита не очень-то и влезает. Всё, что не влезло в один кусочек результата, нам нужно будет перетащить в следующий. Но для этого нужно знать, что, собственно, перетаскивать. А Lua старшие биты не отдаёт, и узнать не получится.
Поэтому если нужно умножение, то куски должны быть размером не больше 32 бит. Тогда даже самое большое произведение — (232 − 1)2 = 264 − 233 + 1 — влезет в 64 бита.
Пораскинув мозгами, я порезал ещё мельче — длиной в 30 бит. Вот как это выглядит.
Итого число представляется 13 кусками, которые я дальше называть буду словами. В таблице они лежат в порядке номеров бит. Это означает, что нулевой бит исходного числа находится в нулевом бите первого слова, а, например, 333-й — в третьем бите 12-го слова. Такой расклад называется little-endian. Так как 384 на 30 не делится, последнее, тринадцатое слово содержит всего 24 бита.
Почему выбрал 30, а не 32, я объясню, когда будем пилить умножение. А пока начнём с малого — единицы и нуля.
4.1. Нейтральные элементы
Это такое умное название для нуля и единицы. Смотрим на схему выше и соображаем код:
local function fieldZero(a) a = a or {} for i = 1, 13, 1 do a[i] = 0 end return a end local function fieldOne(a) a = a or {} for i = 2, 13, 1 do a[i] = 0 end a[1] = 1 return a end Так как вся криптография — это сплошные числодробилки, забивать память таблицами не очень хочется. Поэтому все наши функции будут результат записывать в таблицу, переданную по ссылке первым аргументом.
Ноль и один — всему голова, но пора бы перейти к чему-то поинтереснее.
4.2. Сложение
Складывать числа будем по-простому: от младших слов к старшим.
Правда, нужно ещё помнить о переносах:
00111111 11111111 11111111 11111111 + 00111111 11111111 11111111 11111111 =================================== 01111111 11111111 11111111 11111110 ^ Тут я сложил два каких-то 30-битных числа. А получил 31-битное. Этот лишний бит надо изъять и прибавить к следующему слову. В результате код выглядит вот так:
local carry = 0 for i = 1, 12, 1 do local word = a[i] + b[i] + carry c[i] = word & 0x3fffffff -- ① carry = word >> 30 -- ② end local word = a[13] + b[13] + carry c[13] = word & 0xffffff -- ③ carry = word >> 24 ① — это константа с тридцатью единицами. С помощью неё в c[i] записываем только 30 младших бит.
② — а всё, что не вместилось, закидываем в carry, который на следующей итерации прибавится.
③ — с последним словом делаем всё то же, только там 24 бита, а не 30.
На данном этапе получается вот так:
a + b = c + 2384 * carry.
Это мы сделали сложение по модулю 2384, а у нас он не такой:
p = 2384 − 2128 − 296 + 232 − 1.
Мы работаем по модулю p, так что
2384 − 2128 − 296 + 232 − 1 = 0
(= здесь — это «равно по модулю».)
Перенесём всё, кроме первого, вправо:
2384 = 2128 + 296 − 232 + 1.
Получается, чтобы избавиться от бита переноса (carry), который у нас мог появиться после сложения, надо его прибавить к 128-му, 96-му и нулевому биту, а из 32-го вычесть. Ориентируемся на схему выше и так и поступаем:
c[2] = c[2] - (carry << 2) -- bit 32 c[4] = c[4] + (carry << 6) -- bit 96 c[5] = c[5] + (carry << 8) -- bit 128 for i = 1, 13, 1 do local word = c[i] + carry c[i] = word & 0x3fffffff carry = word >> 30 carry = carry | -(carry & 1 << 63 >> 30) -- ① end В первых трёх строках мы обработали 128-й, 96-й и 32-й биты. Однако это может попортить ранее сделанную работу: например, если c[4] оказался равен 0x3fffffff — 30 единиц, — то после прибавления единицы c[4] будет уже 31-битным. Поэтому в циклике ниже я снова обрабатываю переносы. Заодно такое ловкое построение кода обработает и нулевой бит.
Тут есть одна грабля, которая обозначена ①. Дело в том, что из 32-го бита мы carry вычитали. Что будет, если c[2] был нулевым до вычитания?
00000000 00000000 00000000 00000000 00000000 00000000 00000000 00000000 - 00000000 00000000 00000000 00000000 00000000 00000000 00000000 00000100 ======================================================================== 11111111 11111111 11111111 11111111 11111111 11111111 11111111 11111100 ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~ Внизу оказалось закодировано -4. Когда в цикле дойдём до второго слова, мы в него запишем младщие 30 бит, как и положено. А потом в carry присвоим 34 старших бита (подчёркнуто в схеме выше):
00000000 00000000 00000000 00000011 11111111 11111111 11111111 11111111 Это то, что мы собираемся прибавить к третьему слову. Но нам надо, наоборот, из третьего слова единицу вычесть, потому что во втором было слишком маленькое значение. Для этого вспомним, что x − 1 = x + (−1), то есть если мы хотим именно складывать, то в carry надо записать −1. И если бы carry был не 64-битным числом, а 34-битным, то именно −1 бы там и оказалось:
11 11111111 11111111 11111111 11111111₂ = -1 Нам нужно лишь расширить carry до 64-битного числа, чтобы получилось так:
11111111 11111111 11111111 11111111 11111111 11111111 11111111 11111111₂ = -1 Для этого самый старший бит здесь (тридцать третий) нужно прокопировать влево, пока не наберём, сколько нужно. Эта операция называется знаковым расширением. Именно этим и занимается строка, отмеченная ①:
carry = carry | -(carry & 1 << 63 >> 30) Сначала мы достаём 33-й бит (1 << 63 << 30 — то же, что и 1 << 33):
00000000 00000000 00000000 00000011 11111111 11111111 11111111 11111111 = carry & 00000000 00000000 00000000 00000010 00000000 00000000 00000000 00000000 = 1 << 63 >> 30 ======================================================================== 00000000 00000000 00000000 00000010 00000000 00000000 00000000 00000000 Затем мы берём противоположное число (такое, сумма с которым даёт ноль). Алгоритм его нахождения прост: перевернуть все биты в числе и прибавить единицу.
00000000 00000000 00000000 00000010 00000000 00000000 00000000 00000000 ~: 11111111 11111111 11111111 11111101 11111111 11111111 11111111 11111111 +1: 11111111 11111111 11111111 11111110 00000000 00000000 00000000 00000000 Наконец, проставляем эти биты в carry:
00000000 00000000 00000000 00000011 11111111 11111111 11111111 11111111 = carry | 11111111 11111111 11111111 11111110 00000000 00000000 00000000 00000000 = -(carry & 1 << 63 >> 30) ======================================================================== 11111111 11111111 11111111 11111111 11111111 11111111 11111111 11111111 = -1 Вуаля. Если бы 33-й бит в carry не стоял, то справа от | оказался бы 0, и вернулось бы то же число. Так что в результате цикл автомагически расставит переносы на словах.
После всех этих операций мы получим либо число, полностью сокращённое по модулю p (то есть от 0 до p − 1), либо почти сокращённое — в пределах от p до 2p − 1. Иначе говоря, в c у нас останется лежать либо (a + b) % p, либо p + (a + b) % p. Чтобы довести процедуру до конца, надо проверить, меньше ли результат модуля p, и вычесть p, если это не так. Но мы этим заморачиваться будем только по необходимости, чтоб экономить ресурсы, потому что складывать мы будем очень много, и почти сокращённых чисел в большинстве случаев будет достаточно.
Вот так складываются 2 числа по фиксированному модулю. Трюки, которые отсюда нужно извлечь:
сложение с переносом перенос старших бит каждого слова в последующее сокращение лишних бит результата по модулю p знаковое расширение
4.3. Вычитание
Вычитать число из другого будем кодом, похожим на предыдущий:
-- "p right-shifted by 9" local prs9 = { -- ① 0x3ffffe00, 0x000007ff, 0x00000000, 0x3fff8000, 0x3ffdffff, 0x3fffffff, 0x3fffffff, 0x3fffffff, 0x3fffffff, 0x3fffffff, 0x3fffffff, 0x3fffffff, 0x1ffffffff, } local function fieldSub(c, a, b) local carry = 0 for i = 1, 12, 1 do local word = a[i] - b[i] + prs9[i] + carry -- ① c[i] = word & 0x3fffffff carry = word >> 30 carry = carry | -(carry & 1 << 63 >> 30) -- ② end local word = a[13] - b[13] + prs9[13] + carry -- ① c[13] = word & 0xffffff carry = word >> 24 c[2] = c[2] - (carry << 2) -- bit 32 c[4] = c[4] + (carry << 6) -- bit 96 c[5] = c[5] + (carry << 8) -- bit 128 for i = 1, 13, 1 do local word = c[i] + carry c[i] = word & 0x3fffffff carry = word >> 30 carry = carry | -(carry & 1 << 63 >> 30) end end Отличия я обозначил. Во-первых, числа, с которыми мы работаем, неотрицательные. Но нужно что-то делать, если вычитаем, например, единицу из нуля.
Для этого функция выше вычисляет не a − b, а a + 29p − b. Этим ничего не сломаю, потому что 29p делится на p без остатка, и при сокращении по модулю он уйдёт. Но 29p — это 393-битное число (у него проставлен 392-й бит), а другие числа много меньше его (мы их держим в диапазоне от 0 до 2p - 1, напоминаю). Поэтому результат точно будет неотрицательным. Девятка здесь выбрана от балды: можно было и двойку поставить, например.
А во-вторых, так как мы вычитаем числа, то может возникнуть та же проблема, что и раньше: в carry может понадобиться записать отрицательное число, чтоб вычесть его из следующего слова. Поступаем так же: производим знаковое расширение после каждого возможного вычитания.
В остальном это то же сложение.
4.4. Умножение
Принцип работы умножения в модульной арифметике — сначала сделать само перемножение, а потом мучительно его сократить по модулю. Перемножать будем по словам «в столбик». В десятичной системе счисления мы каждую цифру одного числа множим на каждую цифру другого числа, а потом всё суммируем:
123 × 321 ====== 3 = 1×3 + 2 = 1×2 + 1 = 1×1 + 6 = 2×3 + 4 = 2×2 + 2 = 2×1 + 9 = 3×3 + 6 = 3×2 + 3 = 3×1 ====== + 3 = 3 + 8 = 2 + 6 + 14 = 1 + 4 + 9 + 8 = 2 + 6 + 3 = 3 ====== 39483 Обычно в записи одну цифру умножают сразу на всё число, потому что в мозгу человека хватает кратосрочной памяти, чтоб запомнить переносы. Таким образом промежуточные результаты группируются по-горизонтали — так, как на схеме посередине (группы разделены пустыми строками).
Но никто не мешает группировать и по-вертикали, как после второй черты. Например, первую строку получим, помножив 1-ю цифру одного числа на 1-ю цифру другого. Вторая строка — сумма произведения 1-й цифры одного числа на 2-ю цифру другого и произведения 2-й цифры на 1-ю. В общем случае будет так: если в числах K цифр, то n-ая строка — сумма всех a[i] * b[n + 1 − i], где i меняется от math.max(1, n − K + 1) до math.min(n, K).
Это были десятичные числа, где a[i] с b[j] в пределах от 0 до 9. Ничто не мешает взять числа покрупнее. В частности, a[i] и b[j] можно сделать 30-битными — по выбранному размеру слова. Формула останется той же самой, а промежуточный результат будет группирован по словам, что очень удобно. В коде это выглядеть будет так:
local d = {} for n = 1, 25, 1 do local dn = 0 for i = math.max(1, n - 12), math.min(n, 13), 1 do local j = n + 1 - i dn = dn + a[i] * b[j] end d[n] = dn end Операция, которую мы проделываем, в математике называют свёрткой. В результате свёртки в d будут лежать промежуточные результаты умножения. Остаётся только перенести старшие биты каждого слова в последующее и провести сокращение по модулю.
У себя в коде я этот цикл развернул, потому что умножений в ECDSA надо проделывать десятки тысяч, и вот такие вложенные циклы там ни к чему. Получилось так:
local a1, a2, a3, a4, a5, a6, a7, a8, a9, a10, a11, a12, a13 = table.unpack(a, 1, 13) local b1, b2, b3, b4, b5, b6, b7, b8, b9, b10, b11, b12, b13 = table.unpack(b, 1, 13) local d = { a1 * b1, a1 * b2 + a2 * b1, a1 * b3 + a2 * b2 + a3 * b1, a1 * b4 + a2 * b3 + a3 * b2 + a4 * b1, a1 * b5 + a2 * b4 + a3 * b3 + a4 * b2 + a5 * b1, a1 * b6 + a2 * b5 + a3 * b4 + a4 * b3 + a5 * b2 + a6 * b1, a1 * b7 + a2 * b6 + a3 * b5 + a4 * b4 + a5 * b3 + a6 * b2 + a7 * b1, a1 * b8 + a2 * b7 + a3 * b6 + a4 * b5 + a5 * b4 + a6 * b3 + a7 * b2 + a8 * b1, a1 * b9 + a2 * b8 + a3 * b7 + a4 * b6 + a5 * b5 + a6 * b4 + a7 * b3 + a8 * b2 + a9 * b1, a1 * b10 + a2 * b9 + a3 * b8 + a4 * b7 + a5 * b6 + a6 * b5 + a7 * b4 + a8 * b3 + a9 * b2 + a10 * b1, a1 * b11 + a2 * b10 + a3 * b9 + a4 * b8 + a5 * b7 + a6 * b6 + a7 * b5 + a8 * b4 + a9 * b3 + a10 * b2 + a11 * b1, a1 * b12 + a2 * b11 + a3 * b10 + a4 * b9 + a5 * b8 + a6 * b7 + a7 * b6 + a8 * b5 + a9 * b4 + a10 * b3 + a11 * b2 + a12 * b1, a1 * b13 + a2 * b12 + a3 * b11 + a4 * b10 + a5 * b9 + a6 * b8 + a7 * b7 + a8 * b6 + a9 * b5 + a10 * b4 + a11 * b3 + a12 * b2 + a13 * b1, a2 * b13 + a3 * b12 + a4 * b11 + a5 * b10 + a6 * b9 + a7 * b8 + a8 * b7 + a9 * b6 + a10 * b5 + a11 * b4 + a12 * b3 + a13 * b2, a3 * b13 + a4 * b12 + a5 * b11 + a6 * b10 + a7 * b9 + a8 * b8 + a9 * b7 + a10 * b6 + a11 * b5 + a12 * b4 + a13 * b3, a4 * b13 + a5 * b12 + a6 * b11 + a7 * b10 + a8 * b9 + a9 * b8 + a10 * b7 + a11 * b6 + a12 * b5 + a13 * b4, a5 * b13 + a6 * b12 + a7 * b11 + a8 * b10 + a9 * b9 + a10 * b8 + a11 * b7 + a12 * b6 + a13 * b5, a6 * b13 + a7 * b12 + a8 * b11 + a9 * b10 + a10 * b9 + a11 * b8 + a12 * b7 + a13 * b6, a7 * b13 + a8 * b12 + a9 * b11 + a10 * b10 + a11 * b9 + a12 * b8 + a13 * b7, a8 * b13 + a9 * b12 + a10 * b11 + a11 * b10 + a12 * b9 + a13 * b8, a9 * b13 + a10 * b12 + a11 * b11 + a12 * b10 + a13 * b9, a10 * b13 + a11 * b12 + a12 * b11 + a13 * b10, a11 * b13 + a12 * b12 + a13 * b11, a12 * b13 + a13 * b12, a13 * b13, } Вот тут можно поговорить, наконец, о том, почему слова именно 30-битные. Тринадцатый элемент — вот этот:
a1 * b13 + a2 * b12 + a3 * b11 + a4 * b10 + a5 * b9 + a6 * b8 + a7 * b7 + a8 * b6 + a9 * b5 + a10 * b4 + a11 * b3 + a12 * b2 + a13 * b1, Это сумма 13 произведений 30-битных чисел. Каждое из чисел не превышает 230 − 1, поэтому каждое отдельное произведение не больше 260 − 231 + 1. А значит, сумма меньше либо равна 13 × (260 − 231 + 1) = 0xcffffff98000000d. В этом числе 64 бита, поэтому промежуточный результат вмещается без переполнения. Если же слова делить на 31 бит, то d[13] мог легко переполниться. Поэтому ограничился тридцатью.
Это был первый этап умножения. Дальше предстоит 768-битное число сократить до 384 бит. Начинаем с известных технологий — переноса старших бит:
local carry = 0 for i = 1, 25, 1 do local word = d[i] + carry d[i] = word & 0x3fffffff carry = word >> 30 end Двадцать шестое слово результата будет храниться в carry. Теперь надо избавиться от 13 старших слов, чтобы у нас остались только младшие 13 слов результата. Потому что их мы уже сокращать по модулю умеем (см. сложение, если не умеем). Пользуемся тем же трюком, что и раньше:
2384 = 2128 + 296 − 232 + 1,
поэтому (домножим на 26 + 30k):
2390 + 30k = 2136 + 30k + 2102 + 30k − 238 + 30k + 26.
После небольших перегруппировок получим, наконец, такое:
230(13 + k) = 230(13 + k − 9) + 14 + 230(13 + k − 10) + 12 − 230(13 + k − 12) + 8 + 230(13 + k − 13) + 6.
Страшно, да? А вот как это будет в коде:
local i = 26 repeat d[i - 9] = d[i - 9] + (carry << 14) d[i - 10] = d[i - 10] + (carry << 12) d[i - 12] = d[i - 12] - (carry << 8) d[i - 13] = d[i - 13] + (carry << 6) i = i - 1 carry = d[i] until i == 13 По окончании у нас остаются 13 слов (390 бит), которые нужно свести до 384 бит. Для начала снова проделаем переносы, потому что мы активно в эти слова писали (и активно из них вычитали — не забываем знаковое расширение):
carry = 0 for i = 1, 12, 1 do local word = d[i] + carry c[i] = word & 0x3fffffff carry = word >> 30 carry = carry | -(carry & 1 << 63 >> 30) end local word = d[13] + carry c[13] = word & 0xffffff carry = word >> 24 carry = carry | -(carry & 1 << 63 >> 24) И после этого, наконец, последний раз сократим результат по модулю. Здесь всё абсолютно то же самое, что и в конце сложения:
c[2] = c[2] - (carry << 2) -- bit 32 c[4] = c[4] + (carry << 6) -- bit 96 c[5] = c[5] + (carry << 8) -- bit 128 for i = 1, 13, 1 do local word = c[i] + carry c[i] = word & 0x3fffffff carry = word >> 30 carry = carry | -(carry & 1 << 63 >> 30) end И на этом всё. Произведение a и b лежит в c и не превышает 2p − 1.
4.5. Возведение в квадрат
В эллиптической криптографии умножений в поле делается очень много. Но значительная часть умножений на самом деле является возведением в квадрат (a * a). В этом случае a[i] * b[j] = a[j] * b[i] (индексы меняются местами, а результат тот же, потому что a = b), и из 169 64-битных умножений почти половина (78) излишня, потому что вычисляет уже посчитанное значение.
Поэтому имеет смысл сделать отдельную функцию для возведения в квадрат. Её код отличаться от прошлой будет только в двух местах: там нет распаковывания b, а свёртка (d) считается немного иначе:
local d = { a1 * a1, a1 * a2 << 1, (a1 * a3 << 1) + a2 * a2, a1 * a4 + a2 * a3 << 1, (a1 * a5 + a2 * a4 << 1) + a3 * a3, a1 * a6 + a2 * a5 + a3 * a4 << 1, (a1 * a7 + a2 * a6 + a3 * a5 << 1) + a4 * a4, a1 * a8 + a2 * a7 + a3 * a6 + a4 * a5 << 1, (a1 * a9 + a2 * a8 + a3 * a7 + a4 * a6 << 1) + a5 * a5, a1 * a10 + a2 * a9 + a3 * a8 + a4 * a7 + a5 * a6 << 1, (a1 * a11 + a2 * a10 + a3 * a9 + a4 * a8 + a5 * a7 << 1) + a6 * a6, a1 * a12 + a2 * a11 + a3 * a10 + a4 * a9 + a5 * a8 + a6 * a7 << 1, (a1 * a13 + a2 * a12 + a3 * a11 + a4 * a10 + a5 * a9 + a6 * a8 << 1) + a7 * a7, a2 * a13 + a3 * a12 + a4 * a11 + a5 * a10 + a6 * a9 + a7 * a8 << 1, (a3 * a13 + a4 * a12 + a5 * a11 + a6 * a10 + a7 * a9 << 1) + a8 * a8, a4 * a13 + a5 * a12 + a6 * a11 + a7 * a10 + a8 * a9 << 1, (a5 * a13 + a6 * a12 + a7 * a11 + a8 * a10 << 1) + a9 * a9, a6 * a13 + a7 * a12 + a8 * a11 + a9 * a10 << 1, (a7 * a13 + a8 * a12 + a9 * a11 << 1) + a10 * a10, a8 * a13 + a9 * a12 + a10 * a11 << 1, (a9 * a13 + a10 * a12 << 1) + a11 * a11, a10 * a13 + a11 * a12 << 1, (a11 * a13 << 1) + a12 * a12, a12 * a13 << 1, a13 * a13, } Я тут ещё упоролся и вместо умножения на 2 побитово сдвигаю влево.
4.6. Обращение чисел
Как уже говорил, деление мы заменяем умножением на обратный элемент в поле. Для поиска обратного элемента по модулю простого числа у нас есть 2 варианта: либо применять расширенный алгоритм Евклида, либо воспользоваться малой теоремой Ферма. Проще всего второе, к тому же время его работы не зависит от входного числа, что важно в криптографии.
Теорема это выглядит так: если a ≠ 0 и p простое, то a p − 1 = 1. Так как a p − 1 — это a × a p − 2, получается, что a p − 2 — обратный элемент. Значит, нам нужно всего лишь умножить a с собою p − 2 раз.
Тут возникает проблема, потому что p − 2 — это слегка большое число со 116 знаками в десятичной записи. Если влоб работать на процессоре Рутмастера, который делает 10 миллиардов умножений в секунду, то обращение займёт всего 125 унтригинтиллионов лет. Достойная работа для вычислительного устройства.
Ну на самом деле нет, есть варианты сильно лучше. Посмотрим пока на более приземлённые примеры: например, как мы можем возвести число в 45-ую степень? А вот, например:
b = a -- b = a b = b^2 -- b = a² b = b^2 * a -- b = a⁵ b = b^2 * a -- b = a¹¹ b = b^2 -- b = a²² b = b^2 * a -- b = a⁴⁵ Вместо 45 умножений у нас всего 5 возведений в степень и 3 умножения. Алгоритм для такого называется square-and-multiply, и мы его подробно изучим позднее. Пока стоит отметить, что если мы позволим раскошелиться на временные переменные, то от одного умножения можно избавиться:
b = a -- b = a b = b^2 -- b = a² b = b^2 * a -- b = a⁵ c = b^2 -- c = a¹⁰ c = c^2 -- c = a²⁰ b = c^2 * b -- b = a⁴⁵ Для экспонент на сто порядков больше тоже есть варианты сэкономить до приличия. Алгоритм square-and-multiply вышеупомянутый позволит возвести в степень 39402006196394479212279040100143613805079739270465446667948293404245721771496870329047266088258938001861606973112317 всего за 383 возведения в квадрат и 318 умножений.
Но мы вместо этого воспользуемся тулзой addchain. Получится что-то такое:
local function fieldInvert(b, a) local t1, t2, t3, t4 = {}, {}, {}, {} fieldSq(b, a) fieldMul(b, a, b) fieldSq(b, b) fieldMul(t2, a, b) fieldSq(b, t2) fieldSq(b, b) fieldSq(b, b) fieldMul(b, t2, b) fieldSq(t1, b) for i = 1, 5, 1 do fieldSq(t1, t1) end fieldMul(t1, b, t1) fieldSq(t3, t1) for i = 1, 11, 1 do fieldSq(t3, t3) end fieldMul(t1, t1, t3) for i = 1, 6, 1 do fieldSq(t1, t1) end fieldMul(b, b, t1) fieldSq(t1, b) fieldMul(t3, a, t1) fieldSq(t1, t3) fieldMul(t1, a, t1) fieldSq(t4, t1) for i = 1, 30, 1 do fieldSq(t4, t4) end fieldMul(t3, t3, t4) fieldSq(t4, t3) for i = 1, 62, 1 do fieldSq(t4, t4) end fieldMul(t3, t3, t4) fieldSq(t4, t3) for i = 1, 125, 1 do fieldSq(t4, t4) end fieldMul(t3, t3, t4) for i = 1, 3, 1 do fieldSq(t3, t3) end fieldMul(t2, t2, t3) for i = 1, 33, 1 do fieldSq(t2, t2) end fieldMul(t1, t1, t2) for i = 1, 94, 1 do fieldSq(t1, t1) end fieldMul(b, b, t1) fieldSq(b, b) fieldSq(b, b) fieldMul(b, a, b) end Функция слегка гигантская, но ничего ракетоподобного в ней не обнаруживается. Вместо 318 умножений здесь их всего 15 (возведений в квадрат осталось столько же), так что заморачивались не зря. Хотя всё равно дорогая операция.
Итак, теперь мы можем складывать и вычитать числа по модулю p, а ещё умножать и даже делить. Славно поработали. Теперь, наконец, можно вернуться к тому, с чего начинали, то есть к эллиптическим кривым.
5. Реализация операций на эллиптической кривой
Непосредственно операций в поле, которые мы муторно реализовывали всю прошлую главу, ECDSA вообще не юзает. Но она выполняет операции над точками на эллиптической кривой, а они уже задействуют наше поле.
Формулы товарища Вейерштрасса, которые мы вывели в третьей главе, хороши, но нам не очень подходят. Причин на то две.
Во-первых, там надо дофига делить. Складывать точки нужно будет сотни раз, и каждый раз высчитывать обратный элемент — очень спорная затея. Это всё же самая дорогая операция, что у нас есть.
Во-вторых, нет нуля. Тут я даже комментировать не буду.
В общем, половину третьей главы вы читали зря, потому что юзать мы будем что-то совсем другое. Вместо того, что хранить координаты точки (x, y), мы будем таскать триплет (X, Y, Z). Каждая из больших букв — это 384-битное число из поля (то есть работаем по модулю p). Соотносятся эти представления так:
x = X / Z 2,
y = Y / Z 3.
То есть, чтобы получить тройку по (x, y), можно просто засетить X = x, Y = y и Z = 1. Чтоб вернуться, нужно обратить Z, взять квадрат, домножить на X (получаем x) и на Z −1Y (получаем y).
В принципе, можно увидеть, что для одной и той же точки (x, y) координаты (X, Y, Z) мы можем взять совершенно разные в зависимости от Z. Преимущество здесь в том, что для сложения точек по (X, Y, Z)-координатам нам вообще не нужно дорогущее деление. Оно там потвикает как-то X, Y, Z, и Z уже может стать не единичным, но пока нам не нужны конкретные x и y, никаких делений проводить не потребуется. Сделаем мы его лишь единожды — когда закончим работу с кривой.
А ещё эта тройка позволяет кодировать точку на бесконечности...
5.1. Нейтральный элемент
Вспомним исходное уравнение кривой:
y2 = x3 + ax + b.
Подставим x = X / Z 2 и y = Y / Z 3 в него:
Y 2 / Z 6 = X 3 / Z 6 + aX / Z 2 + b.
Домножим обе части на Z 6:
Y 2 = X 3 + aXZ 4 + bZ 6.
У точки на бесконечности координата Z равна нулю (из-за чего x и y как раз в бесконечность улетают). Тогда от уравнения останется лишь это:
Y 2 = X 3.
В общем, как и у любой другой точки, точку на бесконечности тоже можно представить кучей способов. Например, (1, 1, 0). Я у себя даже вообще все координаты занулял, но здесь давайте всё же возьмём эту тройку как наш канонический нейтральный элемент.
local function groupJacobianZero(q) q = q or {} q[1] = fieldOne() -- X q[2] = fieldOne() -- Y q[3] = fieldZero() -- Z return q end Хранить координаты эти будем в таблице по порядку: X, Y, Z. В функции написано Jacobian, потому что так называются координаты, которые мы используем.
Вооружившись нулём, пришло время сделать что-нибудь полезное. Например, сложить-таки точки.
5.2. Сложение точек
В интернетах есть целый сайт с формулами для вычисления точек на эллиптических кривых. На нём можно найти вот эту запись. Формула не то чтобы очевидная, но нас особо тонкости её не волнуют.
Что волнует, так это её цена: это 11 умножений, 5 возведений в квадрат, 9 сложений/вычитаний и 4 умножения на константу 2 (сдвиг влево на один бит). Они там отсортированы в порядке увеличения числа умножений как самой дорогой операции, и конкретно эта формула торчит сверху.
Там, конечно, есть и ещё более дешёвые, но у них всех есть какие-то предположения к входным точками, например что Z = 1, и нам не подходят. Поэтому берём эту формулу и переписываем её на Lua.
local function groupJacobianAddUnchecked(d, p, q) local x1, y1, z1 = p[1], p[2], p[3] local x2, y2, z2 = q[1], q[2], q[3] local x3, y3, z3 = d[1], d[2], d[3] local z1z1, z2z2, u1, u2, s1, s2, h, i, j, r, v = {}, {}, {}, {}, {}, {}, {}, {}, {}, {}, {} fieldSq(z1z1, z1) -- Z₁Z₁ = Z₁² fieldSq(z2z2, z2) -- Z₂Z₂ = Z₂² fieldMul(u1, x1, z2z2) -- U₁ = X₁ * Z₂Z₂ fieldMul(u2, x2, z1z1) -- U₂ = X₂ * Z₁Z₁ fieldMul(s1, z2, z2z2) -- Z₂³ fieldMul(s1, y1, s1) -- S₁ = Y₁ * Z₂³ fieldMul(s2, z1, z1z1) -- Z₁³ fieldMul(s2, y2, s2) -- S₂ = Y₂ * Z₁³ fieldAdd(z3, z1, z2) -- Z₁ + Z₂ fieldSub(h, u2, u1) -- H = U₂ - U₁ fieldMul2(i, h) -- 2 * H fieldSq(i, i) -- I = (2 * H)² fieldMul(j, h, i) -- J = H * I fieldSub(r, s2, s1) -- S₂ - S₁ fieldMul2(r, r) -- r = 2 * (S₂ - S₁) fieldMul(v, u1, i) -- V = U₁ * I local t2 = {} fieldMul2(t2, v) -- 2 * V fieldSq(x3, r) -- r² fieldSub(x3, x3, j) -- r² - J fieldSub(x3, x3, t2) -- X₃ = r² - J - 2 * V fieldMul(t2, s1, j) -- S₁ * J fieldMul2(t2, t2) -- 2 * S₁ * J fieldSub(y3, v, x3) -- V - X₃ fieldMul(y3, r, y3) -- r * (V - X₃) fieldSub(y3, y3, t2) -- Y₃ = r * (V - X₃) - 2 * S₁ * J fieldSq(z3, z3) -- (Z₁ + Z₂)² fieldSub(z3, z3, z1z1) -- (Z₁ + Z₂)² - Z₁Z₁ fieldSub(z3, z3, z2z2) -- (Z₁ + Z₂)² - Z₁Z₁ - Z₂Z₂ fieldMul(z3, z3, h) -- Z₃ = ((Z₁ + Z₂)² - Z₁Z₁ - Z₂Z₂) * H return fieldZeroFlag(r) end Здесь ещё юзаются 2 новые функции. fieldMul2 производит умножение на 2. Реализовано это на основе сложения, где a[i] + b[i] просто заменяется на (a << 1).
И другая функция — fieldZeroFlag. Она возвращает 1, если в переменной r записан ноль, либо 0 в противном случае. Так как все числа у нас в промежутке от 0 до 2p − 1, узнать, равно ли число нулю, можно, например, так:
local function fieldReduceQuick(b, a) -- a minus p local amp = {} local borrow = 0 for i = 1, 13, 1 do local word = a[i] - p[i] - borrow amp[i] = word & 0x3fffffff borrow = word >> 63 end fieldCopy(b, a) fieldCmov(b, amp, borrow ~ 0x1) return borrow end local function fieldZeroFlag(a) local b = {} fieldReduceQuick(b, a) local bits = 0 for i = 1, 13, 1 do bits = bits | b[i] end return -bits >> 63 ~ 0x1 end Сначала приводим число к диапазону от 0 до p − 1, для чего по-простому вычитаем из числа p. Выбираем из результата вычитания и исходного числа неотрицательное, а другое выкидываем. Затем проходимся по всем битам в выбранном числе, склеивая их трубою (|) в переменную bits. Если bits нулевой, то −bits тоже будет нулевым; если же bits ненулевой, то −bits будет отрицательным. Поэтому если извлечь знаковый бит и инвертировать его (при помощи XOR: ~ 0x1), получим как раз то, что и нужно было.
Могли заметить, что код несколько громоздкий. Например, всё, что делает fieldZeroFlag, — проверяет, является ли какое-либо из слов ненулевым. Это можно сделать сильно проще вот так:
for i = 1, 13, 1 do if b[i] ~= 0 then return 0 end end return 1 Но тогда время исполнения этого цикла будет зависеть от исходного числа. Так, конечно, можно делать, но с точки зрения практической криптографии это нехорошо. Время работы криптографических алгоритмов не должно зависеть от значения секретных данных, потому что никто не мешает замерить это время и получить информацию о содержимом секретов. Если же оно не зависит, то замеры никакой информации раскрыть не могут. Такие алгоритмы называются работающими за константное время.
В посте у меня есть только алгоритм проверки подписи, а не генерации её. Проверка же работает только с публичными значениями. Поэтому конкретно в этом случае, действительно, можно было сделать и по-простому. Но когда я реализовывал, я планировал ещё и генерацию сделать, которой нужен уже секретный ключ. Поэтому все приводимые алгоритмы работают за константное время, хотя на реализацию генерации подписи я к концу подзабил.
Вернёмся к сложению, однако. Давным-давно, в другой главе у нас были две различные формулы для сложения двух точек. Даже в координатах (X, Y, Z) от этого мы никуда не ушли: формула выше, которую реализовали для сложения двух точек, выдаёт правильный результат только в том случае, если P и Q — две различные ненулевые точки.
Проверить, что точки ненулевые, весьма просто: у них должна быть ненулевая координата Z. Но как можно убедиться в том, что они различны, ведь одну и ту же точку можно записать в координатах (X, Y, Z) огромным числом способов? Для этого придётся снова позаниматься алгеброй.
Сначала поразмышляем в обратном направлении: что будет, если P = Q? Тогда, вполне очевидно, x1 = x2 и y1 = y2. Выразим это через X, Y, Z:
x1 = x2 ⇒ X1 / Z12 = X2 / Z22,
y1 = y2 ⇒ Y1 / Z13 = Y2 / Z23.
Вследствие этого в формуле выше получим, что
S2 = Y2Z13 = Z13Z23Y1 / Z13 = Y1Z23 = S1.
Тогда r = S2 − S1 = 0. Аналогично можно получить, что U1 = U2, из-за чего H обнуляется, и поэтому Z3 также будет нулевым. То, что Z3 нулевое, значит, что результат суммы — точка на бесконечности.
Теперь попробуем изменить порядок мысли. Пусть мы провели сложение двух ненулевых P и Q и получили точку на бесконечности и r = 0. Можем ли тогда сказать, что P = Q, и формулу на самом деле применять было нельзя?
Из того, что r = 0, мы получаем, что S1 = S2. Z1 и Z2 не равны нулю (входные точки ненулевые), поэтому:
Y1Z23 = Y2Z13 ⇒ Y1 / Z13 = Y2 / Z23,
то есть y1 = y2.
Теперь разберёмся с иксами. Формула для Z3 такая: ((Z1 + Z2)2 − Z12 − Z22) × H. Так как мы работаем в поле, можно спокойно раскрыть скобки внутри:
Z3 = (Z12 + 2Z1Z2 + Z22 − Z12 − Z22) × H = 2Z1Z2H.
Опять вспоминаем, что Z1 и Z2 ненулевые, а результат нулевой. Можем ли тогда сказать, что H = 0? Да, можем, потому что работаем в поле, где это выполняется.
Тут я б хотел отметить, что необходимо именно поле, а не кольцо. Допустим, мы работали бы по модулю составного числа 15. Но тогда могли бы взять два ненулевых числа 3 и 5, перемножить и словить ноль в ответе! В поле такого невозможно. К счастью, у нас все операции по модулю простого числа p, поэтому это именно поле.
Итак, получили, что H = 0. Это означает, что U1 = U2. Подставим формулы для U1 и U2:
X1Z22 = X2Z12 ⇒ X1 / Z12 = X2 / Z22.
Слева x1, справа x2, и они равны. Победа: мы доказали, что r = 0 и Z3 = 0 означает, что x1 = x2 и y1 = y2.
Алгоритм сложения тогда выглядеть будет так:
Если P = O, то вернуть Q. Если Q = O, то вернуть P. Провести сложение по формуле выше. Если функция вернула 1 (то есть r = 0) и результат — точка на бесконечности, то результат неверный и нужно пересчитать его по формуле удвоения точки. Иначе вернуть результат (он корректен).
В целом всё закономерно, но нужна формула удвоения. Давайте реализуем и её.
5.3. Удвоение точек
Для удвоения (вычисления P + P) формулы тоже приведены на том же сайте. Юзать будем вот такую. В коде она выглядит так:
local function groupJacobianDouble(q, p) local x1, y1, z1 = p[1], p[2], p[3] local x3, y3, z3 = q[1], q[2], q[3] local delta = {} fieldSq(delta, z1) -- delta = Z₁² local gamma = {} fieldSq(gamma, y1) -- gamma = Y₁² local beta = {} fieldMul(beta, x1, gamma) -- beta = X₁ * gamma local x1mdelta = {} local alpha = {} fieldSub(x1mdelta, x1, delta) -- X₁ - delta fieldAdd(alpha, x1, delta) -- X₁ + delta fieldMul(alpha, x1mdelta, alpha) -- (X₁ - delta) * (X₁ + delta) fieldMul3(alpha, alpha) -- alpha = 3 * (X₁ - delta) * (X₁ + delta) fieldAdd(z3, y1, z1) -- Y₁ + Z₁ fieldSq(z3, z3) -- (Y₁ + Z₁)² fieldSub(z3, z3, gamma) -- (Y₁ + Z₁)² - gamma fieldSub(z3, z3, delta) -- Z₃ = (Y₁ + Z₁)² - gamma - delta fieldMul4(y3, beta) -- 4 * beta fieldMul2(beta, y3) -- 8 * beta fieldSq(x3, alpha) -- alpha² fieldSub(x3, x3, beta) -- X₃ = alpha² - 8 * beta fieldSq(gamma, gamma) -- gamma² fieldMul8(gamma, gamma) -- 8 * gamma² fieldSub(y3, y3, x3) -- 4 * beta - X₃ fieldMul(y3, alpha, y3) -- alpha * (4 * beta - X₃) fieldSub(y3, y3, gamma) -- Y₃ = alpha * (4 * beta - X₃) - 8 * gamma² end Здесь функции нужны ещё fieldMul3, fieldMul4 и fieldMul8. Их всех можно определить на основе fieldMulWord, которая умножает 384-битное число на одно 30-битное слово по модулю p:
local function fieldMulWord(c, a, b) local carry = 0 for i = 1, 12, 1 do local word = a[i] * b + carry c[i] = word & 0x3fffffff carry = word >> 30 end local word = a[13] * b + carry c[13] = word & 0xffffff carry = word >> 24 c[2] = c[2] - (carry << 2) -- bit 32 c[4] = c[4] + (carry << 6) -- bit 96 c[5] = c[5] + (carry << 8) -- bit 128 for i = 1, 13, 1 do local word = c[i] + carry c[i] = word & 0x3fffffff carry = word >> 30 carry = carry | -(carry & 1 << 63 >> 30) end end По сути это то же сложение, только вместо a[i] + b[i] пишем a[i] * b. Вместо 169 умножений, как в fieldMul, здесь их нужно всего 13, поэтому в этих заморочках имеется смысл.
Формула удвоения никаких особых случаев, как при сложении, не имеет и без проблем удваивает любую входную точку. Она ещё и дешевле, чем сложение. Но в данном случае формулу удвоения нужно было реализовывать для обеспечения корректности, а не простой оптимизации.
Люди кучу раз пытались реализовывать криптографию на эллиптических кривых и постоянно спотыкались о то, что формула сложения выдаёт неверный результат, если P = Q. Ситуация эта не то чтобы очень частая, поэтому найти проблему трудно, и она проскакивает незамеченной, пока какой-нибудь заинтересованный человек не сломает из-за неё всю криптосистему.
На помощь приходят «полные» формулы сложения, которые, хотя и дорогие по сравнению с обычными формулами, выдают корректный результат на всех входных данных. Таким образом, с ними все эти проблемы избегаются. Найти их можно, например, в этой статье: https://eprint.iacr.org/2015/1060.pdf. В качестве упражнения можно попробовать реализовать их самому.
5.4. Вычитание точек
Последняя базовая операция, которая нам нужна в группе эллиптической кривой, — это вычитание одной точки из другой. Его можно реализовать при помощи сложения с противоположной точкой. У противоположной точки другой знак у координаты y — и, следовательно, Y. В модульной арифметике «поменять знак» значит просто вычесть y из нуля. Поэтому вычитание выглядит следующим образом:
local function groupJacobianSub(d, r, q) local mq = {q[1], {}, q[3]} fieldSub(mq[2], fieldZero(), q[2]) groupJacobianAdd(d, r, mq) end Теперь, когда все основные операции у нас реализованы, наконец можно поговорить о ECDSA.
6. Умножение на скаляр. Приватные ключи
С точками прямо сейчас мы можем совершать две вещи: обнулять и суммировать. Инструменты скудные, поэтому попробуем что-нибудь изобразить помощнее.
Для этого возьмём какую-нибудь точку P на кривой и будем её прибавлять к самой себе несколько раз: P + P, P + P + P и т. д. Чтоб не взрываться с записи, обозначать такое мы будем дальше так: [k]P, где k — какое-то число, — это k точек P, сложенные вместе. Операция эта, как выясняется, не бессмысленна, потому что затрагивает много вопросов о жизни и смерти циклических групп.
Вот что я имею в виду. В нашей группе точек явно не бесконечно. Как минимум, различные точки должны иметь различные координаты, а различных координат существует только p2. Причём большинство из них не будут решать уравнение, так что точек на кривой ещё меньше.
Получается, когда мы будем увеличивать k в [k]P, рано или поздно мы начнём повторяться. Причём, если подумать, то первой повторится точка P, потом [2]P и так далее. Более того, повторению P будет предшествовать получение точки на бесконечности. Значит, есть какое-то такое ненулевое N, что [N]P = O, и потом снова [N + 1]P = P, [N + 2]P = [2]P, ... То есть все эти сложения цикличны.
Минимальное такое N в математике называется порядком элемента. Равно ли N количеству точек на кривой (это количество называется порядком группы)? То есть будем ли мы повторяться с безысходности, когда пробежим вообще все существующие точки?
В общем случае — не факт. Но группа эллиптической кривой secp384r1 является циклической. Это значит, что на кривой есть такая точка G, порядок которой таки равен порядку группы. Она называется порождающим элементом либо генератором группы. Прибавляя генератор к самому себе, мы каждый раз будем получать разные точки, пока не обойдём всю кривую.
Обзовём, как и в стандарте ECDSA, порядок группы буквой n. В случае secp384r1 это простое число, слегка меньшее, чем p. Координаты генератора G на кривой тоже даются в стандарте.
Ключевая вещь, благодаря чему эллиптические кривые вообще несут какую-то пользу в криптографии, такая. Для любой точки P можем набрать достаточно генераторов в сумме, чтобы её получить, то есть существует k такое, что [k]G = P. По k найти P очень просто, и мы этим дальше будем как раз заниматься. А вот обратно, — зная P, найти k — дико сложно. Единственные варианты, которые придумали, лишь слегка умнее перебора втупую.
Поэтому в ECDSA рандомное k от 1 до n − 1 является приватным ключом и никому не раскрывается. А по нему вычисляется Q = [k]G, координаты которой будут уже публичным ключом, известными всему миру.
Как уже писал ранее, на генерацию ключей и подписей я подзабил, когда дошёл до этого, и реализовал только проверку подписи. А там вместо вычисления [k]G нужно считать [u]G + [v]Q. Можно было бы, конечно, реализовать отдельно алгоритм для вычисления [k]P с любых k и P, потом применить его дважды и сложить результаты, но я не стал. Заместо этого я реализовал отдельный специализированный алгоритм, который вычисляет сразу всё выражение [u]G + [v]Q для заданных u, v и Q.
6.1. Алгоритмы умножения на скаляр
Итак, нужно точку P скалярно множить на скаляр. Скаляр увесистый — 384 бит в ширину. Пока будем одно число так влоб суммировать, перед глазами новемвигинтиллион вселенных успеют до текущего состояния с нуля развиться. Хм, где-то мы уже с такой проблемой боролись...
В общем, я бы мог тут упороться с мультипликативными группами поля и дискретными логарифмами. Но проще будет сказать, что все те методы, с помощью которых возводили числа в гигантские степени по модулю, могут быть применены и для умножения точек кривой на скаляр. Например, я там говорил про алгоритм square-and-multiply («взять квадрат и помножить»). Здесь это будет double-and-add («удводить и сложить»). Давайте всё-таки изучим, что это за алгоритм такой.
Пусть наше k не 384, а 5 бит в ширину и равен двадцати двум — либо 101102 в двоичной записи. Для начала мы обнуляем результат, после чего начинаем обходить биты в числе k от старшего ко младшему. Каждый раз при переходе к следующему биту результат удваиваем, после чего, если бит попался единичный, ещё и прибавляем к нему P.
Вот что это вытворит для нашего k:
k = 10110 result = O 1. удваиваем результат: result = [2]result = [2]O = O 2. получаем следующий бит: k = 10110 ↳ 1 3. он единичный, поэтому прибавляем к результату P: result = result + P = O + P = P result = [1₂]P 1. удваиваем результат: result = [2]result = [2]P = [2]P 2. получаем следующий бит: k = 10110 ↳ 0 3. он нулевой, так что ничего не делаем result = [10₂]P 1. удваиваем результат: result = [2]result = [2][2]P = [4]P 2. получаем следующий бит: k = 10110 ↳ 1 3. он единичный, поэтому прибавляем к результату P: result = result + P = [4]P + P = [5]P result = [101₂]P 1. удваиваем результат: result = [2]result = [2][5]P = [10]P 2. получаем следующий бит: k = 10110 ↳ 1 3. он единичный, поэтому прибавляем к результату P: result = result + P = [10]P + P = [11]P result = [1011₂]P 1. удваиваем результат: result = [2]result = [22]P 2. получаем следующий бит: k = 10110 ↳ 0 3. он нулевой, так что ничего не делаем result = [10110₂]P После третьего шага я каждый раз показываю скаляр в результате в двоичной записи. Битики в скаляре всегда совпадают с битиками, которые мы до этого момента прочли. В конце алгоритма, соответственно, битики в скаляре совпадают с k, и мы получаем корректный результат.
Алгоритм в целом уже неплох, но сейчас надо делать столько сложений, сколько в скаляре бит, а их там до 384 может быть. Хотелось бы уменьшить объём работы. Идея состоит в том, чтобы вместо прибавления одного бита каждый раз, докидывать их сразу пачкой в результат. Сделать это можно так:
Выбираем размер окна — w бит. Множим точку P на нечётные скаляры от 3 до 2w − 1: получим [3]P, [5]P, [7]P, ..., [2w − 1]P. Обнуляем результат. Обходим биты в скаляре k от старшего к младшему и каждый раз: Удваиваем результат. Читаем следующий бит из k, и, если он ненулевой, делаем следующее: Среди последующих w − 1 битов находим последний единичный и прочитываем все биты до него (включительно). Все прочитанные биты, включая самый первый из 4.2, объединяем в число x. Удваиваем результат столько раз, сколько дополнительных бит скушали в 4.2.1. И прибавляем к результату [x]P.
Возьмём какой-нибудь более амбициозный k = 6717 = 11010001111012, а размер окна выберем w = 3. Алгоритм отработает следующим образом:
удваиваем P, запоминаем [2]P вычисляем: [3]P = [2]P + P [5]P = [2]P + [3]P [7]P = [2]P + [5]P k = 1101000111101 result = O 1. удваиваем результат: result = [2]result = [2]O = O 2. читаем следующий бит: k = 1101000111101 ↳ 1 он единичный, поэтому: 1. среди последующих 2 битов находим последний единичный и считываем до него: k = 1101000111101 ^· 2. x = 11₂ 3. удваиваем результат 1 раз: result = [2]result = O 4. прибавляем [x]P к результату: result = result + [11₂]P = [11₂]P 1. удваиваем результат: result = [2]result = [2][11₂]P = [110₂]P 2. читаем следующий бит: k = 1101000111101 ↳ 0 1. удваиваем результат: result = [2]result = [2][110₂]P = [1100₂]P 2. читаем следующий бит: k = 1101000111101 ↳ 1 3. он единичный, поэтому: 1. среди последующих 2 битов находим последний единичный и считываем до него: k = 1101000111101 ·· 2. x = 1₂ 3. удваиваем результат 0 раз 4. прибавляем [x]P к результату: result = result + [1₂]P = [1101₂]P 1. удваиваем результат: result = [2]result = [2][1101₂]P = [11010₂]P 2. читаем следующий бит: k = 1101000111101 ↳ 0 1. удваиваем результат: result = [2]result = [2][11010₂]P = [110100₂]P 2. читаем следующий бит: k = 1101000111101 ↳ 0 1. удваиваем результат: result = [2]result = [2][110100₂]P = [1101000₂]P 2. читаем следующий бит: k = 1101000111101 ↳ 0 1. удваиваем результат: result = [2]result = [2][1101000₂]P = [11010000₂]P 2. читаем следующий бит: k = 1101000111101 ↳ 1 3. он единичный, поэтому: 1. среди последующих 2 битов находим последний единичный и считываем до него: k = 1101000111101 ~^ 2. x = 111₂ 3. удваиваем результат 2 раза: result = [2]result = [2][11010000₂]P = [110100000₂]P result = [2]result = [2][110100000₂]P = [1101000000₂]P 4. прибавляем [x]P к результату: result = result + [7]P = [1101000111₂]P 1. удваиваем результат: result = [2]result = [2][1101000111₂]P = [11010001110₂]P 2. читаем следующий бит: k = 1101000111101 ↳ 1 3. он единичный, поэтому: 1. среди последующих 2 битов находим последний единичный и считываем до него: k = 1101000111101 ~^ 2. x = 101₂ 3. удваиваем результат 2 раза: result = [2]result = [2][11010001110₂]P = [110100011100₂]P result = [2]result = [2][110100011100₂]P = [1101000111000₂]P 4. прибавляем [x]P к результату: result = result + [5]P = [1101000111101₂]P В данном случае пришлось сделать 14 удвоений (1 в подготовительной фазе и 13 в цикле) и 7 сложений (3 + 4 соответственно). С прошлым алгоритмом double-and-add бы потребовалось 13 удвоений, но 8 сложений. Так как сложения дороже, то алгоритм даёт чистый выигрыш в стоимости операций.
В целом он обычно работает не хуже, чем double-and-add, и особенно хорошо, если в скаляре несколько единиц подряд. Но если единичек мало, то подготовительные вычисления могут насчитать точки, которые я потом нигде не заиспользую, и время на вычисления потрачу зря.
Более того, нам необязательно именно складывать постоянно. Можно и вычитать! Несколько странно видеть, конечно, вычитание при вычислении суммы, но смысл в нём имеется, если в скаляре очень много единиц идут подряд.
Например, как в k = 1023 = 11111111112, где их все десять. Вместо 10 удвоений и 10 сложений (double-and-add) или 11 удвоений и 7 сложений (улучшенный алгоритм) можно поступить так. Посчитать [1024]P — на это нужно 11 удвоений. А затем вычесть P. И выйдет это сильно дешевле.
6.2. Генерация цепи сложений-вычитаний
Вот ровно этим код ниже и занимается:
local function getChain(k) local r = {} for i = 0, 383, 1 do -- ① r[i + 1] = k[i // 30 + 1] >> i % 30 & 0x1 end r[385] = 0 -- ⑪ for i = 1, 384, 1 do -- ② if r[i] == 1 then for b = 1, 5, 1 do -- ③ if i + b > 384 then -- ④ break elseif r[i + b] == 1 then -- ⑤ local bit = r[i + b] << b local factor = r[i] + bit -- ⑥ if factor <= 31 then -- ⑦ r[i] = factor r[i + b] = 0 else factor = r[i] - bit -- ⑧ if factor < -31 then -- ⑨ break end r[i] = factor for j = i + b, 385, 1 do -- ⑩ if r[j] == 0 then r[j] = 1 break end r[j] = 0 end end end end end end return r end В данном случае я выбрал w = 5 (и при подготовке нужно будет вычислить от [3]P до [31]P). Цикл ① разбивает 384-битный скаляр k — табличку из 13 30-битных слов — на отдельные биты в r. Эта таблица r кодирует действия улучшенного второго умножения на скаляр, который мы рассматривали. Работать будет так. Проходимся справа налево по r (начиная с первого ненулевого числа) и каждый раз:
удваиваем результат если значение x положительно, то прибавляем [x]P к результату если значение x отрицательно, то вычитаем [−x]P из результата (минус нужен, чтоб отрицательное число стало положительным)
После ① там каждый x — либо 0, либо 1. Поэтому следующий гигантский цикл их сгруппирует вместе для оптимизации. Работает он так.
② — проходимся по всем битам в r. Делаем это в порядке, обратном тому, как потом будем в double-and-add применять: от младшего к старшему.
③ — когда встречаем единицу, считываем следующие биты, пока они есть (④). Среди них интересуют только единичные (⑤).
⑥ — сначала я пытаюсь подтянуть этот единичный бит в r[i], чтобы мы вместо двух сложений сразу прибавили предпосчитанное значение со скаляром побольше (⑦). Например, если r[i] и r[i + 3] оба были равны единице, r[i + 3] я обнулю, а сложение перенесу в r[i], который теперь будет равен 1 + (1 << 3) = 9: будет один прибавляться [9]P, а не два раза по P.
Для b ≤ 4 условие ⑦ будет всегда выполняться, но когда дойдём до b = 5, минимальный скаляр уже будет 33, который мы не считали. Поэтому попробуем вместо этого вычитание (⑧). Если проблем не возникает (⑨), то дальше надо «скомпенсировать» в r то, что мы вычли.
Например, был у нас скаляр 1011111000012 и всё то же окно w = 5. На первом же этапе мы запишем r[1] = 1 − (1 << 5) = −31. Поэтому вместо сложений будем вычитать [31]P. Но тогда и исходное число должно быть на 31 больше. Поэтому 31 = 111112 надо к скаляру прибавить.
В результате у нас ещё одна единичка будет перенесена в старшие разряды, из-за чего единичные разряды будут обнуляться, пока не дойдём до нуля. Дальше будем уже выстраивать цепочку операций для 1100000000002.
Вот в цикле ⑩ это и происходит. А так как у нас старших ненулевых бит среди 384 может не оказаться, надо предусмотреть резервный 385-й (⑪).
Как-то так мы получаем цепочку операций для эффективного вычисления умножения точки на скаляр. Отмечу, что последовательность операций в цепочке зависит от скаляра, который подаётся на вход. Если он секретный, то все эти оптимизации использовать нежелательно. Но так как мы собираемся только проверять подписи, а там все данные публичны, проблемы в неконстантности не возникает.
6.3. Алгоритм вычисления [u]G + [v]P
Собственно, вот сам алгоритм, с помощью которого будем вычислять [u]G + [v]P:
local g = { -- ① { 0x32760ab7, 0x295178e1, 0x355296c3, 0xbc976f, 0x142a3855, 0x1d078209, 0x39b9859f, 0x0ed8a2e9, 0x2d746e1d, 0x1c7bcc82, 0x1378eb1c, 0x08afa2c1, 0xaa87ca, }, { 0x10ea0e5f, 0x290c75f2, 0x17e819d7, 0x182c7387, 0x30b8c00a, 0x28c44ed7, 0x2147ce9d, 0x076f4a26, 0x1c29f8f4, 0x22fe4a4b, 0x06f5d9e9, 0x12a5898b, 0x3617de, }, fieldOne(), } local gWindow = groupDoScalarMultPrecomputation(g) -- ② for i, p in ipairs(gWindow) do -- ③ groupJacobianToAffine(p, p) end local function groupJacobianDoubleBaseScalarMulAdd(d, p, u, v) local pWindow = groupDoScalarMultPrecomputation(p) -- ④ local gChain = getChain(u) -- ⑤ local pChain = getChain(v) -- ⑤ groupJacobianZero(d) local i = 385 while i > 1 and gChain[i] == 0 and pChain[i] == 0 do -- ⑥ i = i - 1 end for i = i, 1, -1 do -- ⑦ groupJacobianDouble(d, d) if gChain[i] > 0 then groupJacobianMixedAdd(d, d, gWindow[1 + (gChain[i] >> 1)]) -- ⑨ elseif gChain[i] < 0 then -- ⑧ groupJacobianMixedSub(d, d, gWindow[1 + (-gChain[i] >> 1)]) -- ⑨ end if pChain[i] > 0 then groupJacobianAdd(d, d, pWindow[1 + (pChain[i] >> 1)]) elseif pChain[i] < 0 then -- ⑧ groupJacobianSub(d, d, pWindow[1 + (-pChain[i] >> 1)]) end end end ① — вне функции определяю g — координаты генератора.
② — здесь же делаю подготовительные вычисления (считаю [3]G, [5]G, ..., [31]G).
③ — для ещё большей эффективности все эти точки я перевожу в координаты (x, y, 1), для чего нужно обращать элемент Z. К счастью, проводится это лишь единожды, а результаты запоминаются. В принципе, можно сохранить результаты вычислений в файл и просто загружать оттуда, а не вычислять здесь.
④ — внутри функции уже выполняю подготовительные вычисления для точки P, которая заранее мне неизвестна.
⑤ — затем составляю цепочки операций для вычисления [u]G и [v]P в отдельности.
⑥ — нахожу первую ненулевую операцию в какой-либо из цепочек.
⑦ — и от неё начиная, иду по остальным операциям с шагом −1.
⑧ — стоит обратить внимание, что тут 2 ветки: на > 0 и < 0. Если в цепочке записан ноль, ничего делать не нужно.
⑨ — благодаря упорной работе в ③ формулы для сложения (вычитания) точек можно заметно упростить. Эти формулы я реализовал в функциях groupJacobianMixedAdd и groupJacobianMixedSub по аналогии с предыдущими, поэтому останавливаться на них не буду.
Сам алгоритм умножения я уже показывал. Здесь он адаптирован, чтоб работать сразу с двумя точками. Так как [2]([uʹ]G + [vʹ]P) = [2uʹ]G + [2vʹ]P, биты в результате не перемешиваются, и алгоритм будет работать, как и ожидалось. В результате получим именно [u]G + [v]P, который нам и нужен.
Это был последний алгоритм на эллиптических кривых, который нужен для проверки подписей в ECDSA. Поэтому можем, наконец, перейти к его реализации. А начнём с декодирования публичного ключа из байтового представления.
7. Декодирование точек эллиптической кривой
Как писал ранее, публичный ключ в ECDSA — это (x, y)-координаты точки [k]G. О том, как именно их нужно кодировать в байты, написан целый стандарт SEC 1. В общем-то, варианта два:
Несжатое представление: пишем байт \x04, потом координату x и координату y как 48-байтовые big-endian-числа. Либо сжатое представление, в котором координата y упускается в принципе, за исключением самого младшего его бита.
Дело в том, что, как я говорил, если (x, y) решает уравнение, то и (x, −y) тоже, то есть один x точку определяет неоднозначно. Но один из этих игреков будет чётным, а другой — нет. Собственно, именно это и показывает самый младший бит. Поэтому в сжатом представлении точка записывается так: если y чётный, пишем \x02, а если нет, то \x03; после этого дописываем координату x как 48-байтовое big-endian-число.
Если нам дали точку в сжатом представлении, то придётся решать уравнение y2 = x3 − 3x + b (на secp384r1 параметр a = −3). Для этого надо уметь извлекать квадратный корень по модулю — число, квадрат которого даст правую часть. Этот квадратный корень — совсем не та вещь, которая для обычных чисел. Например, у числа 2 по модулю p квадратный корень существует и равен
0x83da5b0a808954ec4da354c3d9a24626426da013408406c410ecf78879b99d711354f9802677568aa327717821c98efb.
В целом, если модуль простой, то существуют алгоритмы нахождения квадратного корня. В случае, если модуль этот ещё и даёт остаток 3 при делении на 4 (выполнено для secp384r1), то алгоритм упрощается до возведения в степень (p + 1) / 4. Правда, нужно потом результат проверить (возвести в квадрат и сравнить с исходным), потому что не у каждого числа имеется корень.
На практике же часто этим не заморачиваются и используют несжатое представление. Из него мы получаем два числа: x и y. Но, вообще говоря, не факт, что они взяты с кривой, а не просто рандомные какие-то. Поэтому нужно обязательно проверить, что выполняется y2 = x3 − 3x + b. Тут уже можно просто y возвести в квадрат вычесть правую часть, никаких квадратных корней не требуется. Если получили ноль, то всё прекрасно.
Само собою, отсутствие этой проверки — тоже неплохой такой источник рофловых дыр в безопасности, как и кривые формулы сложения.
8. Арифметика по модулю n
Самое весёлое в ECDSA — это тот момент, когда, реализовав арифметику в одном поле, ты вчитываешься в стандарт и понимаешь, что нужна ещё и в другом. Напомню кратко содержание одной из прошлых глав.
На кривой secp384r1 есть такой генератор G — точка, умножение которой на разные скаляры даёт все точки кривой. По крайней мере, пока эти скаляры 0 ≤ k < n, потому что [n]G = O, и последующие скаляры начинают цикл снова. И вот по простому модулю n-то как раз арифметику ECDSA и хочет. А до этого у нас везде был p, который тоже простой, но другой.
Вот они для сравнения:
p = 0xfffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffffeffffffff0000000000000000fffffffd = 2³⁸⁴ - 2¹²⁸ - 2⁹⁶ + 2³² - 1 n = 0xffffffffffffffffffffffffffffffffffffffffffffffffc7634d81f4372ddf581a0db248b0a77aecec196accc52973 = 2³⁸⁴ - 0x389cb27e0bc8d220a7e5f24db74f58851313e695333ad68d Так что да, придётся запилить ещё и операции по модулю n. Вот какие именно операции нужны в ECDSA для проверки подписи:
вычитание умножение возведение в квадрат обращение по модулю
Так как модуль n в длину такой же, как и p, — 384-битный, то резать его будем так же: тринадцатью словами по 30 бит. Благодаря этому алгоритмы в большей части остаются те же, надо только поменять, как именно мы будем сокращать результаты по модулю n. Потому что в данном случае уже будет выполняться другое равенство:
2384 = 0x389cb27e0bc8d220a7e5f24db74f58851313e695333ad68d.
Эту колбасу порежем на 30-битные куски:
2384 = 0x389 × 2180 + 0x32c9f82f × 2150 + 0x8d220a7 × 2120 + 0x397c936d × 290 + 0x34f58851 × 260 + 0xc4f9a54 × 230 + 0x333ad68d.
Тогда, например, scalarSub, который производит вычитание, выглядеть станет уже так:
-- "order right-shifted by 9" local orderrs9 = { 0x0a52e600, 0x20cb5666, 0x14ef5d9d, 0x06d92458, 0x1bbeb034, 0x2c0fa1b9, 0x3ff8ec69, 0x3fffffff, 0x3fffffff, 0x3fffffff, 0x3fffffff, 0x3fffffff, 0x1ffffffff, } local function scalarSub(c, a, b) local carry = 0 for i = 1, 12, 1 do local word = a[i] - b[i] + orderrs9[i] + carry c[i] = word & 0x3fffffff carry = word >> 30 carry = carry | -(carry & 1 << 63 >> 30) end local word = a[13] - b[13] + orderrs9[13] + carry c[13] = word & 0xffffff carry = word >> 24 c[1] = c[1] + carry * 0x333ad68d c[2] = c[2] + carry * 0xc4f9a54 c[3] = c[3] + carry * 0x34f58851 c[4] = c[4] + carry * 0x397c936d c[5] = c[5] + carry * 0x8d220a7 c[6] = c[6] + carry * 0x32c9f82f c[7] = c[7] + carry * 0x389 carry = 0 for i = 1, 13, 1 do local word = c[i] + carry c[i] = word & 0x3fffffff carry = word >> 30 end end Тут ещё можно заметить, что последний цикл не делает знакового расширения, потому что мы не производим вычитания при сокращении по модулю.
Эти огромные куски приносят дварф-фортрессовское веселье при реализации умножения, когда мы избавляемся от старших 13 слов. Дело в том, что если мы сразу все 13 старших слов сократим в одном цикле, то промежуточные результаты вычислений могут легко переполниться, из-за чего алгоритм насчитает какую-то чепуху. Поэтому приходится посередине операции останавливаться на то, чтобы перенести старшие биты в следующие слова, после чего снова возвращаться к сокращению. Получается вот так:
local i = 26 repeat d[i - 7] = d[i - 7] + carry * 0xe272 -- ① d[i - 8] = d[i - 8] + carry * 0x327e0bc8 d[i - 9] = d[i - 9] + carry * 0x348829f9 d[i - 10] = d[i - 10] + carry * 0x1f24db74 d[i - 11] = d[i - 11] + carry * 0x3d62144c d[i - 12] = d[i - 12] + carry * 0x13e69533 d[i - 13] = d[i - 13] + carry * 0xeb5a340 i = i - 1 carry = d[i] until i == 19 carry = 0 for i = 20 - 13, 19, 1 do local word = d[i] + carry d[i] = word & 0x3fffffff carry = word >> 30 end d[20] = carry i = 20 repeat d[i - 7] = d[i - 7] + carry * 0xe272 -- ① d[i - 8] = d[i - 8] + carry * 0x327e0bc8 d[i - 9] = d[i - 9] + carry * 0x348829f9 d[i - 10] = d[i - 10] + carry * 0x1f24db74 d[i - 11] = d[i - 11] + carry * 0x3d62144c d[i - 12] = d[i - 12] + carry * 0x13e69533 d[i - 13] = d[i - 13] + carry * 0xeb5a340 i = i - 1 carry = d[i] until i == 13 ① — это n, сдвинутый влево на 6 бит. Мы аналогично поступали и с p.
Почему останавливаться надо именно на девятнадцатом слове? Ну, на самом деле, взято это с потолка (как номер посередине между 13 и 26). Спустя год после того, как я этот код написал, я запилил скрипт для Z3, который вычисляет всё то же. Если я не накосячил в нём, то контрпримера, когда бы промежуточные результаты снова переполнились, он найти не смог. Так что должно работать.
И, разумеется, в функции вычисления квадрата числа тоже нужно таким образом переделать сокращение по модулю. Ctrl+C, Ctrl+V.
Содержание функции scalarInvert, как и в прошлый раз, генерируется тулзой addchain. Цепочка операций, правда, получается совсем другой, потому что и степень поменялась с p − 2 на n − 2. Получается вот так.
local function scalarRepeatedSq(b, a, n) scalarSq(b, a) for i = 2, n, 1 do scalarSq(b, b) end end local function scalarInvert(b, a) local t1, t2, t3, t4, t5, t6, t7, t8, t9, t10, t11 = {}, {}, {}, {}, {}, {}, {}, {}, {}, {}, {} scalarSq(t4, a) scalarMul(t2, a, t4) scalarMul(t1, t4, t2) scalarMul(t3, t4, t1) scalarMul(t5, t4, t3) scalarMul(b, t4, t5) scalarMul(t6, t4, b) scalarMul(t4, t4, t6) scalarSq(t7, t4) scalarMul(t7, a, t7) scalarSq(t9, t7) scalarSq(t9, t9) scalarSq(t10, t9) scalarSq(t8, t10) scalarRepeatedSq(t11, t8, 5) scalarMul(t8, t8, t11) scalarRepeatedSq(t11, t8, 10) scalarMul(t8, t8, t11) scalarRepeatedSq(t11, t8, 4) scalarMul(t10, t10, t11) scalarRepeatedSq(t10, t10, 21) scalarMul(t8, t8, t10) scalarRepeatedSq(t10, t8, 3) scalarMul(t9, t9, t10) scalarRepeatedSq(t9, t9, 47) scalarMul(t8, t8, t9) scalarRepeatedSq(t9, t8, 95) scalarMul(t8, t8, t9) scalarMul(t8, t4, t8) scalarRepeatedSq(t8, t8, 6) scalarMul(t8, t3, t8) scalarRepeatedSq(t8, t8, 3) scalarMul(t8, t2, t8) scalarRepeatedSq(t8, t8, 7) scalarMul(t8, t6, t8) scalarRepeatedSq(t8, t8, 6) scalarMul(t8, t6, t8) scalarSq(t8, t8) scalarMul(t8, a, t8) scalarRepeatedSq(t8, t8, 11) scalarMul(t8, t7, t8) scalarSq(t8, t8) scalarSq(t8, t8) scalarMul(t8, a, t8) scalarRepeatedSq(t8, t8, 8) scalarMul(t8, t6, t8) scalarSq(t8, t8) scalarSq(t8, t8) scalarMul(t8, t2, t8) scalarRepeatedSq(t8, t8, 6) scalarMul(t8, b, t8) scalarRepeatedSq(t8, t8, 4) scalarMul(t8, t3, t8) scalarRepeatedSq(t8, t8, 6) scalarMul(t7, t7, t8) scalarRepeatedSq(t7, t7, 5) scalarMul(t7, b, t7) scalarRepeatedSq(t7, t7, 10) scalarMul(t7, t6, t7) scalarRepeatedSq(t7, t7, 9) scalarMul(t6, t6, t7) scalarRepeatedSq(t6, t6, 4) scalarMul(t6, b, t6) scalarRepeatedSq(t6, t6, 6) scalarMul(t5, t5, t6) scalarRepeatedSq(t5, t5, 3) scalarMul(t5, a, t5) scalarRepeatedSq(t5, t5, 7) scalarMul(t5, b, t5) scalarRepeatedSq(t5, t5, 7) scalarMul(t5, t1, t5) scalarRepeatedSq(t5, t5, 5) scalarMul(t5, t3, t5) scalarRepeatedSq(t5, t5, 5) scalarMul(t4, t4, t5) scalarRepeatedSq(t4, t4, 5) scalarMul(t4, b, t4) scalarRepeatedSq(t4, t4, 4) scalarMul(t4, b, t4) scalarRepeatedSq(t4, t4, 5) scalarMul(t3, t3, t4) scalarRepeatedSq(t3, t3, 3) scalarMul(t3, t2, t3) scalarRepeatedSq(t3, t3, 7) scalarMul(t3, t2, t3) scalarRepeatedSq(t3, t3, 6) scalarMul(t3, b, t3) scalarRepeatedSq(t3, t3, 4) scalarMul(t3, t1, t3) scalarRepeatedSq(t3, t3, 3) scalarMul(t3, t2, t3) scalarRepeatedSq(t3, t3, 4) scalarMul(t3, t2, t3) scalarRepeatedSq(t3, t3, 4) scalarMul(t2, t2, t3) scalarRepeatedSq(t2, t2, 6) scalarMul(t2, t1, t2) scalarRepeatedSq(t2, t2, 5) scalarMul(t1, t1, t2) scalarRepeatedSq(t1, t1, 6) scalarMul(b, b, t1) scalarSq(b, b) scalarMul(b, a, b) scalarRepeatedSq(b, b, 4) scalarMul(b, a, b) end Чтоб проверить, что в лапше этой я не накосячил, написал на питоне скриптик, который её парсит и символически производит те же операции. А потом сравнил степень, в которую функция возводит фактически, с ожидаемой. В общем, она совпадает, так что за счёт 11 временных переменных мы получаем правильный результат.
Вот так быстренько пробежались по операциям в поле скаляров. Теперь-то точно всё, что для сверки подписей в ECDSA требовалось, у нас наконец-то имеется.
9. Проверка подписей ECDSA
Что, вообще, такое подпись и где там эллиптические кривые? Подпись ECDSA — это просто пара чисел, которые названы прозрачно понятными именами в лучших традициях математики — r и s. Сначала посмотрим, как они получаются:
Пусть d — приватный ключ (рандомное число от 1 до n − 1), а Q = [d]G — публичный. Генерируем рандомное число k от 1 до n − 1. Вычисляем R = [k]G и разбираем её на координаты R.x и R.y. Затем высчитываем r = R.x % n — это нужно потому, что R.x был по модулю p, а нужен он по модулю n. Наконец, s = k −1(h + rd) % n, где h — это 384-битный хэш сообщения, вычисленный при помощи SHA-384, распаршенный как big-endian-число, взятое по модулю n. Если r или s оказался нулевым, возвращаемся на шаг 2.
Проверяются подписи потом по такому алгоритму:
Проверяем, что 0 < r < n, 0 < s < n, а также валидность публичного ключа. Считаем хэш SHA-384 от сообщения и 48-байтовый результат декодируем как скаляр h (сокращённый по модулю n). Вычисляем u = hs −1 % n, v = rs −1 % n. Находим R = [u]G + [v]Q. Убеждаемся, что R — это не точка на бесконечности, и расщепляем по координатам R.x и R.y. Наконец, чекаем, что R.x % n = r.
Если какая-либо из проверок не удалась, подпись признаётся невалидной. Проверки на r, s и R нужны, чтобы мы не пытались проверять нулевую подпись: без проверки такая подпись «действительна» для любого сообщения. И, само собою разумеется, есть люди, которые забывают про все три этих проверки.
Как видно, алгоритмы эти — просто какие-то числодробилки. На самом деле даже вдаваться в суть преобразований необязательно, чтобы реализовать. Так как все функции, которые в алгоритме проверки подписей юзаются, мы, разлившись на восемь тысяч слов, написали, остаётся их только собрать воедино:
function lib.ecdsaVerifyDecoded(message, hash, r, s, q) if scalarCanonicalFlag(r) == 0 or scalarCanonicalFlag(s) == 0 then -- ① return nil, "invalid signature" end if scalarIsZero(r) or scalarIsZero(s) then return nil, "invalid signature" end local messageHash = hash():update(message):finish():sub(1, 48) local e = scalarFromBytes(messageHash) -- ③ local sInv = {} scalarInvert(sInv, s) local u, v = {}, {} scalarMul(u, e, sInv) scalarReduceQuick(u, u) scalarMul(v, r, sInv) scalarReduceQuick(v, v) local p = groupJacobianZero() groupJacobianDoubleBaseScalarMulAdd(p, q, u, v) if groupJacobianZeroFlag(p) == 1 then return nil, "invalid signature" end groupJacobianToAffine(p, p) -- ② fieldReduceQuick(p[1], p[1]) -- r minus x local rmx = {} scalarSub(rmx, r, p[1]) if not scalarIsZero(rmx) then return nil, "invalid signature" end return true end ① — функции scalarCanonicalFlag просто проверяют, находится ли число в диапазоне от 0 до n − 1. Как это делать, я показывал ранее вместе с реализацией fieldReduceQuick.
② — не забываем, что точки все у нас представлены координатами (X, Y, Z), а нужны (x, y). Это единственное место, где нужно инвертировать Z.
Тут ещё есть одна функция, отмеченная ③, которая загружает скаляр из 48-байтной строки в табличку. Выглядит она так:
local function scalarFromBytes(s) local a13 = (">I3"):unpack(s, 1) local a12 = (">I4"):unpack(s, 4) >> 2 local a11hi = (">I2"):unpack(s, 7) & (1 << 10) - 1 local a11lo = (">I3"):unpack(s, 9) >> 4 local a10 = (">I5"):unpack(s, 11) >> 6 & 0x3fffffff local a9hi = (">I2"):unpack(s, 15) & (1 << 14) - 1 local a9lo = (">I2"):unpack(s, 17) local a8 = (">I4"):unpack(s, 19) >> 2 local a7hi = (">I3"):unpack(s, 22) & (1 << 18) - 1 local a7lo = (">I2"):unpack(s, 25) >> 4 local a6 = (">I5"):unpack(s, 26) >> 6 & 0x3fffffff local a5hi = (">I3"):unpack(s, 30) & (1 << 22) - 1 local a5lo = (">I1"):unpack(s, 33) local a4 = (">I4"):unpack(s, 34) >> 2 local a3hi = (">I4"):unpack(s, 37) & (1 << 26) - 1 local a3lo = (">I1"):unpack(s, 41) >> 4 local a2 = (">I5"):unpack(s, 41) >> 6 & 0x3fffffff local a1 = (">I4"):unpack(s, 45) & 0x3fffffff return { a1, a2, a3hi << 4 | a3lo, a4, a5hi << 8 | a5lo, a6, a7hi << 12 | a7lo, a8, a9hi << 16 | a9lo, a10, a11hi << 20 | a11lo, a12, a13, } end Тут не сказать чтобы дикая магия творится. Просто раскрыл блокнот, нарисовал там 48-байтный прямоугольник, разбил на куски по 30 бит, посчитал смещения и выразил свои чувства в коде.
Кстати, о чувствах.
10. Декодирование подписей
Подпись — это пара 384-битных чисел r и s. Как представить их в байтах? Если у вас всё в порядке с головой, вы, наверное, просто бы закодировали r (в little-endian или big-endian) как 48-байтовую строку, затем приклеили бы так же кодированный s и получили вполне человеческую запись на 96 байт. Кодировать в одну сторону тривиально, в другую ошибиться при декодировании практически негде — мир да согласие.
Но, конечно же, дыр в реализациях ECDSA разработчикам алгоритма показалось слишком мало, и вместо этого подписи кодируются при помощи такого изысканного деликатеса, как ASN.1. Это... ну, я бы назвал это JSON-ом или XML, какими бы они были, если бы их спроектировал большой неповоротливый комитет в 1984 году.
Правда, ASN.1 расшифровывается как abstract syntax notation (one), то есть нотация в абстрактном синтаксисе, что, скажем так, само себе не самое оригинальное имя. Означает оно то, что как именно кодировать (и декодировать) значения, ASN.1 не определяет. Он только даёт язык описания схемы (как XSD), в которой наша подпись выглядит вот так:
ECDSA-Sig-Value ::= SEQUENCE { r INTEGER, s INTEGER } Все эти ключевые слова капсом и пучки укропа (::=) вместо человеческих равенств — настоявшиеся пряности, волнующие язык.
Кодировать эту субстанцию можно в куче разных форматов, включая даже XML и JSON. В мобильной связи, слышал, ASN.1, как майонез, заливают во все дыры и протоколы, и JSON-кодирование ASN.1-схем разработчиков софта радует каждый раз, когда им приходится с этим сталкиваться. Меня пока жизнь уберегла от таких психических травм, к великому счастью.
В криптографии применяется формат DER. В двух местах: в сертификатах и вот этих самых подписях. Так как рисовал формат DER тот самый неповоротливый комитет, который предусмотрел футганы на все случаи жизни, читать стандарт и реализовывать его — утончённое удовольствие для истинных гурманов.
Так что знакомьтесь — я гурман. DER-декодер реализовывал на Lua дважды. Для сертификатов в первую очередь, но здесь вот, в ECDSA, тоже пригодилось. В принципе, ничего не мешает сделать функцию, которая декодирует чисто эту структуру в формате DER, но раз уж у меня был декодер готовый, то его и заюзал.
Были и другие гурманы. На багах в декодере DER проверяльщики подписей тоже сыпались...
11. Тестирование
Как понять, что 2 тысячи строчек кода, которые я написал, делают ровно то, что написано, если код почти весь состоит из неочевидных числопревращений? Тестировать такие низкоуровневые вещи всегда дико весело. К счастью, нам сильно повезло с предметной областью.
Если бы такие же числодробили считали какую-то математическую (физическую) модель, как это отлаживать, я не сильно представляю. В криптографии же это достаточно просто. Во-первых, если где-то накосячить, то почти наверняка результат не сойдётся с ожидаемым. Во-вторых, все вычисления детерминистичны. В-третьих, алгоритмы стандартизованы, и имеется куча имплементаций, с которыми можно сравниваться.
Прежде всего я запилил реализацию операций на эллиптической кривой на Python. В отличие от Lua, в языке Python числа могут быть произвольной длины, и напрямую работать с 384 битами там не представляет никаких трудностей. С помощью неё нагенерил несколько тестовых примеров, с которыми мог сравнивать код на Lua.
Также у корпорации Google на GitHub завалялся один такой большой репозиторий. Там куча JSON-файлов с тестовыми примерами для самых разных криптографических алгоритмов. Достаточно просто распарсить JSON и заюзать данные в своих тестах.
Кроме того, авторы криптографических алгоритмов часто описание самого алгоритма снабжают тестовыми векторами — примерами получаемых результатов при конкретных входных данных. Это те же тестовые данные, по сути. ECDSA не исключение, и тестовые векторы из стандарта тоже проверяю в своих тестах.
Для ряда операций применял более изощрённые методы. Например, когда показывал scalarInvert, я упоминал, что корректность цепочки проверил, её проинтерпретировав символически скриптом на Python. Слова-то сложные, но по сути я просто спарсил регексами кусок кода, который обращал скаляр, и проследил за изменениями экспонент на каждой производимой операции.
Наконец, термоядерный вариант — формальная верификация программ. Формальные методы грубо можно поделить на две группы, которые я назову так: автоматические и ручные. В автоматических формулируется, что именно доказывается, а система потом определяет, верно ли доказываемое утверждение или нет (ищется контрпример). В ручных системах доказательство нужно писать самому, а система лишь проверяет, что доказательство не нарушает логических правил.
Опыт использования последних у меня уже был. Как-то, например, корректность программы на 12 строчек доказывал — определение условия корректности и его доказательство на языке Coq (сейчас он как-то иначе называется) потребовало почти в 100 раз больше строк кода, и он не то чтобы вышел сильно читабельным. Отчасти потому, что в этих вещах я несколько нуб всё ещё, но и в целом что-то объёмное так доказывать — титанический труд.
А вот с автоматическими системами именно для доказательства корректности своих программ сталкиваться не приходилось. Поэтому для проверки, что в scalarMul никакие промежуточные результаты не переполняют 64-битные числа, попробовал реализовать скрипт на Python для SMT-решателя Z3. Условия корректности я там, правда, в полуручном режиме сгенерировал, но всё равно, когда спустя 2.5 часа Z3 написал, что контрпример не существует, было такое стойкое ощущение магии. Крутая вещь, короче говоря.
Так что могу сказать, что код в статье, или по крайней мере в libtls13, делает примерно то, чего от него и требуется.
12. Завершение сеанса
Пакет libtls13 можно скачать с oppm. Вот ссылка на пакет. А вот на сам файл с реализацией ECDSA.
Статья изначально была написана в поджатом изложении на английском языке год назад (но не опубликована), и стена текста выше — это её расширенный перевод.
Аддендум 1. Зачем всё это
Если вы не разделяете моих интересов в математическую сторону вопроса и ориентированы на результат («хочу запилить ECDSA»), большое количество математики может вызвать вопросы, неужели она там действительна нужна в таком объёме. Поэтому решил добавить к статье ещё перечень глав с описанием, зачем я их написал.
§1: эллиптические кривые. Глава определяет уравнение эллиптической кривой и точки на них. Представление точки как пары координат нужно, чтобы было понятно, как устроен приватный ключ (§6), и для реализации операций на кривой (§5). Уравнение нужно, чтобы проверить (§7), что публичный ключ валидный, а не просто рандомная пара чисел. §2: кольца и поля. Глава знакомит с арифметикой по модулю, операциями сложения, умножения и деления, а также объясняет, почему модуль берётся именно простой. Понятие кольца поясняет, почему я могу без опаски манипулировать уравнениями при сложении точек (§5.2). Понятие поля позволяет переиспользовать формулы для вычисления суммы точек кривой (§3), используется для обращения элементов по модулю (§4.6), а отсутствие делителей нуля в нём — критичный факт, на который опирается доказательство корректности алгоритма сложения точек (§5.2). §3: группа эллиптической кривой. Глава определяет две базовые операции на эллиптической кривой: сложение точек и взятие противоположного элемента. Также показывает, почему сложение устроено так необычно. Наконец, вводит нейтральный элемент (нуль) — точку на бесконечности. То, что нуль — это бесконечность, применяется при реализации кривой (§5.1). Алгоритм проверки подписей требует проверки точек (публичного ключа и результата вычислений) на нуль (§9). Сложение точек в алгоритме ECDSA повторяется многократно (§9) — это основная вещь, из-за которой подпись нельзя подделать. Сложность формулы сложения обосновывает выбор альтернативного представления точек в реализации (§5). Взятие противоположной точки используется для реализации вычитания (§5.4). Элементы теории групп (порядок, цикличность) вводятся при описании, почему приватный ключ и публичный ключ именно такие (§6). §4: арифметика в полях. Глава описывает представление элементов поля в Lua, показывает модуль p и параметры кривой a, b. Вид модуля используется во всех операциях для сокращения результатов по модулю (§4.1–§4.6). Значение параметров нужно для проверки публичного ключа (§7). Представление очевидным образом юзается в подглавах (§4.1–§4.6), а также переиспользуется при реализации другого поля (§8). Вся глава нужна для реализации операций на эллиптической кривой (§5). §4.1: нейтральные элементы. Подглава показывает устройство нуля и единицы. Также приводится принцип возврата результата вычислений. Очевидным образом юзается в последующих подглавах (§4.2–§4.6). Нуль и единица в поле нужны для определения нуля на кривой (§5.1). §4.2: сложение. Подглава показывает, как складывать большие числа по большому модулю p. Описывает технику сложения с переносом, сокращение результата по модулю и реализацию знакового расширения на Lua. Сложение в поле применяется при сложении точек кривой (§5.2) и удвоения (§5.3). С небольшими изменениями тот же алгоритм используется для умножения на 30-битное слово (§5.2, §5.3) и вычитании (§4.3). Сокращение по модулю — финальная часть других операций в поле (§4.3–§4.5). §4.3: вычитание. Подглава описывает алгоритм вычитания по модулю, который даёт неотрицательный результат. Как и сложение, используется в формулах для сложения (§5.2) и удвоения (§5.3) точек. Также нужно для вычитания точек, чтобы инвертировать знак у координаты y (§5.4). Аналогичный алгоритм применяется для вычисления по модулю n (§8). §4.4: умножение. Подглава рассказывает об альтернативном представлении умножения «в столбик», вводит операцию свёртки и показывает, как она совмещается с модульной редукцией для вычисления умножения двух чисел по модулю p. Используется при обращении чисел по модулю (§4.6). Этот же алгоритм специализируется для случая a = b, чтобы получить алгоритм возведения в квадрат по модулю (§4.5). Применяется в формулах сложения (§5.2) и удвоения (§5.3) точек. Сложность умножения обосновывает выбор формул. Усложнение для модуля n описывается в §8. §4.5: возведение в квадрат. Подглава описывает видоизменённый алгоритм умножения, чтобы вычислять квадрат числа по модулю. Нужен, чтобы остальные операции (§4.6, §5.2, §5.3) работали заметно быстрее. Применяется в формулах сложения (§5.2) и удвоения (§5.3) точек. Аналогичный алгоритм применяется для вычисления по модулю n (§8). §4.6: обращение чисел. Подглава показывает, как с помощью малой теоремы Ферма обратить элемент по модулю p. Впервые упоминается инструмент addchain, с помощью которого аналогичная операция проделывается по модулю n (§8). Обращение нужно для деления чисел по модулю (§5). Непомерная сложность алгоритма обращения обосновывает выбор альтернативного представления точек на кривой (§5). В этом представлении обращение — финальная операция при вычислении [u]G + [v]Q, которое встречается в алгоритме проверки подписи (§9). §5: реализация операций на эллиптической кривой. Глава вводит более эффективное представление точек на кривой и описывает, как перейти к нему из обычного представления в виде пары (x, y) и обратно. Представление нужно, чтобы операции сложения (§5.2) и удвоения (§5.3) работали гораздо быстрее. Также это представление позволяет выразить точку на бесконечности (§5.1). Описание перехода в обычные координаты используется для ускорения (§6.3) алгоритма вычисления [u]G + [v]Q, применяющегося в ECDSA (§9). §5.1: нейтральный элемент. Подглава описывает устройство точки на бесконечности, а также показывает, как точка хранится в памяти программы. Сравнение с нейтральным элементом — это одна из проверок, которые необходимо выполнить в ECDSA для результата вычислений (§9), и защищает от подделки подписи. Вид нейтрального элемента (нулевая координата Z) используется для составления полного алгоритма сложения точек и доказательства его корректности (§5.2). Представление очевидным образом используется во всём последующем коде. §5.2: сложение точек. Подглава описывает ключевую операцию в эллиптической криптографии — сложение точек, показывает необходимость отдельной формулы удвоения. Сложение в ECDSA повторяется многократно (§9), чему посвящена вся глава §6. Неправильный результат, выдаваемый формулой для сложения точки с самой собой, требует реализации отдельной формулы удвоения (§5.3). В главе также определяется функция fieldReduceQuick, которая применяется в ECDSA для полного сокращения координаты точки-результата по модулю p (§9). На её основе же сделана функция scalarCanonicalFlag, реализующая первую проверку алгоритма ECDSA (§9). §5.3: удвоение точек. Подглава приводит алгоритм вычисления удвоения точки кривой. Необходимо для корректной работы сложения точек (§5.2). Благодаря дешевизне в сравнении со сложением применяется в алгоритме многократного сложения (§6, в частности §6.3). §5.4: вычитание точек. Подглава показывает, как с помощью сложения реализовать вычитание точек кривой. Применяется для ускорения работы многократного сложения (§6.1–§6.3). §6: умножение на скаляр. Приватные ключи. Глава вводит элементы теории групп: понятие порядка, генератора. Вводится обозначение [k]P для многократного сложения (умножения на скаляр) точки кривой. Определяется число n. Описывается устройство приватных и публичных ключей, а также неявно формулируется задача дискретного логарифмирования. Глава показывает, зачем в принципе в криптографии применяются эллиптические кривые. Устройство публичного ключа используется для его декодирования из байтового представления (§7). В главе описывается, почему вычисления в ECDSA производятся именно по модулю n, а не p. Само значение n определяется при реализации поля вычетов по этому модулю (§8). Умножение на скаляр — единственная операция на эллиптической кривой, явно используемая в ECDSA (§9). §6.1: алгоритмы умножения на скаляр. Подглава описывает алгоритм double-and-add и его модификацию. Без этих алгоритмов время работы ECDSA бы выражалось экзотическими числами с сотней знаков в десятичной записи. На модификации алгоритма основывается фактически реализованный алгоритм для вычисления [u]G + [v]Q (§6.3). Принцип работы алгоритмов необходимо знать, чтобы понимать, что вообще происходит при генерации цепи сложений-вычитаний (§6.2). §6.2: генерация цепи сложений-вычитаний. Подглава приводит алгоритм генерации цепи операций для умножения точки на большой скаляр. Используется в конечном алгоритме для вычисления [u]G + [v]Q (§6.3). Позволяет ещё сильнее ускорить этот алгоритм. §6.3: алгоритм вычисления [u]G + [v]Q. Подглава показывает, каким образом выполняется операция на эллиптической кривой, используемая в ECDSA для проверки подписи. Используется в ECDSA для проверки подписи (§9). §7: декодирование точек эллиптической кривой. Глава описывает два формата представления точки в байтах и рассказывает о взятии квадратного корня по модулю. Точка — это публичный ключ. Собственно, программе обычно он поступает в байтовом массиве, поэтому необходимо сначала декодировать, чтобы что-то можно было с ним делать, например проверять подписи. Также в главе описываются важнейшие проверки, которые следует сделать имплементации для публичного ключа, используемого в ECDSA (§9). Проверка публичных ключей обязательна в алгоритме обмена ключами ECDH. Отсутствие может привести к компрометации секретных данных (invalid curve attack). §8: арифметика по модулю n. Глава описывает модификацию алгоритмов из §4 для вычисления по модулю n. В частности, описывается, как производить сокращение по модулю и обращение элементов. Вычитание, умножение и обращение по модулю используются непосредственно в ECDSA (§9). Обращение по модулю использует функцию возведения в квадрат для ускорения работы. §9: проверка подписей ECDSA. Глава приводит реализацию алгоритма проверки подписи ECDSA на Lua. Также она содержит функцию преобразования 48-байтового массива в представление скаляра в виде 13 30-битных слов. Собственно, ради этого вся статья и писалась. ECDSA используется, чтобы убедиться в подлинности сообщения или документа. В частности, алгоритм цифровой подписи — одна из составляющих TLS. Цифровая подпись на основе ECDSA применяется в правительстве США. Схожий с ECDSA алгоритм используется и в России. §10: декодирование подписей. Глава вводит форматы ASN.1 и DER, с помощью которых кодируется подпись. Чтобы подпись, поступившую как набор байт, вообще можно было проверить, её надо декодировать из этих арканных форматов данных. §11: тестирование. Глава приводит список методов тестирования, с помощью которых обеспечивается корректность реализации криптографических алгоритмов. Безопасность основывается на корректности этих алгоритмов. Если в них бага, то либо сервис будет отвергать корректные данные, либо, что ещё хуже, принимать некорректные. Последнее может привести к компрометации всего сервиса. §12: завершение сеанса. Глава заключительная, приводит ссылки на полную реализацию. Небольшая часть функций, в ней используемых, в статье не описаны или описаны лишь вкратце. По ссылке при желании можно ознакомиться с полным кодом программы, скриптами, тестами.
- Подробнее...
- 0 комментариев
- 1 166 просмотров