Неполное LU-разложение

Предобусловливание основанное на классических методах

Ранее было показано, что методы Якоби и Гаусса-Зейделя, а также SOR и SSOR могут быть представлены в виде (2.8). Другой формой этого выражения является

x k +1=Gx k +f,

(2.44)

которое получается из (2.14) домножением левой и правой частей на K ^–1. Проведя аналогию с (2.1) и (2.2), нетрудно увидеть, что (2.44) можно рассматривать как итерационную схему, построенную для решения СЛАУ

(I–G)x=f,

(2.45)

в которой

f=K^–1b; G=K^–1R=K^–1(K–A)=I–K^–1A.

Таким образом, система (2.51) эквивалентна системе

K^–1Ax=K^–1b,

т.е. предобусловленной СЛАУ (2.1) с матрицей предобусловливания K. Для методов Якоби, Гаусса-Зейделя, SOR и SSOR эти матрицы соответственно, равны (K заменено на M):

M_J=D;

M_GS=D–E;

M_SOR (D–wE);

M_SSOR (D–wE)D^–1(D–wF).

Скалярные множители и при практических реализациях SOR и SSOR-предобусловливания обычно не учитываются, так как дают лишь масштабирование, практически не влияющие на скорость сходимости.

При w=1 получается симметричный метод Гаусса-Зейделя (SGS)

M_SGS=(D–E)D^–1(D–F).

Заметив, что (D – E) соответствует нижнетреугольной части матрицы A, а (D – F) верхнетреугольной части, то

M_SGS=LU,

где

L=(D–E)D^–1=I–ED^–1, U=D–F.

Таким образом для решения систем вида w = M _SGS^–1 x необходимо последовательно решить

(I–ED^–1)z=x, (D–F)w=z.

Исторически предобусловливание, как собственно и итерационные методы вообще, использовалось для решения разреженных СЛАУ большой размерности. Поэтому использование предобусловливания для решения СЛАУ с плотными матрицами необходимо использовать с некоторыми добавлениями о которых речь пойдет позже.

Для систем с разреженными матрицами хорошим предобусловливанием является неполное LU-разложение. Данное разложение заключается в представлении матрицы в виде

M=LU+R,

(2.46)

где L и U как и прежде нижне- и верхнетреугольные матрицы, соответственно, а матрица R является матрицей ошибки. Тогда приближенное представление M» LU называется неполным LU-разложением матрицы M или коротко её ILU-разложением. Как хорошо известно, в результате реализации стандартного LU-разложения на месте исходной матрицы находится результат разложения, т.е. матрицы L и U, одна из которых имеет единичную диагональ.

Наиболее простым вариантом является ILU(0) разложение. Оно заключается в применении LU-разложения к матрице M, но если m_ij = 0, то полагается l_ij = 0 или u_ij = 0. Поскольку в рассматриваемой задаче матрица M плотная, то применение ILU(0) не целесообразно, без предварительного процесса разрежевания.

С учетом вышесказанного представление матрицы предобусловливания в виде удобном для использования итерационного метода с неявным преобусловливанием (M» LU) требует осуществления двух этапов. На первом этапе определяется структура разреженности матрицы М (путем предфильтрации матрицы А), т.е. m_ij = a_ij, если пара индексов (i, j) соответствует предписанной структуре и m_ij = 0 в противном случае (данный этап более подробно рассмотрен ниже). На втором этапе производится собственно само ILU-разложение, в частности ILU(0).

Алгоритм ILU(0) разложения

Для i =1,… N

w= a_{i *}

Для k =1,… i –1и w_k ¹0

w_k = w_k / u_kk

Для j = k +1,… N

w=w– w_k × u_{k *}

Увеличить j

Увеличить k

l_ij = w_j, для j =1,… i -1

u_ij = w_j, для j = i,… N

Увеличить i

Очевидно, что построенная таким образом матрица LU удовлетворяет всем требованиям к матрице предобусловливателя. Действительно, она близка к матрице A; она легко вычислима по приведенному выше не сложному алгоритму; наконец, она легко обратима, так как является произведением двух треугольных матриц. Таким образом, выбор M = LU является достаточно хорошим способом предобусловливания.

В данной версии системы реализовано нахождение матрицы предобусловливания как с помощью ILU(0), так и с помощью стандартного LU-разложения.

Дата добавления: 2014-12-27