Метод прогону
Більшість технічних задач зводиться до розв’язування систем лінійних алгебраїчних рівнянь, в яких матриці містять багато нульових елементів, а ненульові елементи розміщені за спеціальною структурою (стрічкові квазітрикутні матриці).
Задачі побудови інтерполяційних сплайнів, різницевих методів розв’язування крайових задач для диференціальних рівнянь зводяться до розв’язування систем лінійних алгебраїчних рівнянь з трьохдіагональною матрицею А. В матриці А всі елементи, що не лежать на головній діагоналі і двох сусідніх паралельних діагоналях, дорівнюють нулю.
В загальному вигляді такі системи записують так:
аіхі-1 + bixi + cixi+1 = di (1)
1 ≤ i ≤ n ; a1 = 0 ; cn = 0
або в розгорнутому вигляді:
b1x1 + c1x2 = d1
a2x1 + b2x2 + c2x3 = d2
a3x2 + b3x3 + c3x4 = d3
∙∙∙ ∙∙∙ ∙∙∙ ∙∙∙ ∙∙∙ ∙∙∙ ∙∙∙ ∙∙∙ ∙∙∙ ∙∙∙ ∙∙∙ ∙∙∙ ∙∙∙ ∙∙∙ ∙∙∙ ∙∙∙ ∙∙∙ ∙∙∙ ∙∙∙ ∙∙∙ (2)
aixi-1 + bixi + cixi+1 = di
an-1xn-2 + bn-1xn-1 + cn-1xn = dn-1
anxn-1 + bnxn = dn
Вибір найбільшого елемента при виключенні невідомих за методом Гауса в таких системах робити не можна, оскільки перестановка рядків руйнує структуру матриці. Найчастіше для розв’язку системи з трьохдіагональною матрицею використовують метод прогону, який є частковим випадком методу Гауса.
Прямий хід прогону (алгоритм прямого ходу методу Гауса).
Кожне невідоме хі виражається через хі+1 з допомогою прогоночних коефіцієнтів Аі та Ві
хі = Аіхі+1 + Ві ; (3)
Наприклад, з першого рівняння системи (2) знайдемо
, звідки (4)
З другого рівняння системи (3) виразимо х2 через х3, замінюючи х1 формулою (3) або (4)
a2х1 + b2x2 + c2x3 = a2(A1x2 + B1) + b2x2 = d2
Звідси знайдемо
або х2 = А2х3 + В2
, беручи за е2 = а2А1 + b2
Запишемо:
Аналогічно для кожного і прогоночні коефіцієнти з рівняння хі = Аіхі+1 + Ві мають вигляд:
(5)
еі = аіАі-1 + bі ;
При цьому враховуючи, що а1 = сn = 0 приймаємо
Ао = 0 ; Во = 0. (6)
В розгорнутому вигляді формула (5) буде мати вигляд формули (7). Значення прогоночних коефіцієнтів можна одержати і таким шляхом. В рівнянні (3) понизимо індекс на одиницю та підставимо значення хі-1 в і-е рівняння системи (1)
хі-1 = Аі-1∙хі + Ві-1
аіхі-1 + bixi + cixi+1 = di ai(Ai-1xi + Bi-1) + bixi + cixi+1 = di
(7)
Обернений хід прогонки (аналог оберненого ходу методу Гауса).
Він полягає в послідовному обчисленні невідомих хі. Спочатку знаходять хn. Для цього формулу (7) запишемо при і = n (враховуючи, що Cn =0)
Долі використовуючи формулу (3) знаходимо послідовно всі невідомі xn-1, xn-2, … , x1.
Майже у всіх задачах, що приводять до розв’язку системи (2) з трьохдіаго-нальною матрицею, забезпечується умова переважання діагональних коефіцієнтів
Це забезпечує існування єдиного розв’язку та достатню стійкість методу прогону відносно похибок заокруглення.
Для запису коефіцієнтів ai, bi, та прогоночних коефіцієнтів Аi-1, Ві-1 використати один і той же масив.
Метод LU – розкладу
(метод розкладу на трикутні матриці
або метод LU – факторизації )
Алгоритми цього методу досить близькі до методу Гауса. Основна перевага методу LU – факторизації в порівнянні з методом Гауса є можливість більш простого одержання розв’язку для різних векторів b системи лінійних алгебраїчних рівнянь
A · x = b, (1)
A – матриця розміру nn із сталими коефіцієнтами;
b – n-мірний вектор вільних членів;
х – n-мірний вектор невідомих.
Матриця системи А представляється у вигляді добутку двох трикутних матриць L та U
А =L · U, де
,
L – нижня трикутна матриця;
U – верхня трикутна матриця.
Зауважимо, що на головній діагоналі матриці U стоять одиниці. Це в свою чергу означає, що визначник матриці А дорівнює добутку діагональних елементів lii матриці L. (det U = 1; det L = l11l22…lnn; det A = det L).
Отже, систему (1) можна записати у вигляді:
L · U · x = b. (2)
Введемо допоміжний вектор у як
U · x = у (3)
Тоді можна записати, що
L · у = b. (4)
Розв’язування системи (1) або (2) відбувається в два етапи: спочатку розв’язується система L · у = b, а потім U · x = у. Така послідовність розв’язку дає перевагу методу (в порівнянні з методом Гауса) – якщо потрібно розв’язати декілька систем з однією і тією ж матрицею А, задача суттєво спрощується, оскільки зберігаються матриці L та U.
Розв’язок системи L · у = b називається прямим ходом, а системи U · x = у – оберненим ходом.
Розглянемо прямий хід. Завдяки спеціальній формі матриці L вектор у можна легко визначити. Для цього (4) перепишемо у вигляді системи рівнянь
Звідси можна одержати, що в загальному вигляді
(5)
для
При оберненому ході вектор х визначається з системи рівнянь (3)
x1 + u12x2 + u13x3 + … + u1nxn = y1
x2 + u23x3 + … + u2nxn = y2
x3 + … + u3nxn = y3
………………………………
xn = yn.
Починаючи з останнього рівняння, можна послідовно знайти компоненти вектора х. В загальному вигляді вони визначаються за формулами
для (6)
Тепер розглянемо LU – розклад матриці А.
У варіанті методу, що називається методом Крута, використовується наступна послідовність знаходження елементів матриць L та U () – де k – крок.
Для всіх
де (7)
де (8)
Домовимось, як звичайно, вважати значення суми рівним нулю, якщо верхня границя (межа) сумування менша від нижньої.
Алгоритм
1) k =1
Обчислюємо
для всіх
для всіх .
2) k = k + 1
де
де .
Продовжуємо до тих пір, доки k = n.
Зауваження щодо алгоритму побудови програми:
Необхідні в цих співвідношеннях значення елементів матриць L та U обчислюються на попередніх етапах процесу.
Кожний елемент aij матриці А потрібен тільки для обчислення відповідних елементів матриці L та U ( тобто, в подальшому процесі aij не потрібні). Оскільки нульові елементи матриць L та U, а також одиничну діагональ матриці U запам’ятовувати не потрібно, тобто в процесі обчислення матриці L та U можуть бути записані на місце матриці А. Причому матриця L розташована в нижньому трикутнику (i ≥ j), U – відповідно в верхньому трикутнику (i < j) матриці А.
Розклад матриці А на матриці L та U , як правило, об’єднують з прямим ходом. Може бути використана наступна послідовність обчислень: спочатку обчислюється перший стовпець матриці L , потім перший рядок матриці U та перший елемент вектора у; далі – другий стовпець матриці L , другий рядок матриці U та другий елемент вектора у і т.д.