CC BY
32
УДК 519.688
ЭКСПЕРИМЕНТЫ С ПАРАЛЛЕЛЬНЫМ АЛГОРИТМОМ ВЫЧИСЛЕНИЯ
ПРИСОЕДИНЕННОЙ МАТРИЦЫ
© А.А. Бетин
Ключевые слова: параллельный алгоритм вычисления присоединенной матрицы; присоединенная матрица; эксперименты; кластер.
Приводятся результаты экспериментов с параллельным алгоритмом, проведенных на кластере МВС-100К в МСЦ РАН.
Проблема создания блочных алгоритмов рассматривалась Е. Kaltofen [1] и D. Wiedemann [2], и блочных рекурсивных алгоритмов — S. Watt [3]. Основная трудность таких алгоритмов заключается в появлении необратимых диагональных блоков. Известно, что если матрица A обратима, то все диагональные миноры матрицы ATA тоже обратимы. S. Watt [3] предложил обращать матрицу ATA и использовать тождество A-1 = (ATA)-1AT . Здесь AT транспонированная матрица к матрице A, у нее обратимы все главные миноры, в этом случае требуется два дополнительных матричных умножения.
Рассматривается другой подход к построению блочных рекурсивных алгоритмов. Присоединенная матрица — это транспонированная матрица алгебраических дополнений. Если детерминант матрицы обратим, то обратная матрица может быть вычислена как присоединенная матрица, деленная на детерминант [4].
Известно следующее разложение на множители обратной матрицы. Если ( A C \
A = ( ) — обратимая матрица и A — ее обратимый блок, то можно разложить
V B D ;
на множители ее обратную матрицу A-1 :
(I 1 їь. 1 C I О I О ' A-1 О '
О I О (D - BA-1C)-1 I B - О I
Применение детерминантных тождеств позволяет вычислять присоединенную матрицу с помощью аналогичного разложения присоединенной матрицы [4].
Пусть К - коммутативная область, А = (а^-) - квадратная матрица порядка п над областью К . Будем использовать следующие обозначения. - - определитель подматрицы А, полученной окаймлением углового минора матрицы А порядка к — 1 строкой і и столбцом і . Введем матрицу, построенную на этих минорах Ак = (ак -). Обозначим через А ^ ’ * подматрицу матрицы Ак , стоящую на пересечении строк в + 1,8 + 2, ••• , і и столбцов р + 1, р + 2, ••• , д , где 0 < в < і < п , 0 < р < д < п . Справедливо следующее разложение.
Пусть квадратная матрица А порядка п — к над кольцом К разбита на блоки
A
k+1;k , n k , n
AC
BD
1748
А — квадратный блок порядка в — к, (1 < к < в < п) , ак и а5 не равны нулю, тогда присоединенная матрица (ак)-га+к+1Ак+11’к’га* раскладывается на множители.
' (а5)-1ага1 —(а5)-1^С ' 0 1—1 0 1—1 Е 0 '
0 I 0 (ак)-1С —В а51 0 I
где ^ = (ак)-5+к+1Ак;+1;к’"* , С = (ак^^А^1^*, I — единичная матрица, и справедливо тождество
А5+1;5,п* = (ак )-1(а^ — ),
при этом а0 = 1 и Ад”’”* - присоединенная матрица для матрицы А. Здесь символ * обозначает присоединенную матрицу А* = А-Ме£А .
Существует обобщение этого разложения на случай вырожденного блока А.
Пусть Лга — модуль над Л размерности п, М € лгахга — оператор, действующий в Л”. Ищется обобщенная присоединенная матрица А, такая, что матрица 5 = АМ имеет ступенчатый вид, в частности, если М обратима, то 5 — единичная матрица.
Пусть полугруппа Рп образована матрицами, у которых число единичных элементов совпадает с рангом матрицы, а остальные элементы равны нулю. Полугруппа образована диагональными матрицами порядка п с элементами 0 и 1 на диагонали, |^га| = 2” и единичная матрица I — единица в ив Рп . На введем частичный порядок:
I < Р ^ .] — I € и чертой будем обозначать инволюцию на : I = I — I. Для каждой
матрицы Е € Р” определим диагональные матрицы: 1е = ЕЕ1 € ^га, Ре = ЕтЕ €
Для матрицы Е матрица 1е является левым аннулятором, а матрица Ре - правым анну-
лятором.
Обозначим множество ступенчатых матриц порядка п: 5” = {5 |3Е € Р”, 3^ € Л\0 : 5 = 1^5, ^Е = 5Ре}. Другими словами, для каждой ступенчатой матрицы 5 € 5” найдется матрица Е € Р” такая, что множества нулевых строк у матриц 5 и Е совпадают, и каждый ненулевой столбец матрицы ^Е совпадает с одноименным столбцом матрицы 5. Будем записывать: Е = Е^, 5 € 5га .
Введем определение. Отображение АеХ : Лгахга х (Л\0) ^ (Лгахга)3 х (Л\0)
(А,5,Е5 ,Л) = АеЖ*(М,4),
при п = 2к назовем расширенным присоединенным отображением пары ( М, ^о ), если оно определено рекурсивно следующим образом.
При М = 0 : Аеж*(0, ^0) = (^01, 0, 0, ^0).
При к = 0 и М = а = 0 : Аех*(а, ^0) = (^0, а, а, а).
При к > 0 и М = 0 разобьем матрицу М на четыре равных блока М = (М^), г, ,7 € {1, 2}. Пусть
(Ац, 5ц, Ец, ^п) = Аеж*(Мц, ^0), (1)
обозначим
М12 = А11М12/^, М21 = М21 1/^0, М22 = (М22 — М21Е11М12)/^0-
Пусть
(A12, Б12, E12, ^12) = Aeж^(I11Ml12, ^11), (2)
(А21, 521, Е21, ^21) = АеХ(^^21, ^11) • (3)
Обозначим
М22 = A21M^212Y12/(d11)2, = ^21^12/^11.
1749
Пусть
(4)
Обозначим
M21 — SiiY2i/dii, M22 — Si2^2i + (IiiM112diid2i — SiiE;f1A2iM212)Yi2/(dii)2.
M?2 — Mi22Yi2/(dsdii), M232 — S22 + /2iM222Yi2/ds,
L — (Idii — TiiM^Ev,) ^12^11^22/(dii)2,
Q — (Id12d21 — T21M22E22 d11)A22/ds j F — S11E21d22 + Ml2E22A22/dsj G — A21(M212E12A12 d0 + M21End12d1i) A11/(d11d0)j
Для произвольной матрицы M € Rnx” расширенное присоединенное отображение (A, S, E, d) — Aext(M, 1) определяет расширенную присоединенную невырожденную матрицу A, ступенчатую матрицу S и матрицу Es , такие, что AM — S и dE — SJe , определитель d.
В определении расширенного определенного отображения содержится алгоритм вычисления присоединенной матрицы.
На основе приведенного алгоритма был разработан и реализован параллельный алгоритм вычисления присоединенной матрицы. Опишем граф алгоритма.
Алгоритм содержит четыре рекурсивных вызова, кроме этого используются рекурсивные алгоритмы умножения матриц, умножения и деления матрицы на число. Умножение и деление матрицы на число эффективнее делать сразу в рекурсивном умножении матриц. Каждому алгоритму соответствует тип вершины.
Граф алгоритма состоит из вершин пяти типов:
1 - главная вершина — корневая вершина дерева алгоритма (A, S, Es, d) — Aext(M, d0) ;
2 - вершина типа A * B ;
3 - вершина типа AB d0 ;
4 - вершина типа AB + CD ;
5 - вершина типа d0 + CDd0 , где A , B , C , D — матрицы, d0 , d1 — числа.
Для хранения результатов вычисления, промежуточных результатов и результатов вычисления вершин используется массив объектов results, для хранения входных данных используется массив объектов params.
Все деревья, исходящие из вершины, разбиты на пучки. Деревья в одном пучке вычисляются параллельно, пучки пронумерованы и вычисляются в соответствии со своими номерами. Пучок из одного дерева — это вычислительный блок в вершине.
A
(L + FG/(diid2i)/di2 — FA21/(diid2i) —QG/(d11d21d12) QA21/(d 11 d21)
S
Тогда
(A,S,E,d22) — Aext(M,d0).
Здесь использованы обозначения
A? — Eij Elj j — EljEij j — dij1 — E'2j Sij j j € {1, 2}-
1750
Описание вершин будем вести следующим образом: выражения в круглых скобках будут обозначать пучки деревьев, которые вычисляются параллельно, каждое число в скобках будет обозначать тип вершин, числа между скобками - номера вычислительных блоков. Для обозначения входных и выходных данных для вершины будем использовать запись {Ау, Ву, Еу, } = Рг(М, ^0) , где переменные в фигурных скобках — результат вычисле-
ния вершины, Рг(М, ^0) — означает тип вершины и входные данные.
Приведем строение каждой вершины, начиная с простых.
Вершина второго типа вычисляет матричное произведение Ах В . Она имеет строение 0(4444)1.
У этой вершины один пучок из четырех поддеревьев. Каждое поддерево — это вершина типа АВ + СД .
(4444) = ({50} = Р4(А0, В0, АЬВ2) {51} = Р4(А0, Вь А1, В3) {52} = Р^, В0, А3, В2) {5з} = Р4(А2,В1,Аз,Вз) );
матрицы Аг,Ву получены в нулевом вычислительном блоке.
Вычислительные блоки во второй вершине:
0 : разбиваем матрицы А и В на четыре равных блока А =
Ао
А3
Аі
А4
в =
Во
В2
Ві
Вз
1 : объединяем матрицы , 5 =
Вершина третьего типа вычисляет
50 5і
5з 54
ЛхБ
^о- Она имеет строение 0(5555)1.
У этой вершины один пучок из четырех поддеревьев. Каждое поддерево — это вершина типа Ар4 + Ар4 .
(5, 5, 5, 5) = ({50} = Р5(Ао, Во, Аі, В2, ^о, ^і) {5і} = Р5(Ао, Ві, Аі, Вз, ^о, ^і)
{52} = B0, ^3 В2, ^ ^і) {53} = Р5( А2, Ві, ^3 B3, 4, ^і) );
матрицы А^,В^ получены в нулевом вычислительном блоке.
Вычислительные блоки в третьей вершине:
0 : разбиваем матрицы А и В на четыре равных блока А =
Ао
А3
Аі
А4
В=
Во
В3
Ві
В4
1 : объединяем матрицы , 5 =
5о
53
5і
54
Вершина четвертого типа вычисляет АВ + СД . Имеет строение (22)1.
У этой вершины один пучок из двух поддеревьев. Каждое поддерево - это вершина типа А * В.
(2.2) = ({Я0} = Р5(А, В) {Я1} = Р5(С, Д)).
Вычислительные блоки в четвертой вершине:
1 : Я = Я0 + Я1;
Вершина пятого типа вычисляет Ар^0 + С|_о ^0 . Имеет строение (33)1.
У этой вершины один пучок из двух поддеревьев. Каждое поддерево — это вершина типа Ар^0 .
(3.3) = ({Я0} = Рз(А, В, 4, й) {Я1} = Рз(С, Д 4, Й1)).
Вычислительные блоки в пятой вершине:
1 : Я = Я0 + Я1.
Вершина первого типа имеет следующее строение:
0(1)1(33)2(211)3(22333)4(33323)5(12)6(323)7(33)8(32)9.
а
1751
Пучки деревьев для вершин первого типа следующие:
(1) = ({Alb Slb Elb dll} = P1(M, do));
(33) = ({M^} = P3(All,M 1,1,do) {M^} = P3(M2,Yll,Xl));
(211) = ({Y1} = P2(A11,M1) {A21, S21, E21, d21} = P1(M211,d11) {A12,S12,E12,d12} = Pl(*4,dll));
(22333) = ({A1} = P2(Al2,An) {Y>} = P2(Sn,Y2l, 1, -dll) {M^} = P2(A2l^2) {Y3} =
P3(S11, X6, 1 d11) {Y16} = Pl(M2,X13, 1 d12));
(33323) = ({Yis} = P3(M22,Xl0,1,dn) {Y7} = P3(Ху,Хю, 1,dn) {M222} =
P3OY2,Yl2,Xy,Xl0, -dll) {Y4} = P2(y3,M212) {Y14} = P3(Sll,Xll,dll,do));
(12) = ({A22, S22, E22, d22} = P1 (X5, ds) {Y5} = P2(X7, Y12));
(323) = ({Ml32} = P3(M22, Y22,1, -ds) {A2} = P2(A22, A21) {Y6} = P3(Xs, Y22,1, -ds)); (33) = ({Y13} = P3(X24,Xl2,1, ds) {Y15} = P3(M22,Xl2,1, ds));
(32) = ({R2} = P3(P, G, 1, dl2) {Y19} = P3(F, G)).
Вычислительные блоки в первой вершине:
0 і n = M-fe2ngth; разбиваем входную матрицу M на четыре равных блока,
Mo M1
M=
1
M=
V M3 M4 xl = -do; хіз = eiiaii; d21 = dll * dll
2 : X2 — E11M112; X4 — mvolution(En) * M^;
3: M22 * sii — Yi; Y12 — E12S12 — d12T12; X6 — E21A2i; X11 — E21A2i;
Y21 — E21S21 — d2iT2i;
4 : X9 — TiiMi2; X10 — Ei2Ai;
5 : X5 — involution^) * M222; ds — — ^-fai-?4; G — Yi-+?8;
6: Det — d22; Mi22 — Y5 + (Si2 * d2i); Y22 — E22S22 — Det * /22; L — (Ai~DlfWn;
X — /21MI2; M31 — ; M3i — ^;
7: M232 — Y6 + S22; X12 — E22A2; X14 — /2iM222;
8: P — Ad?3; F — ?-+2f-; R1 — F; R3 — P;
9: R0 — Lt?9; S — jom(M131,M132,M231,M232); A — join(R0,R1,R2,R3);
E — jom(Eii, E12, E21, E22) •
Здесь обозначены: join — объединяет четыре блока в одну матрицу, involution — под инволюцией понимается операция / — I — / .
17Б2
Таким образом граф параллельного алгоритма будет иметь следующей вид:
С параллельным алгоритмом были проведены эксперименты на кластере МВС-100К в МСЦ РАН. Кластер состоит из 1184 вычислительных модулей, каждый модуль состоит из двух четырехядерных процессоров 1п1е1Хеоп 3Ггц и 8Гб оперативной памяти.
В экспериментах использовались матрицы размеров 512 х 512 , 1024 х 1024 , 2048 х 2048 , 4096 х 4096 и 8192 х 8192 . Присоединенная матрица вычислялась в конечном поле Ър , где р < 232 . Результаты экспериментов приведены в табл. 1-3. За эффективность вычислений на к процессорах по сравнению с вычислением на п процессорах примем число
= (Т"~" - 1)/(-1)-
п
Для вычисления эффективности в процентах будем умножать ап,к на 100%. Отметим, что
£ к
при эффективность равна 0, а при — = — эффективность равна 100%.
п
Результаты экспериментов с матрицами с плотностью 100%.
Таблица 1
Время вычислений присоединенной матрицы, сек
Количество процессоров 2 4 8 16 32 64
Размер матрицы
512x512 11.4 11.571 12.521 13.512
1024x1024 73.7 56.8 55.1 52.6 53.2 60.3
2048x2048 515.9 443.9 332.8 333.3 282.2 367.0
4096x4096 3878.8 2889.4 2305.8 1708.5 1936.2 2555.4
1753
Таблица 2
Эффективность вычислений присоединенной матрицы
Количество процессоров a2,4 ^4,8 a8,16 a16,32 a32,64
Размер матрицы
512x512 - - -0.9% -7.6% -7.9%
1024x1024 29% 3% 4% -1% -11%
2048x2048 16% 33% 0% 18% -23%
4096x4096 34.0% 25% 34% -11% -24%
Таблица 3
Время вычислений и эффективность вычислений присоединенной матрицы
с плотностью 1 %
Количество процессоров 2 4 8 16 a2,4 a4,8 a8,16
Размер матрицы
512x512 3.737 4.572 5.098 5.500 -18.3% -10.2% -7.4%
1024x1024 44.311 36.153 27.919 30.507 22.0% 29.8% -2%
2048x2048 356.581 263.753 197.893 188.559 35.0% 33.0% 5.0%
4096x4096 3855.991 1951.726 1618.777 1546.455 100% 20.0% 4.0%
8192*8192 5499.869 4848.779 4329.546 2915.142 13% 12% 49%
Результаты экспериментов показали, что с увеличением количества процессоров эффективность вычислений снижается, но тем не менее время вычислений уменьшается. Эффективность вычислений достигает 100%, что говорит об эффективном распараллеливании алгоритма.
Увеличение размеров задачи приводит к увеличению эффективности вычислений, но для больших задач возникает проблема нехватки памяти. Для решения этой проблемы было решено хранить матрицы не в оперативной памяти а в виде файлов на жестком диске. Проводятся эксперименты с такими матрицами.
ЛИТЕРАТУРА
1. Kaltofen E. Analysis of Coppersmith’s block Wiedemann algorithm for the parallel solution of sparce linear systems // Applied algebra, algebraic algorithms and error-correcting codes. 1993. V. 673. P. 195-212.
2. Wiedemann D.H. Solving sparce linear equations over finite fields // IEEE Trans. Inform. Theory. 1986. V. 32 (1). P. 54-62.
3. Strassen V. Gaussian Elimination is not optimal // Numerische Mathematik. 1969. V. 13. P. 354-356.
4. Малашонок Г.И. Матричные методы вычислений в коммутативных кольцах: Монография. Тамбов: Изд-во ТГУ им. Г.Р. Державина, 2002. С. 78-82.
БЛАГОДАРНОСТИ: Работа выполнена при финансовой поддержке АВЦП «Развитие научного потенциала высшей школы (2009-2010 годы)» (проект № 2.1.1/1853); темплана 1.12.09.
Поступила в редакцию 28 августа 2010 г.
Betin A. A. Experiments with parallel algorithm for calculation of the adjoint matrix.
Experiments with parallel algorithms for computing the characteristic polynomials of integer and polynomials matrices on the MVS100k clusters are stated and discussed.
Key words: calculation of an adjoint matrix, parallel algorithm, cluster.
1754
