Материал из eSyr's wiki.

Предыдущая лекция | Следующая лекция

Содержание 1 Глава 1. Численные методы линейной алгебры 1.1 Параграф 2. Связь метода Гаусса с разложением матрицы на множители (продолжение) 1.1.1 Связь метода Гаусса с разложением A = B × C 1.2 Параграф 3. Обращение матрицы методом Гаусса-Жордана 1.3 Параграф 4. Метод квадратного корня

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

Параграф 2. Связь метода Гаусса с разложением матрицы на множители (продолжение)

Решение нелинейной системы (3), (4), полученной на предыдущей лекции.

• c_ij = ^{(a_ij − Σ_l=1i−1 b_il×c_lj)}/_bii, i < j ... (3)
• b_ij = a_ij − Σ_l=1^j−1 b_il×c_lj, i ≥ j ... (4)

Первый шаг: первая строка C

• b₁₁ = a₁₁
• c_1j = ^a_1j/_b11, j = 1…m

Второй шаг: первый столбец B

• b_i1 = a_i1, i = 1…m

Третий шаг: вторая строка C

• b₂₂ = a₂₂ − b₂₁ × c₁₂
• c_2j = ^{(a_2j − b₂₁×c_1j)}/_b22, j = 2…m

Четвёртый шаг: второй столбец B

• b_i2 = a_i2 − b_i1×c₁₂, i = 2…m

Таким образом, мы решили нелинейную систему, грамотно организовав алгоритм.

Факторизация A = B × C (2) осуществляется при b_ii ≠ 0

Утверждение. Пусть все главные угловые миноры A отличны от 0:

A1 = a11 ≠ 0

;

A2 = ( a11 a12 ) ≠ 0 a21 a22

;

…;

Ai = ( a11 ... a1i ) ≠ 0 ... ai1 ... aii

Тогда разложение матрицы (2) существует и единственно.

Доказательство.

∃! A_i = B_i C_i = b₁₁ × b₂₂ × … × b_{i − 1,i − 1} × b_ii

b₁₁ × b₂₂ × … × b_{i − 1,i − 1} = A_i − 1

b_ii = ^A_i/_Ai − 1, i = , i = 2…m

ч. т. д.

Связь метода Гаусса с разложением A = B × C

В методе Гаусса для того, чтобы свести матрицу к верхнетреугольной матрице с единицами на главной диагонали, требуется ^(m3 − m)/₃ операций.

Задача. Показать, что для реализации (вычисления) по формулам (3), (4) требуется точно такое же число действий: ^(m3 − m)/₃ (решение).

В результате факторизации получаем систему BСx = f. Система Ax = f (1) сводится к решению двух уравнений By = f (5), Cx =  y (6)

Число действий для решения (5):
b_i1y₁ + b_i2y₂ + … + b_iiy_i = f_i, i = 1..m
y_i = ^{(f_i - Σ_{l = 1}i − 1b_ily_l)}/_bii
i − 1 умножение + 1 деление — i действий
Итого Σ_{i = 1}^mi = ^m(m+1)/₂ действий.
В точности то же самое число действий, что требуется для вычисления правой части в методе Гаусса.

Число действий для решения (6):
x_i + C_{i, i + 1}x_i+1 + … + C_imx_m = y_i, i = 1..m
x_i = y_i − Σ_{l = i + 1}^mc_ilx_l
m − i умножений, i = 2…m
Итого ^m(m − 1)/₂
В точности то число действий, которое требуется для обратного хода Гаусса.

Весь метод Гаусса требует ^(m3 − m)/₃ + m² = ^m/₃ × (m² + 3m − 1). В чём же выигрыш? Он будет показан при нахождении обратной матрицы.

Параграф 3. Обращение матрицы методом Гаусса-Жордана

Пусть есть матрица A размера m × m, det A ≠ 0. Тогда для неё существует A⁻¹: A⁻¹ × A = A × A⁻¹ = E.

Один из способов нахождения матрицы — нахождение алгебраических дополнений. Один из недостатков в том, что при малом определителе гигантские погрешности. Второй способ в приписывании единичной матрицы и проведения над ней действий трёх видов до получения единичной матрицы слева. Тоже не лучше.

Метод Гаусса-Жордана:

• Обозначим A⁻¹ = X. Тогда AX = E
• Запишем в явном виде:
• ∑_{l = 1}^ma_ilx_lj = δ_ij, i, j = 1..m (2)
• m² неизвестных, поэтому действий m⁶

Сейчас мы покажем, что разумно организовав алгоритм, получим порядка m³ действий.

Решение этой системы распадётся на решение m систем, у которых будет меняться только правая часть.

Введём вектор x^(j) = (x_1j, x_2j, …, x_mj)^T

Вектор δ^(j) = (0, 0, …, 0, 1 (j-я позиция), 0, …, 0)

Тогда Ax^(j) = δ^(j), j = 1..m (3)

Факторизуем A = B×C. Для этого потребуется ^(m3 − m)/₃ действий

Далее решаем 2 системы:

BCx^(j) = δ^(j)

Cx^(j) = y^(j)

• By^(j) = δ^(j) (4), j = 2…m
• Cx^(j) = y^(j) (5)

Решение пары таких систем требует m² действий, а всего таких систем m. То есть, итого ^(m3 − m)/₃ + m³ = ⁴/₃ m³ − ^m/₃

Теперь покажем, что число действий можно свести к m³. За счёт чего мы можем снизить число действий:

В уравнении (4) правые части имеют очень специальный вид. Мы в первый раз ненулевую компоненту встретим только на j-м шаге.

Первое уравнение:

• b₁₁y₁^(j) = 0 => y₁^(j) = 0
• b₂₁y₁^(j) + b₂₂y₂^(j) = 0 ⇒ y₂^(j) = 0
• …
• y_j − 1^(j) = 0 (6*)
• b_jjy_j^(j) = 1 ⇒ y_j^(j) = ¹/_bjj
• b_ijy_j^(j) + b_{i, j + 1}y_{j + 1}^(j) + … +b_iiy_i^(j) = 0, i = j + 1..m
• y_i^(j) = ^{(Σ_{l = j}i − 1b_ily_l(j))}/_bii, i = j+1..m (6)

Пусть i, j фиксированы. Тогда 1 деление и i − j умножений. Отпустим один из индексов, i. Получим:

• m − j + (m − j − 1) + … + 2 + 1 = ^{(m − j + 1)(m − j)}/₂ умножений
• m − j делений (6) + 1 (6*)
• m − j + 1 делений

Всего при фиксированном j: ^{(m − j + 1)(m − j)}/₂ + ^{2(m − j + 1)}/₂ = ^{(m − j + 1)(m − j + 2)}/₂ умножений и делений.

Отпускаем j, j = 2…m. Получим: ∑_j = 1^m(^{(m − j + 1)(m − j + 2)}/₂) ... (δ)

Задача. Показать, что (δ) = m(m+1)(m+2)/6

Для (5) мы выигрыша получить не можем. Поэтому для (5) при фиксированном j получаем ^m(m − 1)/₂ действий. Для решения всех систем (5) требуется ^m2(m− 1)/₂ действий. Итого получаем:

^(m3 − m)/₃ + ^{m(m + 1)(m + 2)}/₆ + ^{m2(m − 1)}/₂ = ^m/₆ × (2m² − 2 + m² + 3m + 2 + 3m² − 3m) = m³

Таким образом, разумно организовав вычисления, мы получим обратную матрицу за m³ операций.

Параграф 4. Метод квадратного корня

Другое название: метод Холецкого

Решаем уравнение Ax = f (1)

• |A| ≠ 0 (m × m)
• A = A* ⇔ a_ij = a_ji
• A = A^T — для вещественных матриц

Суть метода:

Факторизуем матрицу A в виде A = S*DS (2)

	s11	s12	s13	...	s1m
	0	s22	s23	...	s2m
S =	0	0	s33	..	s3m
	...
	0	0	0	...	smm

s_ii > 0, i = 2…m

D = diag (d₁₁, …, d_mm) d_ii = ±1, i = 1..m

Тогда система (1) решается следующим образом:

• S*DSx = f
• S*Dy = f (3)
• Sx = y

Эти две системы решаются явным образом, ибо матрицы треугольные

рассмотрим вещественную матрицу A = ((a₁₁, a₁₂), (a₁₂, a₂₂)) a_ij — вещественное. Тогда S = ((s₁₁, s₁₂), (0, s₂₂)), S* = ((s₁₁, 0), (s₁₂, s₂₂)), D=((d₁₁, 0), (0, d₂₂))

Умножим DS = ((d₁₁, 0), (0, d₂₂)) ((s₁₁, s₁₂), (0, s₂₂)) = ((d₁₁s₁₁, d₁₁s₁₂), (0, d₂₂s₂₂)). Она сохранила верхнетреугольную форму (позже покажем, что умножение верхнетреугольной на почти верхнетреугольную сохраняет форму). Домножим на S*. S*DS = ((s₁₁, 0), (s₁₂, s₂₂))((d₁₁s₁₁, d₁₁s₁₂), (0, d₂₂s₂₂)) = ((d₁₁s₁₁², d₁₁s₁₂s₁₁), (d₁₁s₁₂s₁₁, d₁₁s₁₂² + d₂₂s₂₂²)) = A. Поэлементно приравниваем компоненты.

• d₁₁s₁₁² = a₁₁
• d₁₁s₁₂s₁₁ = a₁₂
• d₁₁s₁₂² + d₂₂s₂₂² = a₂₂

Посмотрим на первое уравнение. Из него понятно, что d₁₁ = sign a₁₁. Тогда отсюда получаем s₁₁ = sqrt(|a₁₁|). Дальше смотрим на второе уравнение: s₁₂ = a₁₂/d₁₁s₁₁. Теперь последнее уравнение рассматриваем. Перепишем его в виде: d₂₂s₂₂² = a₂₂ − d₁₁s₁₂². Отсюда получим d₂₂ = sign(a₂₂ − d₁₁s₁₂²), s₂₂ = sqrt(|a₂₂ − d₁₁s₁₂²|). Таким образом, мы показали, что для матрицы 2×2 в таком виде осуществима.