Аппроксимация матрицы с положительными элементами матрицей единичного ранга

Большинство современных математических методов решения задач естествознания, техники, экономики требуют решения линейных задач большой размерности. Для понижения вычислительной сложности используется специальная структура матриц, соответствующих этим задачам. За данной структурой лежит следующий общий физический смысл: (1) строки и столбцы таких матриц ассоциированы с некоторыми элементами в пространстве, (2) задана некоторая функция взаимодействия этих элементов, (3) если функция взаимодействия асимптотически гладкая, то взаимодействие разнесенных в пространстве групп элементов можно приблизить малым числом параметров [1] (критерий разделения). Таким образом, блоки, ассоциированные с хорошо разделёнными в пространстве группами элементов, обладают малым рангом.

Другой известный пример блочно-малоранговых матриц связан с матрицами, полученными при дискретизации дифференциальных уравнений. Известно [2, 3], что если матрица A получена при конечно-элементной дискретизации дифференциального уравнения, удовлетворяющего некоторым ограничениям [2-4], то обратная к ней приближается блочно-малоранговой матрицей.

Блочно-малоранговые матрицы представляют из себя приближение с хорошей точностью плотных матриц в малопараметрическом формате. Блоки малого ранга представляются в виде произведения матриц меньшего размера. Это позволяет значительно экономить машинную память. Например, в отличие от плотной ( m^n ) матрицы, для представления которой требует m • n элементов, матрица единичного ранга требует m + n элементов. Другой характерной особенно-

Панюков А.В., Чалуб Х.З., Аппроксимация матрицы с положительными Мезал Я.А. элементами матрицей единичного ранга стью малопараметрического представления является быстрая процедура умножения такой матрицы на вектор ( O ( m + n ) операций вместо O ( mn ) ). Быстрая процедура умножения матрицы на вектор позволяет эффективно применять итерационные методы к решению систем с малопараметрическими матрицами.

Методы приближенной факторизации блочно-малоранговых матриц могут быть применены для приближенного решения и предобуславливания систем с плотными матрицами [5–14] в задачах аэро-, гидро- и электродинамики, а также в прикладной статистике и логистике.

Известны алгебраические методы построения малопараметрических представлений матриц, основанные на малоранговых аппроксимациях отдельных блоков [1, 2]. В данной работе рассмотрен эффективный способ аппроксимации блоков матрицы с положительными элементами матрицей единичного ранга, т. е. в виде произведения столбца на строку.

В первом разделе дана постановка аппроксимирующей задачи. Во втором разделе предложен способ сведения аппроксимирующей задачи к транспортной задаче в матричной постановке. В третьем разделе предложен алгоритм решения аппроксимирующей задачи и доказано, что его вычислительная сложность равна O ( mnlog ( mn ) ) .

Постановка аппроксимирующей задачи

Пусть дана матрица Л = {Xj > 0: i е I, j е J}. Рассматриваемая задача состоит в нахождении таких матриц А = {ai> 0: i е I} и В = {Pj > 0; j е J}, что Л = A • BT , т. е. разложения матрицы в произведения столбца и строки. Данная задача эквивалентна нахождению решения системы алгебраических уравнений Xj = aiPj, i е I, j е J. Понятно, что при произвольных положительных значениях Xj данная система уравнений может оказаться несовместной.

Для построения аппроксимирующей задачи воспользуемся методом наименьших модулей [15–17]. Введем функцию

Fл(a, в ) =

1 pwi

I iе I, jе J

1 ав log

X

значение которой равно среднему значению модулей логарифмов отношения приближенного значения элемента матрицы к точному значению этого элемента.

Очевидно, что inf F ^ = 0 тогда и только тогда, когда система уравнений совместна. Из неотрицательности функции F _Л следует, что значение inf F _Л можно рассматривать как степень несовместности системы. При X >  0 функция F ( 2 ) является непрерывной в окрестности любого минимума, поэтому инфимум достигается, а оптимальным приближенным решением системы с минимальной степенью несовместности можно считать

( а ⁰ , в ⁰ ) е Arg

min

{Pj > 0: j е J } { а- > 0: ге I }

V 1 ^a i ^P j

I log

- е I , j е J       ^X j

Легко заметить, что из оптимальности решения (αo,βo) следует оптимальность множества ре- шений D = {(ao • c, Po /c): c > 0}. Мы будем считать решением аппроксимирующей задачи

(а, в ) =

arg

min (а,в ) е D

II^{( а} , ^{p )}l L

Заметим, что если ( а, в ) е D , то

* α

a k /max a i • max в j

V -е I      jе J      ,

—------------:k е I max α iе I i

> , в * =

в k /max a i • max в j

\ i е I      j е J

<---------------------------------

max β j 'е J j

.

Таким образом, корректная постановка аппроксимирующей задачи является двухуровневой, но для ее решения достаточно найти любое решение задачи (1) нижнего уровня.

Математика

Сведение аппроксимирующей задачи к транспортной задаче в матричной постановке

Поскольку логарифмическая функция монотонна и

( log

I

ав

^X i

^

= 0 ^

(-log Xy + logoi + log Pj = 0)

^ ( - a iy + x i - y j = 0, l y = log X j . X = log a . У у = log j . i e I . j e J .

то далее рассматриваем задачу аппроксимации в терминах x, y :

( x⁰ ,y⁰ ) = arg min Z |^{- l} y ⁺ x i ^- У у |. ^{X e}RJ i e I , je J y e R | ^J |

Данная задача эквивалентна задаче линейного программирования

Z w ij ^ min.

ie I , je J         ^x , ^y , w

^- w ij - ^{- l} y ⁺ x i ⁺ У j - w y , w y - ⁰, i^e I , ^{j e J} , которая в стандартной форме имеет вид

Z wij^min, ie I, je J         x,y, w xi+ Уj+ wy - ly,-xi- Уj+ wy - -lij, wy - 0, ie I,j e J.

Построим задачу, двойственную задаче (18)–(19):

Z ^l ij ⁽ f ij - f ji ) e I, j e J

Z ( f j - f j^ = 0.

j e J

Z (fy-f,.) = 0.

eI fy + fji - 1, fy . fji - 0,

Сделав в задаче (5)–(7) замену переменных

^ max , f	(5)
i e I .	(6)
j e J .	(7)
i e I . j e J .	(8)

g ij = f ij ^- fji . i^e I . ^{j e J} .

получим задачу линейного программирования транспортного типа в матричной постановке

Z ^l g ie I . je J

Z g ij = 0. i^{e 1} . Z g ij = 0. ^{j e J} .

je J

ie I

^ max, ^g

- 1 - g y - 1, i e I . j e J .

Для решения подобных задач большой размерности известно программное обеспечение, которое легко инкапсулируется в систему MS Office. Данное программное обеспечение вместе с решением прямой задачи находит решение соответствующей двойственной задачи:

TAr . s . t ) = Z ( tij ⁺ t ji Hmin, i e I . j e J                     r ’ s ’ t

ri ⁺ s j ⁺ t ij ^- t ji = ^l ij . t ij . t ji ^- 0. i^e I . ^{j e J} .

Сравнивая систему ограничений задачи (1)–(2) с системой ограничений задачи (9), легко заметить. что из допустимости базисного решения ( r,s,t ) задачи (9) следует допустимость решения

R = (x = r. у = s. w = {wy = ty + ty.: i e I. je J])

задачи (1)-(2). Более того. если ( r,s,t ) - оптимальное решения задачи (10). то R - оптимальное решение задачи (1)–(2), так как двойственные им задачи (5)–(8) и (9) имеют соответствующие оптимальные решения.

Алгоритм решения аппроксимирующей задачи

Изложенное выше позволяет предложить алгоритм аппроксимации матрицы Л в виде произведения А • В T . Как отмечалось ранее, задача (9) является транспортной задачей в матричной постановке. По ее решению легко найти решение задачи (10), являющейся двойственной к ней. Это позволяет предложить следующий алгоритм решения аппроксимирующей задачи.

Алгоритм Reduction

Вход: I , J , Л = { А у : i е I , j е J } ;

Выход: A = { a i : i е I } , В = { P _j: j е J } , F _A( а , в ) ;

Шаг 1. По матрице Л = { Х у : i е I, j е J } вычислить матрицу L = { l j = ln A ij : i е I, j е J } .

Шаг 2. Найти решения ( r,s,t ) и g пары взаимно двойственных задач (9) и (10).

Шаг 3. Положить R = ( x = r , y = s , w = { w ij = t ij + t ji : i е I , j е J } ) .

Шаг 4. Вернуть a = { a = exp ( r i' ) : i е I } , в = { в _} = exp ( s j ) : j е J } , F _Л ( a , в ) = Т _Л ( r , s , t ) ;

Конец описания алгоритма Reduction

Для решения на Шаге 2 пары взаимно двойственных задач (8) и (9) транспортного типа известно множество эффективных алгоритмов на основе симплекс-метода [18]. Однако специальный вид ограничений позволяет предложить более эффективный алгоритм. Приведем его описание.

Алгоритм LD_Finder

Вход: I , J , L = { l j : i е I , j е J } ;

Выход: r = { r_j:j е I } , s = { s j :J е J } , t = { ( t_ij,t_jj ) : i е I, J е J } , Т _Л ( r , s , t ).

Шаг 1 . (Построение матрицы L ). Для каждой строки i е I матрицы L выполнить шаги 1.1, 1.2 и 1.3, затем перейти на шаг 2.

Шаг 1.1 . Построить отсортированную строку

^{( k )}-^-12         ik^k 'e.J   /⁽¹⁾⁽²⁾< ^(I J Щ

^L L i J =         k = 1,², . ,l J I, ^{J е} J , ^l - ^l - — - ^l .

[ i j ^{( k} )                                                ij ⁽¹⁾ ij ⁽²⁾              ij ( j J

Шаг 1.2 . Положить k _

= I J I + 1

, k + =

= I J | + 1

7( ^k ± ) , /( ^k - )

= ^l j ⁽ k ± ) ^{+ l} j ⁽ k - )

, ri ~        2

.

Шаг 1.3 . Для k = 1,2,...,

I T I ттпnn^if tjtk /(^) —

I J | положить l /I-, i zix ij ^{( k} )       ij ^{( k} )

^- r .

ɶ

Шаг 2 . (Построение матрицы L ). Для каждого столбца j е J матрицы L выполнить шаги 2.1, 2.2 и 2.3, затем перейти на шаг 3.

Шаг 2.1 . Построить отсортированный столбец

L ^/[ * ][ j ] = { / ( k ^ j : k = 1,2,

...

J I   i ^{( k} ) e J / ⁽¹) < /⁽²) < < /⁽ I^D ]

. ,| J I, J ^е J ,      ^l i (1) j _ ^l i (₂) j-----^l i (! ! !) j } .

Шаг 2.2 . Положить k _

= Ш±1

, k ± =

= I I | ± 1

/ (k ± )      / (k _ )

l ( ^k + ) j ± l ( ^k - ) j

, s ^j          2

.

Шаг 2.3 . Для k = 1,2, . ,

I J I положить 7 _ ⁽ k k ) = ^/( k J

s j .