num-meth

1.1.5.1QR-разложение

Один из широко используемых видов мультипликативного разложения – это разложение исходной матрицы в произведение ортогональной Q и правой верхнетруегольной R матриц. Напомним, что, по определению, особенность ортогональной матрицы Q состоит в том, что Q−1 = QT, т.е. матрица, обратная к ортогональной, совпадает с транспонированной ортогональной.

Один из методов построения QR-разложения состоит в ортогонализации и нормировании ее столбцов. Для системы векторов-строк матрицы A применяется процедура ортогонализации и нормирования. Заметим, что процедура ортогонализации рассчитана на линейно-независимую систему векторов. Линейная независимость строк матрицы A является условием и критерием существования и единственности решения исходной системы.

Представим матрицу A как систему векторов-строк A = (ˆa1, . . . , aˆn) и рассмотрим процесс Грамма-Шмидта ортогонализации этой системы строк. В случае, если система векторов (ˆa1, . . . , aˆn) линейно-независима, то в результате для этой системы будет получен

ˆ

(ˆa2, b1)

Откуда получаем: β21 = − ˆ ˆ .

(b1, b1)

ˆ ˆ

Предположим, что на i-м шаге уже построена система векторов b1, . . . , bi−1. Тогда:

не нулевые т.к. каждый из этих векторов построен как линейная комбинация векторов, взятых из линейно-независимой системы (ˆa1, . . . , aˆn), причем заведомо не все из коэффициентов этой линейной комбинации равны нулю.

Процедура ортогонализации Грамма-Шмидта, с небольшими модификациями, может быть применена и к системе линейно-зависимых векторов. В этом случае на некотором

ˆ

шаге будет получен нулевой вектор bi. Такой вектор отбрасывается, и процесс продолжается дальше. В результате будем получена система попарно ортогональных векторов (ˆa1, . . . , aˆk), где k ≤ n, причем легко видеть, что k совпадает с рангом системы векторов (ˆa1, . . . , aˆn) и, следовательно, с рангом исходной матрицы A.

Описанный процесс дает лишь ортогональный, но не нормированный базис системы векторов (ˆa1, . . . , aˆn). Переход к ортонормированному базису (ˆq1, . . . , qˆn) осуществляется

ˆ

посредством нормирования каждого из векторов bi, i = 1, n:

Теорема 1.1.2. Любую квадратную матрицу порядка n, у которой все угловые диагональные миноры

отличны от нуля, можно представить, и единственным образом, в виде LU-разложения.

Запишем выражения для коэффициентов матриц L и U (вывод приводимых равенств получается в процессе доказательства предыдущей теоремы):

u1j = a1j, j = 1, n

li1 = ai1 , i = 2, n

u11

∑i−1

uij = aij − likukj, i = 2, n, i < j

k=1

1

∑j−1

lij = ujj (aij − likukj), i = 2, n, i > j

k=1

1.1.5.3QRS-разложение

Применение QRS-разложения основывается на идее умножения исходной матрицы справа и слева на одну или несколько ортогональных матриц:

MpMp−1 × · · · × M1AN1 × · · · × Nq−1Nq

Произведение ортогональных матриц – также ортогональная матрица. Поэтому матрица A представима в виде мультипликативного разложения следующего вида:

A = (MpMp−1 × · · · × M1)−1R(N1 × · · · × Nq−1Nq)−1 = QRS

Причем в этом выражении матрицы Q и S – ортогональные.

1.1.5.4Полярное разложение

Полярным разложением квадратной матрицы A называется ее такое мультипликативное разложение

A = SP,

в котором матрица S – симметрическая, а матрица P – ортогональная. Прежде, чем указать, как вычисляются члены разложения, введем несколько определений.

1.1.5.5Каноническое (спектральное) разложение

Определение 1.1.2. Каноническим (спектральным) разложением матрицы M

называется представление этой матрицы в виде произведения:

N = T ΛT −1,

0λn

ния матрицы M, а матрица T составлена из собственных векторов-столбцов матрицы M по соответствующим собственным значениям.

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

Теперь доопределим операцию возведения матрицы в степень для случая, когда показатель степени 0 < α < 1. Арифметическим корнем степени α из неотрицательной матрицы M называется такая матрица N, что N2 = M. Чтобы вычислить матрицу N, нужно найти спектральное разложение матрицы M = T ΛT −1. Затем, полагая по определению:

√

будем считать, что N = T Λ T −1.

Теперь можно написать выражения для сомножителей в полярном разложении матрицы A = SP .

√

S = AAT

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

P = S−1A = S(AT)−1

Матрица P ортогональная и определена единственным образом, если исходная матрица M невырождена.

1.2Решение системы линейных алгебраических уравнений (СЛАУ)

1.2.1Прямые методы решения СЛАУ

Прямыми называются такие методы решения СЛАУ, которые приводят к решению за конечное число арифметических операций и, в случае точности вычисления кажой операции, гарантирует получение точного решения.

1.2.1.1Метод Гаусса и его модификации

Пусть задана СЛАУ:

В матричной форме эта система записывается в виде Ax = h.

16 Глава 1. Численные методы линейной алгебры

Идея метода состоит в том, чтобы последовательно, в ходе серии шагов, исключать переменные. Суть алгоритма состоит в следующем.

Пусть элемент a11 отличен от 0. Тогда умножим первую строку 1.22 на величину a21 ,

a11

а результат вычтем из второго уравнения. Аналогичную процедуру неоходимо выполнить для остальных (n − 2) уравнений (для i-го уравнения нужен первое уравнение умножать

Элемент a11 называется ведущим элементом шага. Заметим, что если в системе 1.23 обнаружатся нулевые строки, то их необходимо вычеркнуть.

Следующий шаг выполняется аналогично. Элемент a′22 ≠ 0 принимают за главный, а из всех уравнений, расположенных ниже второго, вычитают второе уравнение, умноженное

на величину a′23 . И снова, нулевые уравнения вычеркиваются, если они возникают. В

a′22

результате, после не более, чем n, шагов исходная система примет вид:

Матрица этой системы, в общем случае, трапецевидная, а если r = n, то она правая верхнетреугольная. Процесс, в результате которого получается эта матрица, называется прямым ходом метода Гаусса. Очевидно, если матрица треугольная, то исходная система имеет ровно одно решение. В противном случае система имеет много решений, переменные xr+1, . . . , xn принимаются за свободные, а остальные – за главные.

Процесс получения значений главных переменных называется обратным ходом метода Гаусса. В этом процессе на первом шаге вычисляется переменная xr, а результат подставляется в остальные уравнения. Затем то же проделывается с переменной xr−1 и т.д. до тех пор, пока все переменные не получат либо численного значения, либо выражения через свободные переменные.

На практике метод Гаусса используется для решения невырожденных систем с квадратными матрицами коэффициентов при переменных. В процессе решения, для проверки правильности, результат каждого шага заносят в таблицы.

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

1.2.1.1.1Модификации метода Гаусса

. Один из недостатков описанной выше схемы состоит в накоплении ошибок на каждом шаге. Это особенно заметно проявляется, когда ведущий элемент очередного шага прямого хода оказывается близок к 0. Один из способа преодоления этой проблемы заключается в схеме, известной под названием метод Гаусса с выбором главного элемента.

Существуют две разновидности этого метода – с частичным выбором главного элемента, и с полным. Рассмотрим сначала первый вариант.

При частичном выборе на k-м шаге среди всех неприведенных коэффициентов при переменной xk выбирается наибольший по модулю. После этого выполняется перестановка строк матрицы таким образом, чтобы выбранный элемент оказался на позиции главного.

Можно выбирать главный элемент не по столбцу, а по строке (в этом случае потребуется перестановка столбцов матрицы, чтобы выбранный элемент оказался на нужной позиции).

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

Вместе с тем, описанные дополнения к методу Гаусса сами по себе влекут дополнительные расчеты. Поэтому если выполнять перестановку слишком часто, то существенно возрастет общее время расчетов. Этого можно избежать, если ввести параметр ε > 0. Процедуру поиска главного элемента нужно выполнять только если на k-м шаге элемент (k, k) окажется по модулю меньше этого параметра.

1.2.1.2Метод квадратного корня

1.2.1.3Применение мультипликативных разложений матриц в решении СЛАУ

Применение мультипликативных разложений к решению СЛАУ основывается на следующей идее. Пусть задана СЛАУ:

Ax = b

Пусть далее, для матрицы A существует разложение A = A1 ×· · ·×An. Тогда исходная система примет вид:

A1 × · · · × Anx = b

Далее, выполняется одно из двух эквивалентных преобразований системы. Во-первых, можно выполнить последовательное умножение уравнения слева на матрицу A−1 1, A−2 1 . . .

A−k 1 (при этом сомножители выбираются так, чтобы нахождение этой обратной матрицы было возможно более простой задачей). В результате система приобретает вид:

Ak+1 × · · · × Anx = b(k),

где b(k) = A−k 1 × A−k−11 × . . . × A−1 1b– вектор свободных членов, . Если k = 0, то b(k) = b(0) = b. Второй способ заключается в замене переменных. Пусть после k − 1 шагов (если

k = 1, то имеем дело с исходной системой) система приняла вид:

Ak × · · · × Anx = b(k−1)

Тогда произведем замену переменных: y = Ak+1 × · · · × Anx и перейдем к системе:

Aky = b(k−1)

Разрешив эту систему относительно y, получим после k-го шага:

Ak+1 × · · · × Anx = b(k)

При этом b(k) = A−k 1bk−1 = A−k 1 ×· · ·×A−1 1b. Необязательно, что все эти операции выполняются напрямую; эквивалентный результат достигается по возможности более простыми преобразованиями.

Матрицы A1, . . . , An выбираются какого-либо специального вида с таким расчетом, чтобы решение задачи на каждой итерации было возможно проще. В результате вместо решения одной сложной задачи решается ряд простых. В качестве матриц, участвующих в мультипликативном разложении, берутся обычно треугольные, диагональные, двухдиагональные, трехдиагональные, ортогональные и другие виды «хороших» матриц. Далее рассмотрим кратко некоторые из основных разложений.

1.2.1.3.1Применение QR-разложения.

В случае QR-разложения матрица представляется произведением ортогональной и правой верхнетреугольной матриц. При этом исходная система Ax = b запишется в виде:

QRx = b,

а от этой системы можно перейти к следующей:

Rx = Q−1b

Т.к. Q−1 = QT, то последнее равенство можно переписать в виде:

x = QTb,

где матрица R – треугольная.

Поскольку процедура ортогонализации позволяет не только ортонормировать систему линейно-независимых векторов, но и, с некоторыми модификациями процедуры, определить ранг линейно-зависимой системы, в процессе ортогонализации можно определить, имеет ли исходная система единственное решение (если этого не было известно заранее).

Итак, в результате процедуры ортогонализации и нормирования получаем систему векторов:

q1 = u11a1,

q2 = u12a1 + u22a2,

. . .

qn = u1na1 + u2na2 + · · · + unnan

Это равносильно следующему матричному уравнению:

Q = AU,

где Q = (q1, . . . , q)n) – ортогональная матрица, U – правая верхнетреугольная:

00 . . . unn

Отсюда получаем QR-разложение матрицы A в виде A = QU−1 = QR. При этом матрицу U−1 вычислять нет необходимости. В самом деле, выполним цепочку эквивалентных преобразований исходной системы Ax = b. Представив матрицу A в виде ее QR-разложения, получим:

QU−1x = b

Далее, домножим обе части равенства на Q−1 (слева). Т.к. матрица Q – ортогональная, то Q−1 = QT и:

U−1x = QTb

Полученное равенство снова домножим слева на матрицу U и получаем искомое решение исходной системы:

x = UQTb

1.2.1.3.2Применение LU-разложения.

Если матрица A представлена в виде своего LU-разложения, то исходная система преобразуется к виду:

LUx = b

Общий ход решения этой новой системы заключается в выполнении замены z = Ux и последовательном решении двух систем:

{

Lz = b Ux = z

Именно такое LU-разложение и лежит в основе метода Гаусса. Решение первой из указанных систем равносильно прямому ходу в методе Гаусса, а решение второй – обратному ходу в методе Гаусса.

Разновидностью описанного метода является метод Холецкого. Этот метод заключается в том, что исходная матрица также представляется LU-разложением, с тем отличием, что у верхней треугольной матрицы R элементы на главной диагонали должны быть равны 1, а элементы на главной диагонали нижней треугольной матрицы S должны быть отличны от нуля:

Теперь вместо решения исходной системы необходимо решить две более более простых:

{

Sz = b Rx = z

Приведем, без доказательства, формулы для вычисления элементов матриц S и R:

1.2.1.3.3Применение полярного, спектрального и QRS-разложений.

Возьмем QRS-разло -жение матрицы A:

A = (MpMp−1 × · · · × M1)−1R(N1 × · · · × Nq−1Nq)−1 = QRS

Напомним, что в этом выражении матрицы Q и S – ортогональные. Общий ход решения полученной системы сводится к замене переменных y = Sx. В полученной системе:

QRy = b

обе части равенства домножаются слева на матрицу Q−1:

Ry = Q−1b

В силу ортогональности матрицы Q, имеем:

QRS-разложение матрицы A строится таким образом, чтобы матрица R имела какойлибо специальный «удобный» вид, упрощающий решение последней приведенной системы, например, в качестве матрицы R выбирают трехдиагональную матрицу. В таком случае задачу 1.25 можно решать методом прогонки.

В некоторых случаях различные разложения матриц удобно комбинировать. Например, исходную матрицу A можно предварительно представить ее полярным разложением.

После того, как получено полярное разложение матрицы A = SP , матрицу S можно дополнительно подвергнуть QR – разложению. Поскольку матрица S является симметрической по построению, ее можно представить в виде спектрального разложения с ортогональной матрицей T . Тогда получим:

A= T ΛT TP

Вэтом разложении все матрицы, кроме Λ, – ортогональные, а сама матрица Λ – диагональная. Обозначим N = T TP и получим:

A = T ΛN

Это выражение в точности совпадает с QRS-разложением.

1.2.1.4Метод прогонки.

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

Исходная система в таком случае имеет вид:

Выполним переобозначение: ai = ai,i−1, bi = aii, ci = ai,i+1. В этих новых обозначениях имеем: