Задачі лінійної алгебри.
Розв’язування системи лінійних алгебраїчних рівнянь (ЛАР) типу
А · х = b (1)
– матриця коефіцієнтів при невідомих;
– вектор-стовпчик невідомих;
– вектор-стовпчик вільних членів.
В розгорнутому вигляді цю систему можна записати так:
a11х1 + а12х2 + … + а1nxn = b1
a21х1 + а22х2 + … + а2nxn = b2 (2)
… … … … … … … … … … …
an1х1 + аn2х2 + … + аnnxn = bn
Розв’язком системи рівнянь є вектор
Х = (х1, х2, …, хn)T ,
елементами котрого є корені рівняння системи (2), що виражаються у вигляді лінійних функцій.
Методи розв’язування систем
лінійних алгебраїчних рівнянь.
А · х = b
Поділяються на дві групи : прямі та ітераційні.
1) Прямі методи – вони зводяться до кінцевих алгоритмів для обчислення коренів рівнянь (тобто розв’язки шукають за певними формулами). Вони дають розв’язки після виконання відомого для даного n (n – порядок) числа арифметичних операцій.
Іншими словами, прямим методом розв’язування лінійної системи А·x = b називають будь-який метод, котрий дозволяє знайти елементи вектора x з допомогою скінченого числа елементарних математичних операцій: додавання, віднімання, ділення, множення, та, можливо, кореня квадратного.
Оцінити ефективність будь-якого методу можна за допомогою двох – трьох важливих характеристик:
числа операцій, необхідних для реалізації даного методу;
об’єму пам’яті;
чутливості до переносу похибок заокруглення (або обчислювальної стійкості).
Розглянемо поняття стійкості методу. Нехай за значенням вхідної величини х шукається значення вихідної величини y. Якщо х має абсолютну похибку Δх, то розв’язок має похибку Δу. Метод називається стійким за вхідним параметром х, якщо малі похибки у вхідному параметрі х викликають малі похибки і в розв’язку у. Відсутність стійкості (нестійкий метод) означає, що навіть незначні похибки у вхідних даних ведуть до значних похибок в розв’язку, або до зовсім неправильного результату.
Практично всі прямі методи розв’язування систем ЛАР базуються на зведені матриці А до матриці більш простішої структури – діагональної (тоді розв’язок очевидний) або трикутної – та методів розв’язування таких систем.
До групи прямих методів належать:
– метод Гауса та його різновиди:
а) класичний метод Гауса із зведенням матриці А до верхньої трикутної матриці і одержанням розв’язків з допомогою обернених підстановок. Число операцій (вартість методу) – n2/2 операцій сумування, множення та n операцій ділення (можна ними знехтувати в порівнянні з n2/2)
б) метод Гауса з вибором головного елемента (частковим або повним). Число арифметичних операцій при цьому складає ~ n3/3 сумувань та ~ n3/3 множень.
Саме цим визначається повна вартість методу, оскільки вартість розв’язку вже самої трикутної системи незначна в порівнянні з вартістю зведення матриці до трикутного вигляду.
– LU-розклад (lower-upper –нижній-верхній)
А · х = b A = LU Ax = LUx = b
Ly = b
Ux = y
Якщо використовувати алгоритм Краута, то число операцій складе .
З точки зору об’єму обчислень методу LU-розкладу еквівалентний методу Гауса з частковим вибором головного елемента його переваги – це можливість роботи з різними векторами вільних членів b та з транспонованими матрицями АТ (рівняння АТх = b – розв’язок знаходиться за тим же LU-розкладом).
– метод Халецького (схема).
При розкладі симетричних матриць можна зменшити число операцій і об’єм пам’яті. Повна вартість складає половині вартості методу Гауса + n обчислень квадратного кореня. Метод чисельно стійкий.
– метод Жордана (роблять діагональну матрицю замість трикутної). Досить рідко використовується на практиці.
До прямих методів відносяться також методи для кліткових та розріджених матриць.
2) Ітераційні (або наближені) методи – це методи послідовних наближень. В них необхідно задати деякий наближений розв’язок – так зване початкове наближення. Після цього з допомогою деякого алгоритму проводиться один цикл обчислень, котрий називається ітерацією. В результаті ітерації знаходять нове наближення. Ітерації проводять до тих пір, доки не одержать розв’язок із заданою похибкою. Об’єм обчислень при цьому наперед не відомий.
При використанні ітераційних методів виникає проблема збіжності чисельного методу. Вона (збіжність) означає близькість одержаного після проведення ітерації розв’язку до істинного розв’язку. Збіжність ітераційного процесу: якщо в результаті проведення ітерацій одержуємо деяку послідовність х1, х2, … , хn (не важливо, це скалярні чи векторні величини), та якщо ця послідовність збігається до точного розв’язку х = а, тобто існує границя цієї послідовності , то метод є збіжним.
Чому виникла потреба в ітераційних методах? Чим не влаштовують прямі методи?
Основний недолік прямих методів – це нагромадження похибок в процесі розв’язування, оскільки обчислення на будь-якому етапі використовують результати (з похибками) попередніх операцій. Це особливо небезпечно для великих систем (n = 106 і більше) – наростає число операцій, а також для погано обумовлених систем (det A ≈ 0 ) (малі похибки обчислень або вхідних даних можуть породити значні похибки в розв’язку). Тому прямі методи використовують для відносно невеликих ( n < 200 ) систем з густо заповненою матрицею та det A ≠ 0 .
Перевагою ітераційних методів над прямими є те, що окремі похибки, що виникають при обчисленнях, не впливають на кінцевий результат (ітерації закінчуються тільки тоді, коли одержано розв’язок із наперед заданою точністю), а ведуть лише до збільшення числа необхідних ітерацій.
Крім того, розв’язування систем ітераційними методами спрощується ще й тому, що на кожній ітерації розв’язується система з одними і тими ж матрицями.
До ітераційних належить : метод простої ітерації, метод Зейделя, метод верхньої релаксації, та інші.
Прямі методи розв’язування систем лінійних алгебраїчних рівнянь
Класичний метод Гауса.
Розглянемо систему рівнянь четвертого порядку:
А · х = b
(1)
Зауважимо, що елементи вектора-стовпчика вільних членів b занесені в матрицю коефіцієнтів А.
Будемо вважати, що а11 ≠ 0. З першого рівняння знаходимо х1:
х1 = α12х2 + α13х3 + α14х4 + α15 , (2)
де , .
З допомогою рівняння (2) можна виключити х1 з решти рівнянь, для чого достатньо підставити (2) для х1 в друге, третє і четверте рівняння системи. Це і є першим кроком – кроком виключення невідомого х1.
x1 = α12x2 + α13x3 + α14x4 + α15
,
Перехід від початкової системи
до новоствореної
відбувається за такою формулою:
Другий крок – виключення невідомого х2 відбувається аналогічно:
x2 = α23x3 + α24x4 + α25
Третій крок – виключення невідомого х3
х3 = α 34х4 + α 35
,
;
Останнє рівняння можна переписати у вигляді:
або .
Отже, в результаті прямого ходу одержимо систему рівнянь:
х1 = α12х2 + α13х3 + α14х4 + α15 ;
х2 = α23х3 + α24х4 + α25 ;
х3 = α34х4 + α35 ;
х4 = α45 .
Знаходження невідомих проводиться в оберненому ході методу Гауса шляхом зворотніх підстановок.
Якщо п – кількість рівнянь (порядок) системи, то програмування обчислювального процесу проводиться так:
L – кількість кроків виключення ;
j – позначення другого індексу при визначенні α ;
і – номер рядка системи ;
k – номер стовпця.
Можна записати, що для всіх
Обернений хід: , .
Отже, обчислювальна схема прямого ходу методу Гауса
Для L = 1, п – 1
Для j = L + 1, n + 1
alL,j = – aL,j/aLL
Для i = L + 1, n
Для K = L + 1, n + 1
aiK = aiK + aiL∙alLK
піддається спрощенню.
Однак початкове обчислення всіх коефіцієнтів α не є обов’язковим. Це випливає з наступного. Наприклад, перехід від початкової системи коефіцієнтів до наступної відбувається так:
Наприклад, коефіцієнти першого чи другого стовпця нової системи утворюються за правилом
, або
Отже, визначивши, наприклад α12 зразу ж можна переходити до визначення коефіцієнтів нової системи і т.п. Таким чином цикли по J i по K можна об’єднати (оскільки, що J i K змінюються в однакових межах).
Якщо α(L) замінити на а(L) та цикли по J та по K об’єднати в один (тобто J на K), то одержимо загальну форму методу виключення Гауса із стовпцевою формою розкладу матриці А до трикутного вигляду)
В кінці цих перетворень одержимо:
Обчислювальна схема прямого ходу зображується так:
L = 1, n – 1
K = L + 1, n + 1
I = L +1, n
Метод Гауса (рядкова і стовпцева форми
розкладу матриці А до трикутного вигляду).
Розглянемо систему лінійних рівнянь
А · х = b (1)
з невиродженою матрицею А розміру n*n .
Більшу частину всього обчислювального процесу поглинає зведення матриці А до трикутного вигляду.
Можливі дві форми розкладу (зведення) матриці А до трикутного вигляду – рядкова або стовпцева.
Як ми вже знаємо, стовпцева форма розкладу зображується наступною схемою:
Для l = 1 до n – 1
Для k = l + 1 до n
alk = – alk/all (2)
Для і = l + 1 до n
аik = aik + ail*alk
Права частина b системи (1) також може оброблятися в ході зведення А до трикутного вигляду. Тому можна bі приєднати до і-го рядка (член аі,n+1) – як ми й робили раніше (в такому випадку в циклі " k " верхня межа зростає до n + 1). Можна bі залишити на місці не вносячи цього в масив А. В цьому випадку в результаті виконання (якщо не вносити) прямого ходу методу Гауса одержується система рівнянь:
(3)
... ... ... ... ... ... ... ... ... ... ... ... ... ... ... ... ...
.
Стовпчикову форму розкладу можна продемонструвати наступним прикладом:
Перехід, наприклад, від матриці А початкових коефіцієнтів при невідомих при n=4:
(4)
до матриці коефіцієнтів після першого етапу перетворень
(5)
відбувається наступним чином: обчислюються послідовно
;
формуються елементи стовпця
; ; .
Обернений хід при стовпчиковій формі розкладу описується загальною формулою:
, (6)
Розглянемо тепер рядкову форму розкладу матриці А. Вона базується на зведенні системи лінійних рівнянь до трикутного вигляду. Для цього спочатку нормують перше рівняння, ділячи його на а11(0), тобто роблять коефіцієнт при х1 рівним 1. Потім це перше рівняння помножують відповідно на коефіцієнт аі,1(0) при х1 всіх інших рівнянь і послідовно віднімають від усієї решти рівнянь. В результаті х1 буде виключене із всіх рівнянь, крім першого. На другому кроці виключають х2 з третього, четвертого, ..., п –го рівнянь. Цю процедуру повторюють до тих пір, доки вся система не буде зведена до такого трикутного вигляду
Рядкова форма зображається наступною обчислювальною схемою:
Для l = 1 до n – 1
Для і = l + 1 до n
аіl = – aіl/all (7)
Для k = l + 1 до n
аik = aik + ail*alk
Тут перехід від матриці (4) до матриці (5) відбувається так:
Тобто, на відміну від стовпчикової форми, обчислення коефіцієнтів нової матриці відбувається по рядках. Результат же одержується той самий.
Стовпчикова форма:
Рядкова форма:
При (обертанні) обчисленні оберненої матриці доцільно використовувати розклад матриці А до трикутного вигляду за рядковою формою.
Метод Гауса з вибором головного елемента.
При прямому ході потрібно виконувати ділення на . Якщо в системі ЛАР вважати, що аіі = 0, то в цьому випадку метод Гауса використовувати не можна.
Елемент, на який потрібно ділити, називається головним (провідним) або центральним і ділення на нульовий провідний елемент неможливе.
Велики похибки заокруглень також виникають тоді, коли головний елемент не дорівнює нулю, але є малим числом. Наприклад, для системи з трьох рівнянь з трьома невідомими:
1.2357х1 + 2.1742х2 – 5.4834х3 = – 2.0735
6.0696х1 – 6.2163х2 – 4.6921х3 = – 4.8388
3.4873х1 + 6.1365х2 – 4.7483х3 = 4.8755
метод Гауса дає наступні результати:
х1 = 0.99968 х2 = 0.99994 х3 = 0.99991
(5 значущих цифр), що добре узгоджується з точним розв’язком х1 = х2 = х3 =1.
Тепер поміняємо місцями друге і третє рівняння системи і знову розв’яжемо її методом Гауса. В цьому випадку розв’язок буде мати наступний вигляд:
х1 = 2.9021 х2 = 1.4286 х3 = 0.9992
Причина такої відмінності полягає в тому, що в першому випадку після першого етапу провідним елементом є , а в другому (після перестановки) – , що і дає таку розбіжність (при обмеженій розрядності чисел (точності обчислень)).
Для усунення такого роду неприємностей проводять вибір головного елементу таким чином – за головний вибирають елемент з максимальним абсолютним значенням. Причому використовують різні варіанти вибору – по рядку, по стовпцю або по всій матриці. Якщо вибір проводиться по стовпцю (наприклад, при х1) – то в цьому випадку проводиться так званий частковий вибір провідного елемента. Якщо вибирати максимальний елемент з усієї матриці, то такий вибір головного елемента називають повним.
Метод прогону .
Більшість технічних задач зводиться до розв’язування систем ЛАР, в яких матриці містять багато нульових елементів, а ненульові елементи розміщені за спеціальною структурою (стрічкові квазітрикутні матриці).
Задачі побудови інтерполяційних сплайнів, різницевих методів розв’язування крайових задач для диференціальних рівнянь зводяться до розв’язування систем ЛАР з трьохдіагональною матрицею А. В матриці А всі елементи, що не лежать на головній діагоналі і двох сусідніх паралельних діагоналях, дорівнюють нулю.
В загальному вигляді такі системи записують так:
аіхі-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 (враховуючи, що Сп = 0)
Долі використовуючи формулу (3) знаходимо послідовно всі невідомі xn-1, xn-2, … , x1.
Майже у всіх задачах, що приводять до розв’язку системи (2) з трьохдіаго-нальною матрицею, забезпечується умова переважання діагональних коефіцієнтів
Це забезпечує існування єдиного розв’язку та достатню стійкість методу прогону відносно похибок заокруглення.
Для запису коефіцієнтів ai, bi, та прогоночних коефіцієнтів Аi-1, Ві-1 використати один і той же масив.
Приклад.
Розв’язати методом прогонки систему:
4х1 + х2 = 5,6
х1 + 4х2 + х3 = 7,2 х2 + 4х3 +х4 = 7,8
х3 + 4х4 + х5 = 8,4
х4 + 4х5 = 7,4
Розв’язок: х1 =1,1 ; х2 = 1,2; х3 = 1,3; х4 = 1,4; х5 = 1,5
Обчислення визначника квадратної матриці.
Визначником n-го порядку із n2 елементів квадратної матриці А називається алгебраїчна сума всіх можливих добутків n елементів, взятих по одному із кожного рядка та кожного стовпця матриці.
Наприклад,
(1)
Звідки взялися знаки “ – ”?
Спочатку зупинимось на деяких загальних поняттях.
Перестановкою чисел 1, 2, …, n називають будь-яке розташування цих чисел в певному порядку. Число всіх перестановок 1∙ 2∙ …∙ n = n! В одній із перестановок числа розташовані у вигляді натурального ряду. В інших перестановках цей процес порушується.
Розрізняють в розташуванні чисел:
а) інверсію – коли більше число стоїть попереду меншого;
б) порядок – якщо менше число стоїть попереду більшого.
Перестановки можуть бути парні і непарні – з парним та непарним числом інверсій.
Тепер щодо добутків в det A. Знак окремого добутку буде (–1)r , де r – число інверсій, що визначається за другим індексом елементів матриці, при розташуванні перших індексів в порядку зростання (а саме так розташовані перші індекси у всіх добутках (1)). Тому, (1) можна записати у початковому вигляді так:
де r1 = 0; r2 = 2; r3 = 2; r4 = 3; r5 = 1; r6 = 1.
Звідси й випливає формула (1).
Визначник матриці не змінюється при транспонуванні матриці, тобто det AT = det A. При взаємній перестановці двох стовпців або двох рядків визначник не змінює свого абсолютного значення, але змінює знак на протилежний. Якщо в квадратній матриці є два однакових стовпці або два однакових рядки, то det A = 0. Визначник знаходиться тільки для квадратних матриць. Якщо det A ≠ 0, то матриця А називається невиродженою (регулярною, неособливою).
Особлива (вироджена, сингулярна) матриця, для якої det A = 0, оберненої матриці А–1 немає.
Система нормальних рівнянь
Ах = b
має розв’язок, коли матрицю А можна обернути, тобто для неї існує обенена матриця А–1.
Матриця А комутативна з оберненою матрицею, тобто
А–1А = А А–1 = Е
Обчислення матриці А–1 – доволі важкий процес.
Обчислення визначника матриці методом Гауса.
Метод Гауса дає можливість обчислити визначник матриці А системи лінійних рівнянь
Ах = b (1)
Прямий хід метода Гауса не змінює значення визначника матриці А. Якщо в прямому ході рядки матриці не переставляються, то знак визначника не змінюється. Визначник матриці дорівнює добутку провідних (діагональних) елементів в прямому ході методу Гауса
(2)
Для його обчислення здійснюється прямий хід методу Гауса, як і при розв’язку системи, тільки без перетворень вектора b. Після завершення прямого ходу знаходиться детермінант матриці А.
Однак, якщо використовувати описаний метод для обчислення визначника
наступної матриці
,
то це призведе до зупинки з індикацією помилки (ділення на нуль), тому користуються програмою з вибором головного елемента.
Якщо використовують метод Гауса з вибором головного елемента, то в формулу необхідно додати множник (–1)К
,
де К – кількість перестановок рядків.
Обчислення оберненої матриці.
Розглянемо невироджену квадратну матрицю А (det A ≠ 0) і знайдемо елементи її оберненої матиці А–1:
А–1={ αij } , ; ,
використовуючи співвідношення:
А ∙ А–1 = Е , (1)
де Е – одинична матриця.
Співвідношення (1) можна записати в розгорнутому вигляді:
; (2)
– символ Кронекера
В співвідношенні (2) j – зафіксоване, тому j-й стовпчик оберненої матриці αjK можна розглядати як вектор невідомих. Отже, при фіксованому j розв’язок системи (2) буде j-м вектором-стовпчиком матриці А–1.
Покажемо це на конкретному прикладі: n = 4
A ∙ А–1 = Е
Елементи оберненої матриці можна знайти, розв’язуючи n (в даному випадку n = 4) систем алгебраїчних рівнянь:
j = 1
j = 2
j = 3
j = 4
Маємо n систем алгебраїчних рівнянь з однаковою матрицею А та різними стовпцями вільних членів.
Якщо для розв’язування лінійних систем використовувати метод Гауса, то зведення матриці А до трикутного вигляду виконується тільки один раз, а обернений хід виконується n разів (кожний раз для відповідного вектора-стовпчика вільних членів).
Для того, щоб перевірити, чи вірно знайдені елементи оберненої матриці, необхідно перемножити матрицю А на її обернену матрицю А–1, в результаті чого одержимо одиничну матрицю.