Программные продукты и системы / Software & Systems
№ 1 (109), 2015
УДК 511.9+519.612 Дата подачи статьи: 06.10.14
DOI: 10.15827/0236-235X.109.055-062
ПАРАЛЛЕЛЬНЫЕ АЛГОРИТМЫ ВЫЧИСЛЕНИЯ ЛОКАЛЬНЫХ МИНИМУМОВ ЦЕЛОЧИСЛЕННЫХ РЕШЕТОК
О.В. Кузьмин, д.ф.-м.н., профессор, [email protected] (Институт математики, экономики и информатики Иркутского государственного университета, бульв. Гагарина, 20, г. Иркутск, 664003, Россия);
В.С. Усатюк, программист, L@ Lcrypto.com (Братский государственный университет, ул. Макаренко, 40, г. Братск, 665709, Россия)
Сделан обзор методов решения задач поиска кратчайшего базиса (Shortest basis problem) и кратчайшего вектора в решетке (Shortest Vector Problem). Для обобщенного метода приведения базиса решеток - блочного метода Коркина-Золотарева (Block Korkin-Zolotarev) - продемонстрирован метод декомпозиции алгоритмов ортогонализации базиса и поиска кратчайшего вектора в целочисленной решетке. Представлены доказанные ранее оценки точности решения задач поиска кратчайшего вектора и кратчайшего базиса целочисленной решетки в зависимости от метода приведения базиса. Получены экспериментальные оценки точности представления чисел в алгоритмах ортогонализации и поиска кратчайшего вектора на ансамблях случайных решеток и на решетках, сложных по Гольштейну-Майеру. Продемонстрировано падение эффективности слабосвязанных параллельных вычислительных устройств с ростом размерности решетки, обусловленное необходимостью роста точности представления чисел. Полученный параллельный алгоритм позволяет осуществить линейное ускорение алгоритма поиска локальных минимумов в решетке в зависимости от числа вычислительных устройств, допускающих выполнение сколь угодно точной арифметики. Проведено экспериментальное сравнение параллельных реализаций алгоритма поиска кратчайшего вектора на видеокартах CUDAENUM и многоядерных мультипроцессорных конфигурациях parENUM. Представленный алгоритм продемонстрировал сходный результат на решетках малой размерности и существенно лучший в случаях высокой размерности в сочетании с использованием экстремального отсечения ветвей в алгоритме Каннана. Реализация представленного алгоритма заняла первое место среди методов приближенного поиска кратчайшего вектора в случайных ансамблях решеток, построенных на сложных по Голдштейну-Майеру идеалах кольца кругового многочлена.
Ключевые слова: приведение базиса решетки, целочисленные решетки, поиск кратчайшего базиса, поиск кратчайшего вектора, блочный метод Коркина-Золотарева, метод Ленстра-Ленстра-Ловаса.
Математический аппарат теории решеток тесно связан с абстрактной и геометрической алгебрами, числовыми методами приближенного анализа, теорией информации, оптимизацией и криптоанализом, что позволяет применять методы теории решеток для решения ряда задач: линейного программирования, плотной упаковки тел (в частности, сфер в евклидовом пространстве) [1], факторизации многочленов с рациональными коэффициентами [2], поиска диофантовых приближений [3], оптимизации кодирования с использованием адаптивных антенных решеток со слабо коррелированными антенными элементами (MIMO) [4] и многих других [5-9].
Основные понятия и обозначения
Определение 1. Решеткой называется дискретная абелева подгруппа, заданная в пространстве
Rn.
Пусть базис B={bb ..., bn} задан линейно независимыми векторами в Rm. Тогда под решеткой будем понимать множество целочисленных линейных комбинаций этих векторов:
L(b1, •••, bn ) = {е*А : (Х >•••> X ) е ” }>
где m и n - размерность и ранг решетки соответственно, m>n.
Определение 2. Решетки, у которых размерность и ранг равны, называются полными.
Определение 3. Решетки Ь\ и L2, заданные базисами B = {bl,..., bn } , B' = {Ь\,..., Ь'п } , конгру-
энтны; L\(B) = L2(B'), если объемы фундаментальных параллелепипедов, образованных их базиса-
равны: det(Li)= det(L2), где det L =
)det(BT B)
ми
или detL = | detB | для полноразмерных решеток.
Множество конгруэнтных решеток может быть получено в результате умножения базиса B решетки на целочисленные унимодулярные матрицы.
Определение 4. Пусть L - решетка ранга n. i-м соответствующим минимумом, где i=1, 2, ..., n, называется точная нижняя граница радиусов замкнутых шаров с центром в нуле, содержащих i линейно независимых векторов решетки L:
X, (L) = inf{r | dim (span) Lp|B( 0, r))) > i},
где B(0, r) = {x e Rm | x < r} - замкнутый шар ра-
диуса r с центром в нуле.
Длине кратчайшего вектора в решетке L соответствует Xj(L), что приводит непосредственно к
55
Программные продукты и системы / Software & Systems
№ 1 (109), 2015
формулировке задачи поиска кратчайшего и короткого векторов в решетке.
Определение 5. Задача поиска короткого вектора (Д-short vector problem, SBPY(m)): пусть дана m-мерная решетка L(B) ранга n и вещественное Д>1. Найти нетривиальный вектор, в Д раз больший кратчайшего вектора в решетке
Ъ е L :|\b\\ <А-Х1 (L).
При у=1 решается задача поиска кратчайшего вектора в решетке, при у>1 - короткого вектора.
Определение 6. Константой Эрмита называется величина ym, заданная выражением
Y,
m SUP
L
К (L)2
| L - полная
(det L)m
Применяя первую теорему Минковского о соответствующих минимумах, с одной стороны, и лемму о норме кратчайшего вектора в решетке, получим соответственно верхнюю и нижнюю оценки длин кратчайшего вектора в решетке:
min
Ml (L) det(L)»
где Ь - ортогональный вектор, полученный из базиса одним из QR-методов разложения матрицы базиса решетки: Грама-Шмидта, поворота Гивенса или преобразования Хаусхолдера [10].
Определение 7. Задача поиска короткого базиса решетки (Д-short basis problem, SBPY(m)): пусть дан базис полной решетки B и вещественное Д>1. Найти базис
В’ ={b1, Ьbm}: L(B) = L(B'),
ш|| ь:\\ <y nrj bill •
При у=1 полученный базис достигает нижней границы неравенства Адамара и является ортого-
/«nil т
нальным: П bn = det (L) < П||b,|| .
i=i “р i=i р
Определение 8. Базис B={b1, ..., bm} решетки L<zRm приведен по длине, если в результате орто-гонализации решетки методом Грамма-Шмидта выполняется следующее неравенство:
1 2'
где р.у - коэффициенты Грамма-Шмидта.
Определение 9. Упорядоченный по длине базис B={b1, b2, ..., bm} решетки L<zRm приведен блочным методом Коркина-Золотарева с блоком
h, А ’ 1 ^ J <1 < m
Ре[2, m] и точностью
Sell;1
если базис B
приведен по длине или S2 • |bx| <Xj2(Li), i=1, ...,
m, где Xi(L,) - длина кратчайшего вектора в решетке Li, образованной ортогональным дополнением векторного пространства с базисом bi, ...,
bmin(i+p-l, m).
Определение 10. Базис решетки приведен по Ленстра-Ленстра-Ловасу (LLL-алгоритмом), если он приведен блочным методом Коркина-Золотарева с размером блока р=2 и точностью ортого-1
нализации S e l —
12
Определение 11. Локальным минимумом решетки называем вектор b={ bl , b2, ..., bn}, beL, для которого справедливы неравенства: | ci |<| bi |, V b=(c1, c2, ..., cm)eL, i=1, ..., n, при этом хотя бы одно из них строгое.
В теории оптимальных сеток Коробова такие векторы образуют узлы оптимальных сеток. Нетрудно убедиться, что эти векторы являются кратчайшими в решетке, а подмножество локальных минимумов образует кратчайший базис решетки. Число всевозможных конфигураций, которые образуют кратчайшие базисы для m-мерной решетки L объема detL, не превосходит величину
C(n)log”-1 det L
SBP < -
C(n) - некоторая поло-
жительная константа [11].
Таким образом, задача поиска локальных минимумов содержит в себе задачу поиска кратчайшего базиса, известную как алгоритм приведения базиса решеток. На настоящий момент неизвестны полиномиальные алгоритмы приведения базиса по Минковскому, за исключением решеток размерности 2, предложенных Гауссом [12], и 3, предложенных Валле (с кубической сложностью) [13] и улучшенных Семаевым, снизившим сложность до квадратичной [14]. В работе [15] описан алгоритм для приведения базиса по Минковскому до 6 измерений. В работе [16] результат был улучшен до 7-мерных решеток. По причине высокой вычислительной сложности практический интерес имеют алгоритмы поиска локальных минимумов с полиномиальной, псевдо- и субэкспоненциальной сложностью. Для приведения базиса решетки по Ловасу применяется полиномиальный LLL-алгоритм, предложенный в [17]. Для приведения по Коркину-Золотареву применяется блочный метод SE BKZ, предложенный Шнором-Охнером [18]. При изменении размера рассматриваемой подрешетки р и параметра точности ортогональности базиса 5 получаем иерархию методов от приведения по Ленстра-Ленстра-Ловасу до блочного метода Коркина-Золотарева (см. табл.).
Построим, варьируя размер приводимого базиса р и релаксируя норму в блочном методе Коркина-Золотарева, параллельный алгоритм поиска локального минимума целочисленной решетки.
На основании нижней оценки из определения 6 нетрудно убедиться, что поиск локального минимума условно состоит из двух алгоритмов ортого-нализации базиса и поиска кратчайшего вектора в некоторой решетке [19]. Осуществим декомпозицию и распараллеливание этих алгоритмов.
2
n
т
г=1
56
Программные продукты и системы / Software & Systems
№ 1 (109), 2015
Сложность и точность методов приведения решеток в евклидовом пространстве
Complexity and accuracy of lattice reduction methods in Euclidean space
Тип приведения базиса решетки Точность приведения Алгоритмы и их временная сложность
По Минковскому SVP, = 1, т SBP, < (Ym ) ^ Лагранж, 1773 Гаусс (m=2), 1850 Валле (m=3, 0(m3)), 1986 Семаев (m=3, 0(m2)), 2001 Блумер (O(m!x 3m)), 2001
| m Он m 0 щ 1 св Ё 8 1 & ft ь m § О ГО По Коркину-Золотареву, P=m m-1 SVPA <чГ , Коркин-Золотарев, 1872 Лагарис-Ленстра-Шнор, 1990
2<р<от /у" х ЛМ м( М-1) SBPa W6 х (м + 3)! хув^ л 6 х 2М в Шнор-Охнер, 1991
S-приведенный по Ловасу (Lovasz) (Р=2) m-1 SVPA = 2~г , m(m-1) SBPA < а ^ Ленстра-Ленстра-Ловас, 1982
Параллельный алгоритм ортогонализации решеток
Классический метод Грама-Шмидта, применяемый для ортогонализации базиса решетки, требует порядка 2mn2 или 2n3 (при m=n) операций с плавающей запятой. Однако классический метод Грама-Шмидта численно неустойчив вследствие ошибок округления [20]. По этой причине на практике в последовательной реализации алгоритмов применяют модифицированный метод Грама-Шмидта [21]. Последний, как и классический метод Грама-Шмидта, обладает сходной вычислительной сложностью. Однако будучи рекурсивным по своей сути, алгоритм Грама-Шмидта обладает неглубокой декомпозицией данных на уровне компонентов вектора при наличии множества блокировок. Эта особенность приводит к неэффективности параллельных реализаций, которая была продемонстрирована в работе [22], производительность в лучшем из случаев не превосходила 60 % от пиковой производительности на матрицах размера 6 000x6 000 для SCI-кластера на 16 XEON 2 ГГц и 37 % для Beowulf-кластера (CLiC) с 528 Pentium III 800 МГц. Схожий результат был получен для параллельных алгоритмов ортогонализа-ции Грама-Шмидта, реализованных на видеокартах NVIDIA CUDA 295 (1 800 Гфлопс при однократной и 150 Гфлопс при двукратной точности вычисления) на матрицах размера 4 000x4 000 и двойной точности вычислений, 15-кратный прирост по сравнению с одноядерной реализацией [23].
Перейдем от традиционно используемого метода ортогонализации Грама-Шмидта к системе эквивалентных определений на основе QR-методов ортогонализации.
Определение 12. QR-разложением базиса решетки называют представление матрицы B = (bj, b2, ..., b) размера mxn в виде произведения
Q-унитарной и R-верхней треугольной матриц:
= (Ь,b2, ..., bn) = QT х R =
( r 1,1 Г1,2 ‘ ‘ ‘ r ^ r1,n
ObЧ,•••, Чп) 0 Г2,2 * Г2,п
v 0 * 0 r n,n
В свою очередь, верхнюю треугольную матрицу можно представить в виде
f r i,i ri,2 - r ^ r1,n
R = 0 Г2,2 — r2,n
v 0 0 Гп,п y
"r1,1 0 — 0 3 f 1 ^2,1 — ^n-1,1 ^n,1
0 r2,2 0 1 ^3,2 — ^n,2
0
0 — 0 r 0 01 ^n,n-1
v n,n у v0 0 —0 1 J
Таким образом, ортогонализация Грама-Шмидта связана с QR-разложением следующим
г. .
образом: ц. . = —, a* = г.. • q., где q- - i-й вектор’. r . ’
i
столбец ортогональной матрицы Q.
В таком случае традиционные определения методов приведения могут быть переформулированы в форму, использующую параллельные QR-методы ортогонализации.
Определение 13. Базис решетки L(B), образо-
ванный векторами B = {Ьг,bn} с Rm, называют приведенным по длине (ослабл. по Эрмиту) [24], если для коэффициентов верхней треугольной матрицы R, QR-разложения базиса решетки B вы-
полняется у, j
i,
< — r.
21 м
,1 < i < j < n .
57
Программные продукты и системы / Software & Systems
№ 1 (109), 2015
Определение 14. Базис решетки L(B) приведен по Ловасу [17], если он приведен по длине и для коэффициентов R QR-разложения базиса решетки B выполняется условие r2 + R-u - §R-и-i ’ 1<г'-п’ Se(0,25, 1].
Определение 15. Базис решетки L(B) приведен по (Эрмиту)-Коркину-Золотареву (HKZ, KZ) [25], если он приведен по длине и для всех квадратных подматриц R' размера n—i+1, 1—i<n, полученных из R вычеркиванием i первых строк и столбцов; первый вектор-столбец является кратчайшим в решетке L'(R'), Tf(i) = Xf (L').
Опираясь на указанные определения, авторы провели сравнение методов ортогонализации базисов решетки с использованием библиотек линейной алгебры для многоядерных процессоров AMD CORE MATH Library, Intel Math Kernel Library и видеоускорителей NVIDIA CUDA CUB-LAST, AMD OPENCL CLMAGMA, Accelerated Parallel Processing Math Libraries (CLAMDBLAS) (см. рис. 1, 2). Кривые отображены на рисунках сверху вниз в порядке их упоминания в легенде. Видеоускоритель GPU NVIDIA 550 TI (691 Гфлопс при однократной и 58 Гфлопс при двукратной точности вычисления, то есть видеокарте с почти в три раза меньшей производительностью) дает 30-кратный прирост в случае использования метода поворота Гивенса на видеокарте и почти 170-кратный прирост при использовании метода отражений Хаусхолдера на четырех ядрах Phenom X4 965 3.4 ГГЦ [26].
Однако следует отметить, что применение полученных методов затруднено следующими обстоятельствами. С ростом размерности решетки рост необходимой точности представления чисел с плавающей запятой невысок для ансамблей случайных решеток [18]. В случае применения базисов, соответствующих плотным рюкзачным задачам [27] или же решеткам, сложным по Гольдштейну-Майеру [28, 29], необходимая точность арифметики растет линейно (см. рис. 3). Последнее обстоятельство существенно затрудняет использование видеоускорителей. Так, падение пиковой производительности арифметики двукратной точности по сравнению с однократной - 1 к 4 для специализированных вычислительных ускорителей и 1 к 24 для массовых видеокарт. С дальнейшим ростом точности вычислений обменное отношение резко возрастает.
С другой стороны, с ростом размерности решетки время, затрачиваемое на ортогонализацию базиса, становится небольшим по сравнению со временем поиска кратчайшего (короткого) вектора в подрешетке. Однако следует отметить, что в случае применения некоторых вероятностных методов, в частности, основанных на отсечении веток при поиске кратчайшего вектора в подрешетке [30, 31], применение параллельных методов орто-
Рис. 1. Зависимость времени выполнения алгоритма QR-разложения от размера базиса решетки при однократной точности вычислений (float, 7 десятичных цифр)
Fig. 1. Dependence of QR-decomposition algorithm execution time from basis reduction size with single-shot calculation accuracy (float, 7 decimal numerals)
Рис. 2. Зависимость времени выполнения алгоритма QR-разложения от размера базиса решетки при двукратной точности вычислений (double, 15 десятичных цифр)
Fig. 2. Dependence of QR-decomposition algorithm execution time from basis reduction size with double calculation accuracy (double, 15 decimal numerals)
гонализации становится актуальным при приведении решеток высоких размерностей.
58
Программные продукты и системы / Software & Systems
№ 1 (109), 2015
Рис. 3. Рост точности вычислений с ростом размерности решетки, а также соотношение времени, затрачиваемого алгоритмами на поиск кратчайшего вектора и ортогонализацию базиса
Fig. 3. Increasing of calculation accuracy with increasing lattice dimension. In addition, relations of the shortest vector algorithm time and basis orthogonalization
Параллельный алгоритм поиска кратчайшего вектора
Решение задачи поиска кратчайшего вектора в решетке заключается в полном переборе всех линейных комбинаций векторов базиса
У m x,.b,.
i=1 1 1
< A2, x,eZ, где A - норма иско-
мого кратчайшего вектора. В качестве нормы берется верхняя оценка длины кратчайшего вектора,
A = TYT det(L)m , где ym - константа Эрмита, в тех случаях, когда наименьший из векторов в базисе
__ 1
решетки превосходит A, ||bj| > л/Ym det(L)m . По этой причине предварительное приведение базиса решетки позволяет уменьшить пространство перебора Х{. С этой целью распишем базис решетки через ортогональные векторы, для простоты изложения используя метод ортогонализации Грама-
Шмидта. Получим bi = ’^j=^ijbj, 2<i<m, 1<j<i<m,
где р-у - коэффициенты Грама-Шмидта. Легко убедиться, что в этом случае поиск кратчайшего вектора сводится к решению системы неравенств:
*21ЮГ * A2
(*2-1 +^2,2-1 Х2 ) Ъ,
||2 ~ . II | ||2
-2-Л ^ A - X2 И ,
(*1 + £ i кj)21|ъг|2* * a2 - £ i
j=2
где lj = (Xj + ^ )2 ||b|| , и выбору одного из
i=j+1
целочисленных векторов, у которого норма ска-
лярного произведения с базисом решетки минимальна.
Формализуя задачу, получим обход дерева от корня к листу, в каждой из вершин которого решается соответствующее линейное уравнение. Из корня этого дерева выходит
2 •
A = 2 •
У
УгГdet(L)
1Ю1
ветвей,
а в случае предварительного приведения базиса
решетки 2 •
Ю
ветвей. В силу симметричности
дерева (по свойствам нормы) для получения искомого кратчайшего вектора необходимо перебрать только половину его вершин. В результате полного обхода от корня к листу получим предполагаемый кратчайший вектор Х с нормой, меньшей либо равной искомой. Если норма полученного вектора будет меньше заданной ранее, целесообразно обновить ее с целью уменьшения пространства перебора. Алгоритм останавливается, когда завершен обход вершин дерева или в случае поиска короткого вектора получен вектор с достаточной нормой. Каждый из потоков вычисляет свою ветку, исходящую из корня дерева.
Далее приведем параллельный алгоритм поиска кратчайшего вектора в решетке L(B), B={bi , ^
— ,bm} :
Вход: IKII >|| Ь2\\ >->|| К\ , Ргр
Выход: Вектор x=(xb
Х2,
Xm)eZm:
£ xibi
= ^ (ЦБ)) .
1. Вычислить xm-k координатных компонент на k уровнях от корневого узла.
2. Каждый из потоков перебирает свое поддерево xm, ..., xm-k, i<m-k.
h = (Xi + £ jj,i)21|b,X
j=i+1
Если £ lJ < A и i=1, то X = X U xm
j=1
xi = xi+1;
если x =
£ xj bj
j=i
< A, то A = xm
Если £ lj < A и i^1, то i=i-1,
u - Z h - Z (^)-Члт-
j=i+1
2
i=1
2
J-i
m
i=2
X -
59
Программные продукты и системы / Software & Systems
№ 1 (109), 2015
т
Если ^ lj > A , то i=i+1, x,=x,+L
j=i+1
3. Возвращаем вектор с минимальной нормой из множества X.
Данный алгоритм был реализован на основе потоковой модели NPTL (Native POSIX Thread Library [32]) под CentOS 6.3. При осуществлении приведения 103-мерной решетки BKZ-методом при Р=52, 4 потоках, исполняемых на AMD Phenom 965/8 Gb DDR2-800, получено трехкратное ускорение выполнения метода по сравнению c fplll-4.0.1 [33]. Предложенный метод позволяет перебирать порядка 10 миллионов координатных компонент в секунду. Было проведено сравнение с реализациями паралелльного алгоритма поиска кратчайшего вектора на видеокартах CUDAENUM [34] и многоядерных мультипроцессорных конфигурациях parENUM [35]. Для сравнения применялась следующая тестовая конфигурация: fplll-4.0.4, parENUM, pLRT: AMD Phenom II X4 965 3.4Ghz/8Gb DDR2-800|12.5 ГБ/c, ОС CentOS 6.4 (2.6.32-358.6.2.el.x86_64). Пиковая производительность SP 108.8; DP 54.4 Гфлопс. Для ParENUM применялись пакеты mpich2-dev x86_64 L3.2p1-1.fc16, mprf-3.1.2, gmp.x86_64 4.3.1-
7.el6_2.2. Для CUDAENUM Intel Core i7-3930K, 64GB DDR3-2133 (16.5 ГБ/c), GTX 680 2Gb SP 3090; DP 129 Гфлопс (192 ГБ/C), ОС CentOS 6.4 (2.6.32-358.14.1.el.x86_64), NVIDIA CUDA SDK 5.5x64. Приводились решетки, сложные по Голь-штейну-Майеру. Чтобы производительность видеокарты не упала вследствие необходимости применения четырехкратной точности и производительность не упиралась в шину памяти, применялись сравнительно малые размеры подрешетки и значительно более мощные конфигурации рабочей станции. Результаты сравнения представлены на рисунке 4.
Рис. 4. Сравнение времени работы алгоритмов параллельного поиска кратчайшего вектора на решетках, сложных по Гольштейну-Майеру, различной размерности
Fig. 4. Comparison of the shortest vector parallel search algorithm execution time on N-dimensional Goldstein-Mayer lattices
С использованием приведенных алгоритмов было реализовано приложение, занявшее 1-е и 3-е места на международном конкурсе алгоритмов
поиска коротких векторов [36], для норм
i
— (м
A = m • det(L)m и A = 1,05 -Г1—+1
соответственно.
Таким образом, в работе предложен параллельный алгоритм решения задачи поиска локальных минимумов в целочисленных решетках на основе блочного метода Коркина-Золотарева. Приведены точностные характеристики метода при решении задач приближенного поиска кратчайшего вектора и базиса в целочисленной решетке. Выполнено сравнение алгоритма и его аналогов. Реализованный на основе предложенного алгоритма программный комплекс установил мировой рекорд на международном конкурсе поиска кратчайших векторов на ансамблях решеток, построенных на кольце кругового многочлена, сложных по Гольштейну-Майеру.
Литература
1. Conway J.H., Sloane N.J.A. Sphere Packings, Lattices and Groups. Springer. 3rd ed., 1999, 703 p.
2. Dwork C. Lattices and their Application to cryptography [Lecture Notes]. Stanford Univ., 1998, 116 p.
3. Grotschel M., Lovasz L., Schrijver A. Geometric Algorithms and Combinatorial Optimization. Springer-Verlag, 1993, 564 p.
4. Mobasher A. Applications of Lattice Codes in Communication Systems. Canadian theses. Univ. of Waterloo. Dept. of Electrical and Computer Engineering, 2008, 147 p.
5. Bernstein D.J., Buchmann J., Dahmen E. (Eds.) Post Quantum Cryptography. Springer, 2009, pp. 147-191.
6. Darte A., Schreiber R., Villard G. Lattice-based memory allocation Computers, IEEE Transactions, 2005, no. 10 (54), pp. 1242-1257.
7. Lagarias J.C., Odlyzko A.M. Solving low-density Subset Sum problems, J. Assoc. Comp. Mach, 1985, no. 1 (32), pp. 229-246.
8. Коробов Н.М. Приближенное вычисление кратных интегралов с помощью методов теории чисел // Докл. Акад. наук СССР. 1957. Т. 115. № 6. С. 1062-1065.
9. Коробов Н.М. Теоретико-числовые методы в приближенном анализе. М.: Изд-во МЦНМО, 2004. 288 с.
10. Press W.H., Teukolsky S.A., Vetterling W.T. Numerical Recipes: The Art of Scientific Computing. NY: Cambridge Univ. Press, 2007, 1262 p.
11. Быковский В.А. О погрешности теоретико-числовых квадратурных формул // Докл. РАН. 2003. Т. 389. № 2. С. 154-155.
12. Gauss C.F. 1801. Disquisitiones Arithmeticae. Leipzig; Translated into English by A. C. Clarke. New Haven, CT: Yale Univ. Press, 1966.
13. Vallee B. Une Approche Geometrique de la Reduction de Reseaux en Petite Dimension. PhD Thesis, Univ. de Caen, 1986.
14. Semaev I. A 3-dimensional lattice reduction algorithm. Lecture Notes in Computer Science, Springer-Verlag. 2001, no. 2146, pp. 181-193.
15. Рышков С.С. К теории приведения положительных квадратичных форм по Эрмиту-Минковскому. Исслед. по теории чисел. 2, Зап. научн. сем. ЛОМИ, 33. Л.: Наука, 1973. C. 37-64.
16. Afflerbach L., Grothe Н. Calculation of Minkowski-reduced lattice bases. Computing, 1985, no. 3-4 (35), pp. 269-276.
17. Lenstra A., Lenstra H., Lovasz L. Factoring polynomials with rational coefficients. Mathematische Annalen, 1982, no. 261 (4), pp. 515-534.
п2 • det(L)m
60
Программные продукты и системы / Software & Systems
№ 1 (109), 2015
18. Schnorr C.P., Euchner M. Lattice basis reduction: Improved practical algorithms and solving subset Sum Problems. Fundamentals of Computation Theory, 1991, pp. 68-85.
19. Кузьмин О.В., Усатюк В.С. Программный комплекс приведения базиса целочисленных решеток // Программные продукты и системы. 2012. № 4 (100). С. 180-183.
20. Longley J.W. Modified Gram-Schmidt process vs. classical Gram-Schmidt. Communicationsin Statistics - Simulation and Computation. 1981, no. 10 (5), pp. 517-527.
21. Higham N.J. Accuracy and Stability of Numerical Algorithms Society for Industrial and Applied Mathematics, 2002, 711 p.
22. Runger G. and Schwind M. Comparison of Different Parallel Modified Gram-Schmidt Algorithms. Euro-Par 2005, LNCS 3648, pp. 826-836.
23. Milde В., Schneider М. Parallel Implementation of Classical Gram-Schmidt Orthogonalization on CUDA Graphics Cards. TU Darmstadt Fachbereich Informatik Kryptographie und Computeralgebra, 2009, URL: http://www.cdc.informatik.tu-darm-stadt.de/~mischnei/ CUDA_GS_2009_12_22.pdf (дата обращения: 20.06.2014).
24. Hermite C. Extraits de lettres de M. Hermite M. Jacobi sur differents objets de la theoriedes nombres.Reine Angew Math. 1850, no. 40, pp. 279-290.
25. Korkine A., Zolotareff G. Sur les formes quadratiques. Maht. Ann., 6, 1873, pp. 366-389.
26. Усатюк В.С. Реализация параллельных алгоритмов ор-тогонализации в задаче поиска кратчайшего базиса целочисленных решеток. ПДМ. Приложение. 2012. № 5. С. 120-122.
27. Schnorr C.P., Shevchenko T. Solving Subset Sum Problems of Density close to 1 by "randomized" BKZ-reduction. IACR Cryptology ePrint Archive, 2012, 620 p.
28. Goldstein D., Mayer A. On the equidistribution of Hecke points. Forum Mathematicum., 2003, vol. 15, no. 2, pp. 165-189.
29. Кузьмин О.В., Усатюк В.С. Комбинаторная оценка мощности множества, содержащего кратчайший вектор, образованного базисом решетки, приведенной блочным методом Коркина-Золотарева. URL: http://arxiv.org/abs/1302.4062 (дата обращения: 20.06.2014).
30. Schnorr C.-P., Horner H.H. Attacking the Chor-Rivest cryptosystem by improved lattice reduction. Springer-Verlag, LNCS, 1995, vol. 921, pp. 1-12.
31. Gama N., Nguyen P., Regev O. Lattice Enumeration using Extreme Pruning. Advances in Cryptology - EUROCRYPT 2010. 29th Conf. on the Theory and Applications of Cryptographic Techniques, 2010, Springer, pp. 257-278.
32. Kerrisk M. The Linux Programming Interface: A Linux and UNIX System Programming Handbook. No Starch Press, 2010, 1552 p.
33. Application fplll. URL: http://perso.ens-lyon.fr/damien. stehle/fplll (дата обращения: 20.06.2014).
34. Application GPIENUM. URL: http://homes.esat.kuleuven. be/~jhermans/gpuenum/ (дата обращения: 20.06.2014).
35. Application parENUM. URL: http://xpujol.net/ (дата обращения: 20.06.2014).
36. Ideal lattice challenge. URL: http://www.latticechallenge.
org/ideallattice-challenge/index.php (дата обращения:
20.06.2014).
DOI: 10.15827/0236-235X.109.055-062 Received 06.10.14
PARALLEL ALGORITHMS FOR INTEGER LATTICES BASIS REDUCTION
Kuzmin O. V., Dr.Sc. (Physics and Mathematics), Professor, quzminov@ mail.ru (Institute of Mathematics, Economics and Informatics of Irkutsk State University,
Gagarin Blvd. 20, Irkutsk, 664003, Russian Federation);
Usatyuk V.S., Programmer, L@ Lcrypto.com (BratskState University, Makarenko St. 40, Bratsk, 665709, Russian Federation)
Abstract. The article presents a multithread heterogeneous (GPU+CPU) implementation of block Korkin-Zolotarev lattice reduction methods. The paper considers connection and hierarchy of the main three types of lattice basis reduction: size reduced (weak Hermit), Lovasz condition (Lenstra-Lenstr-Lovasz), Korkin-Zolotarev as known as Hermit-Korkin-Zolotarev (BKZ-methods). The authors provide references to applications in: information theory (decoding of coding group in MIMO), calculus (minimize of the positive quadratic form), complexity theory and cryptanalysis of Merkle-Hellman cryptography (solving subset sum problems), algebra and control theory (linear diophantine equation solving system), compiler theory (lattice based memory allocation), Vinogradov-Korobov number-theoretic method in the calculus, methods synthesize cryptographic and cryptanalysis in lattice based cryptography. The authors have estimated experimental bounds of the necessary number precision to solve the shortest basis problem as well as bounds of the shortest basis problem and shortest vector problem solution accuracy under block Korkin-Zolotarev reduction. They decompose Kannan algorithms with extreme pruning techniques for solving the shortest vector problem and QR-decomposition of lattice basis. Implementation is based on POSIX Tread Library, CLAMDBLAS, CUBLAS, CLMAGMA, ACML, IMKL and shows linear speedup with a number of processor admitting variable-precision arithmetic. There is a comparison with several multithread implementations: fpLLL, CUDAEnum and parENUM. The algorithm implementation won the first place in the ideal lattice challenge.
Keywords: lattice basis reduction, integer lattices, shortest vector problem, SVP, shortest basis problem, SBP, Block Korkin-Zolotarev reduction, BKZ, LLL reduction.
References
1. Conway J.H., Sloane N.J.A. Sphere Packings, Lattices and Groups. Springer Publ., 3rd ed., 1999, 703 p.
2. Dwork C. Lattices and Their Application to Cryptography. Lecture Notes, Stanford Univ. Publ., 1998, 116 p.
3. Grotschel M., Lovasz L., Schrijver A. Geometric Algorithms and Combinatorial Optimization. Springer-Verlag, 1993, 564 p.
4. Mobasher A. Applications of Lattice Codes in Communication Systems. Canadian Theses. Univ. of Waterloo Publ., 2008, 147 p.
5. Bernstein D.J., Buchmann J., Dahmen E. (Eds.) Post Quantum Cryptography. Springer Publ., 2009, pp. 147-191.
6. Darte A., Schreiber R., Villard G. Lattice-Based Memory Allocation Computers, IEEE Transactions. 2005, no. 10 (54), pp. 1242-1257.
61
Программные продукты и системы / Software & Systems
№ 1 (109), 2015
7. Lagarias J.C., Odlyzko A. M. Solving Low-Density Subset Sum problems. Journ. Assoc. Comp. Mach. 1985, no. 1 (32), pp. 229-246.
8. Korobov N.M. Approximate calculation of multiple integrals using methods in the theory of numbers. Doklady Akademii naukSSSR [Reports of the USSR Academy of Sciences]. 1957, vol. 115, no. 6, pp. 1062-1065 (in Russ.).
9. Korobov N.M. Teoretiko-chislovye metody v priblizhennom analize [Number-Theoretic Methods in Approximate Analysis]. Moscow, MTSNMO Publ., 2004, 288 p.
10. Press W.H., Teukolsky S.A., Vetterling W.T. Numerical Recipes: The Art of Scientific Computing. NY, Cambridge Univ. Press, 2007, 1262 p.
11. Bykovsky V.A. On the error-theoretic quadrature formulas. Doklady RAN [Reports of RAS]. 2003, vol. 389, no. 2, pp. 154-155 (in Russ.).
12. Gauss C.F. 1801. DisquisitionesArithmeticae. Leipzig (Eng. ed.: by A.C. Clarke. New Haven, CT, Yale Univ. Press, 1966).
13. Vallee B. Une Approche Geometrique de la Reduction de Reseaux en Petite Dimension [A Geometric Approach to Reduction of Networks in Small Dimension]. PhD thesis, Universite de Caen Publ., 1986.
14. Semaev I. A 3-dimensional lattice reduction algorithm. Lecture Notes in Computer Science. Springer-Verlag Publ., 2001, no. 2146, pp. 181-193.
15. Ryshkov S.S. On the theory of positive quadratic forms on Hermite-Minkowski. Issled. po teorii chisel. 2, Zap. nauchn. sem. LOMI [Research on the Theory of Numbers. 2. Proc. of Research Workshops LOMI], 33. Leningrad, Nauka Publ., 1973, pp. 37-64 (in Russ.).
16. Afflerbach L., Grothe Н. Calculation of Minkowski-reduced lattice bases. Computing. 1985, no. 3-4 (35), pp. 269-276.
17. Lenstra A., Lenstra H., Lovasz L. Factoring polynomials with rational coefficients. Mathematische Annalen. 1982, no. 261(4), pp. 515-534.
18. Schnorr C.P., Euchner M. Lattice basis reduction: Improved practical algorithms and solving subset Sum Problems. Fundamentals of Computation Theory. 1991, pp. 68-85.
19. Kuzmin O.V., Usatyuk V.S. Software package for integer lattices basis reduction. Programmnyeprodukty i sistemy [Software & Systems]. 2012, no. 4 (100), pp. 180-183 (in Russ.).
20. Longley J.W. Modified Gram-Schmidt process vs. classical Gram-Schmidt. Communicationsin Statistics -Simulation and Computation. 1981, no. 10 (5), pp. 517-527.
21. Higham N.J. Accuracy and Stability of Numerical Algorithms Society for Industrial and Applied Mathematics. 2002, 711 p.
22. Runger G., Schwind M. Comparison of Different Parallel Modified Gram-Schmidt Algorithms. Euro-Par 2005 Conf LNCS 3648, pp. 826-836.
23. Milde В., Schneider М. Parallel Implementation of Classical Gram-Schmidt Orthogonalization on CUDA Graphics Cards. TU Darmstadt Fachbereich Informatik Kryptographie und Computeralgebra. 2009, 4 p. Available at: http://www.cdc.informatik.tu-darmstadt.de/~mischnei/ CUDA_GS_2009_12_22.pdf (accessed June 20, 2014).
24. Hermite C. Extraits de lettres de M. Hermite M. Jacobi sur differents objets de la theoriedes nombres. Reine Angew Math.1850, no. 40, pp. 279-290.
25. Korkine A., Zolotareff G. Sur les formes quadratiques. Maht. Ann. 1873, no. 6, pp. 366-389.
26. Usatyuk V.S. The implementation of the parallel orthogonalization algorithms in the shortest integer lattices basis problem. Prikladnaya diskretnaya matematika. Prilozhenie [Applied Discrete Mathematics. Application]. 2012, no. 5, pp. 120-122 (in Russ.).
27. Schnorr C.P. Shevchenko T. Solving Subset Sum Problems of Density close to 1 by "randomized" BKZ-reduction. IACR Cryptology ePrint Archive. 2012, 620 p.
28. Goldstein D., Mayer A. On the equidistribution of Hecke points. Forum Mathematicum. 2003, vol. 15, no. 2, pp. 165-189.
29. Kuzmin O.V., Usatyuk V.S. Kombinatornaya otsenka moshchnosti mnozhestva, soderzhashchego kratchayshiy vektor, obrazovannogo bazisom reshetki, privedennoy blochnym metodom Korkina-Zolotareva [Upper and Lower Bound on the Cardinality Containing the Shortest Vectors in a Lattice Reduced by Block Korkin-Zolotarev Method]. Available at: http://arxiv.org/abs/1302.4062 (accessed June 20, 2014).
30. Schnorr C.-P., Horner H.H. Attacking the Chor-Rivest Cryptosystem by Improved Lattice Reduction. Springer-Verlag Publ., 1995, vol. 921 of LNCS, pp. 1-12.
31. Gama N., Nguyen P., Regev O. Lattice Enumeration using Extreme Pruning. Advances in Cryptology -EUROCRYPT 2010. 29th Conf. on the Theory and Applications of Cryptographic Techniques. 2010., Springer Publ., pp. 257-278.
32. Kerrisk M. The Linux Programming Interface: A Linux and UNIX System Programming Handbook. No Starch Press, 2010, 1552 p.
33. Application fplll. Available at: http://perso.ens-lyon.fr/damien.stehle/fplll (accessed June 20, 2014).
34. Application GPIENUM. Available at: http://homes.esat.kuleuven.be/~jhermans/gpuenum/ (accessed June 20, 2014).
35. Application parENUM. Available at: http://xpujol.net/ (accessed June 20, 2014).
36. Ideal Lattice Challenge. Available at: http://www.latticechallenge.org/ideallattice-challenge/index.php (accessed June 20, 2014).
62