CC BY
39
УДК 519.688
ПАРАЛЛЕЛЬНЫЙ ДПФ-АЛГОРИТМ ВЫЧИСЛЕНИЯ ПРИСОЕДИНЕННОЙ МАТРИЦЫ В КОЛЬЦЕ ПОЛИНОМОВ МНОГИХ ПЕРЕМЕННЫХ С ЦЕЛЫМИ КОЭФФИЦИЕНТАМИ
© А. О. Лапаев
Ключевые слова: параллельная компьютерная алгебра, дискретное преобразование Фурье, быстрое преобразование Фурье, полиномиальные алгоритмы, вычисление присоединенной матрицы, параллельные матричные алгоритмы.
Предлагается параллельный алгоритм вычисления присоединенной матрицы над кольцом полиномов многих переменных с целыми коэффициентами, использующий дискретное преобразование Фурье над простым полем.
1 Введение
Одной из задач параллельной компьютерной алгебры является разработка эффективных матричных алгоритмов, в частности для матриц, элементами которых являются полиномы многих переменных.
Классическим методом ускорения вычислений является использование Китайской теоремы об остатках (КТО) [1] для числовых и полиномиальных модулей: задача решается в гомоморфных образах над конечным полем и затем результат интерполируется в кольцо полиномов с использованием схем Ньютона или Лагранжа [1]. Такой подход позволяет избежать сильного роста коэффициентов и степеней полиномов при вычислениях.
В данной статье предлагается вместо КТО для полиномиальных модулей использовать дискретное преобразование Фурье [2]. Преимуществом данного подхода является восстановление результата по полиномиальным модулям за меньшее количество операций. Например, в случае использования т простых полиномиальных модулей сложность восстановления по КТО будет оцениваться как 0(т2), а в случае использования ДПФ -©(тк^2 га).
В параграфе 2 приводится алгоритм вычисления присоединенной матрицы над кольцом Ъ[х\,... ,хп\, использующий дискретное преобразование Фурье (ДПФ) в конечном поле, и вычисляются выражения для сложности приведенного алгоритма. Кроме того, проводится сравнение теоретических оценок для алгоритмов, использующих ДПФ [1, 3, 4] и использующих КТО.
В параграфе 3 рассматривается распараллеливание алгоритма вычисления присоединенной матрицы.
В параграфе 4 приводятся результаты экспериментов с полученным параллельным алгоритмом вычисления присоединенной матрицы, который в конечном поле использует ДПФ.
2 Оценки сложности
Алгоритм для вычисления присоединенной матрицы с использованием ДПФ и его оценки сложности для случая полиномов одной переменной были получены в работе [5]. Теперь рассмотрим случай полиномов многих переменных.
Пусть А = (а^(х)) е Ъп*п[хи... ,ха],
ац(х ъх2,•••,*!<)= X] ... к*х1Хх22хк/-
/с 1=0 к2=0 А:^=0
Пусть шах \а^'к2,'"'к,£\ ^ а и &е%иа^х1,х2,... ,хл) < Зи,и = где degu -
степень полинома по переменной хи. Обозначим 5 = 5152 ... 5и .
Максимальное количество мономов в элементах присоединенной матрицы А* не будет превосходить 771 = . К0ЛИЧеСТВ0 точек ДЛЯ вычисления ДПФ ПО переменной Хк будет
лг* = [4, 6, 7|. Обозначим N = ...N4.
Найдем количество простых 32-битных числовых модулей. Будем использовать оценку абсолютной величины коэффициента определителя матрицы А.
Определитель матрицы А можно вычислить по формуле:
det А = ^ ( — 1) Cllj1(l2j2 • • • ап7тИ (1)
Оь-Лп)
где Ух,...,^п) - перестановка чисел от 1 до гг, ^ - четность этой перестановки. Формула (1) содержит ровно п! слагаемых. Используя формулу (1) и формулу Стирлинга для верхней оценки значения п\, получим верхнюю оценку для максимального числового коэффициента с1е1 А. Тогда количество г простых 32-битных модулей равно
г=Г108л2(^(-) е^"5п-1ап)1, (2)
где /г = 232.
Алгоритм вычисления присоединенной матрицы
1. Выбираются простые числовые модули р\,... ,рг .
2. Устанавливается номер простого модуля £ = 1.
3. Для каждого элемента матрицы А вычисляется дискретное преобразование Фурье ррАА) = №*(<*#)) = ((“«/»по ПР0СТ0МУ модулю р1 [4]. Каждый элемент такой матрицы - одномерный вектор длины N = 21об2т, где и = 1,..., N являются элементами простого поля ЪР1.
4. Для матриц ^(Л(и)) = (а“/),и = вычисляются присоединенные матрицы
^(Л-<и)) = «•').
5. Составляется матрица ^(Л*) = ,... , а- )). Для ее элементов вычисляется
обратное преобразование Фурье А*г = .Р_1(.Р*(Л*)).
6. Увеличиваем £ на 1. Если £ < г, то переходим к шагу 3.
7. По образам А\, А*2,..., А* в кольцах Ър^Ъ^х],... ,ЪРг[х\ вычисляется матрица А*, используя КТО [3] для ЧИСЛОВЫХ модулей РьР2, . • • ,Рг •
Получим выражения сложности алгоритма.
Количество операций над словами при вычислении дискретного преобразования Фурье по алгоритму быстрого преобразования Фурье для п2 полиномов на N точках при использовании г числовых простых модулей равно
п2г (7вГк^Л а] + 9ЛГк^2 -/V).
Каждая кольцевая операция над образами полиномов состоит из rN арифметических операций над словами и гМ вычислений остатка от деления. Будем считать, что вычисление остатка от деления занимает столько же времени, как 7 сложений слов. Тогда количество арифметических операций над образами полиномов исходной матрицы А при использовании алгоритма прямого хода [8] равно
8 п3гЛГ
Для восстановления полинома по его образам необходимо выполнить г обратных преобразований Фурье и восстановить п2 • чисел. Тогда число операций для получения присоединенной матрицы по ДПФ-образу ее элементов будет
тг2(9гЛг к^2 N + 2г2зп(1).
Таким образом, сложность алгоритма будет равна
n2r(7s\\ogh а"| + ЭЛ/- к^2 ЛГ) + 8 п3гЛ^ + п2(9т\/У к^2 N + 2т2 вп11). (3)
Также приведем выражение сложности классического метода, использующего КТО
как для числовых, так и для полиномиальных модулей:
n2r(Js\\ogh а] + т^в2) + 8п3гтгй5 + п2(2(п*з)2г + 2г2вп<*). (4)
Если пй5 равно степени числа 2 и если пренебречь величиной 7\\ogha]n~d по сравнению с гг, то оценки (3) и (4) примут более простой вид (5) и (6):
п2+<*гв( 18 ^2(тгйв) + 8п + 2 г), (5)
п2+аг8(з(2п<1 + 1) + 8п + 2г). (6)
Вычислим выражения (3) и (4) при некоторых значениях параметров а, й, п, с1 и представим результат в виде таблицы. Для полиномов одной переменной с целыми коэффициентами, занимающими 8 бит, получим следующую таблицу, у которой в левом столбце указан порядок матрицы (п), а в верхней строке - размеры исходных полиномов (з).
Таблица 1
Отношение числа операций алгоритмов, использующих КТО и ДПФ, при с1 = 1 и а = 256.
П \ 5 2 4 8 16 32 64
2 0,58 0,57 0,67 0,91 1,39 2,28
4 0,63 0,69 0,87 1,24 1,95 3,30
8 0,75 0,87 1,15 1,69 2,75 4,75
16 0,91 1,09 1,47 2,23 3,71 6,51
32 1,05 1,29 1,79 2,76 4,67 8,36
64 1,17 1,46 2,04 3,20 5,48 9,93
128 1,25 1,58 2,23 3,52 6,06 11,08
Табл. 1 делится побочной диагональю на две области: в верхнем треугольнике более эффективен алгоритм, использующий КТО для полиномиальных модулей, в нижнем -алгоритм, использующий ДПФ.
Приведем аналогичную таблицу для матриц, элементами которых являются полиномы двух переменных.
Из табл. 2 видно, что алгоритм, использующий ДПФ, является эффективным при вычислении присоединенной матрицы для матриц, элементами которых являются полиномы двух переменных при порядке матрицы п > 4.
Таблица 2
Отношение числа операций алгоритмов, использующих КТО и ДПФ, при с1 = 2 и а = 256
n/s 2 4 8 16 32 64
2 0,54 0,62 0,84 1,27 2,08 3,58
4 0,82 1,16 1,82 3,06 5,37 9,69
8 1,63 2,65 4,54 8,13 14,93 27,79
16 3,65 6,41 11,73 21,81 41,25 78,20
32 8,39 15,56 29,45 56,33 108,40 209,35
64 18,99 36,35 70,29 136,68 266,59 520,87
128 41,69 81,24 159,12 312,52 614,74 1210,26
3 Параллельный алгоритм
Пусть имеется параллельная ЭВМ с t процессорами с номерами q — 0,..., t — 1. Пусть на процессоре с номером 0 имеется матрица A Е Z[x\Приведем шаги параллельного алгоритма.
1. Процессор с номером 0 выполняет рассылку матрицы А на остальные процессоры.
2. Процессоры с номерами, отличными от 0, выполняют прием матрицы А.
3. Каждый процессор оценивает количество простых модулей г и максимальную степень каждой переменной Х{ в элементах присоединенной матрицы А*.
4. Для упрощения описания предположим, что г делится на t нацело. В случае, если это не так, первые г mod t процессоров вычисляют присоединенную матрицу по [jJ + 1 модулям, а остальные - по |_yj модулям.
5. На каждом процессоре выполняется разбиение массива простых чисел на t частей согласно условию, описанному на предыдущем шаге.
6. На каждом процессоре выполняется вычисление ДПФ для элементов матрицы последовательно по каждому простому числовому модулю данного процессора. В случае, если г < t, оставшиеся процессоры используются для распараллеливания вычисления
ДПФ.
7. Для полученной матрицы образов ДПФ вычисляется присоединенная Аф по каждому простому числовому модулю на данном процессоре.
8. Для элементов матрицы Adft вычисляется обратное преобразование Фурье по каждому простому числовому модулю на данном процессоре.
9. Полученные матрицы Аф делятся на t частей, и каждый процессор посылает часть с номером q на процессор с номером q.
10. Каждый процессор восстанавливает коэффициенты своей части полиномов по Китайской теореме об остатках.
11. Каждый процессор посылает свою часть восстановленных на шаге 10 данных процессору с номером 0.
12. На процессоре 0 выполняется объединение полученных от других процессоров частей присоединенной к А матрицы.
Приведем псевдокод данного алгоритма.
//Данный код выполняется на всех процессорах кластера. matrixSend(A); //Рассылка матрицы А на все процессоры.
/*Вычисление количества простых модулей и максимальной степени по каждой переменной.*/ getBounds();
^Инициализация структуры, в которой будут хранится промежуточные данные*/ pcs = new PolCoeffStore[primesLocalAmount[rank]];
/*В элементе массива primesLocalAmount с номером j хранится количество простых чисел для процессора с номером j. Переменная rank хранит номер данного процессора.*/ for (int i = 0; i < primesLocalAmountfrankJ; i++) {
//В массиве primes находятся простые модули. int p = primesfij;
//Вычисляем ДПФ-образ матрицы. matrFft = fftComputation(matr, p);
//Вычисляем присоединенную матрицу. adjMatrFft = adjointComputation(matrFft);
//Вычисляем обратное ДПФ для присоединенной матрицы. pcsfi] = fftInverseComputation(fftImRes, р);
}
/*Выполняется восстановление результата по КТО и его сборка на процессоре 0.*/ buildResult();
//Процессор с номером 0 возвращает результат. if (rank == 0) return result; else return null;
4 Эксперименты
Рассмотренный алгоритм был программно реализован. С полученными программами были проведены эксперименты на кластере МВС100К МСЦ РАН. Результаты экспериментов представлены в табл. 3. В таблице приведено время вычислений и эффективность, которая вычислялась по формуле:
hn __ 1
_ Jk_____
Чт,к ~ к_ _ т
Здесь tm, tk ~ время вычислений на т и к процессорах соответственно.
Из табл. 3 видно, что алгоритм обладает хорошей масштабируемостью, т. к. время вычислений уменьшается при увеличении количества процессоров. Например, при при увеличении количества процессоров с 2 до 16 для матриц порядка п = 16 время уменьшилось в 6 раз.
Таблица 3
Время вычисления присоединенной матрицы с использованием к процессоров и эффективность алгоритма для матриц порядка п с полиномами двух переменных степени 4 по каждой с
64-битными коэффициентами
к (процессоры) 2 4 8 16
п (порядок матрицы) Время,с Время,с / (72,4 Время,с / <72,8 ВреМЯ,С / 02,16
4 3 2/50% 2/17% -
8 38 18/111% 11/82% -
16 1457 661/120% 362/101% 230/76%
5 Заключение
В статье описан алгоритм параллельного вычисления присоединенной матрицы в кольце Z[ii,... ,Xd]. Приведены выражения сложности данного алгоритма для случая, когда элементами матрицы являются полиномы многих переменных. Из данных оценок следует, что использование ДПФ имеет преимущество при восстановлении результата. Получена схема распараллеливания данного алгоритма для суперЭВМ с распределенной памятью. Проведены эксперименты на кластере МСЦ РАН, которые показали, что при увеличении количества процессоров в 8 раз время уменьшается приблизительно в 6 раз, что говорит о хорошей масштабируемости рассмотренной схемы распараллеливания.
Следующим этапом исследований является построение и реализация параллельных ДПФ-алгоритмов для умножения матриц, вычисления характеристического полинома и определителя.
ЛИТЕРАТУРА
1. Кнут Д.Э. Искусство программирования. Т.2. Получисленные алгоритмы, 3-е изд. М.: Издательский дом «Вильямс», 2001.
2. Лапаев А.О. Вычислительные алгоритмы для полиномиальных матриц с использованием ДПФ // International conference Polynomial Computer Algebra. St.Petersburg, PDMI RAS, 2009. C. 141-146.
3. Нодеп П., Kumme К. Алгебраическая алгоритмика (с упражнениями и решениями)/ Пер. с франц. М. Мир, 1999.
4. Лапаев А.О. О вычислении многомерного дискретного преобразования Фурье в простом поле // Вестник Тамбовского университета. Сер. Естественные и технические науки. Том 14. Вып. 4. 2009. С. 729-731.
5. Malaschonok G. Valeev Yu. D., Lapaev A. О. On the choice of multiplication algorithm for polynomials and polynomial matrices // Zapiski POMI. 2009. V. 373. P. 157-188.
6. Кормен, Томас X., Лейзерсон, Чарльз И., Ривест, Рональд Л., Штойн, Клиффорд. Алгоритмы: построение и анализ, 2-е издание/ Пер. с англ. М.: Издательский дом «Вильямс», 2005. С. 1296.
7. Лапаев А.О. Параллельное вычисление дискретного преобразования Фурье полинома в простом поле. Материалы 9-й международной конференции-семинара «Высокопроизводительные параллельные вычисления на кластерных системах». Владимир, 2009. С. 272-273.
8. Malaschonok G.I. Algorithms for computing determinants in commutative rings // Diskret. Mat. 1995. Vol. 7. No. 4. P. 68-76.
БЛАГОДАРНОСТИ: Работа выполнена при поддержке программы «Развитие потенциала высшей школы» (проект 2.1.1/10437).
Поступила в редакцию 12 ноября 2010 г.
UDK 519.688
PARALLEL DFT-ALGORITHM FOR COMPUTING OF ADJOINT MATRIX IN THE RING OF INTEGRAL POLYNOMIALS WITH MANY VARIABLES
(c) A. O. Lapaev
Key words: parallel computer algebra, discrete Fourier transform, fast Fourier transform, polynomial algorithms, computation of adjoint matrix, parallel matrix algorithm.
Parallel algorithm for computing of adjoint matrix over Z[xi,...,xn] in terms of discrete Fourier transform over finite field is described.
