Метод Релаксації для розв'язку Слар

Вступ

 

Чисельний розв’язок системи  лінійних алгебраїчних рівнянь (СЛАР) - одне з найбільш часто зустрічаються  завдань у науково-технічних дослідженнях. Така задача виникає в математичній фізиці (чисельне рішення диференціальних  та інтегральних рівнянь), економіці, статистиці тощо. При цьому прикладні задачі часто потребують вирішення великих і надвеликих СЛАР з числом невідомих більше 1000. До таких СЛАР, наприклад, звертаються при чисельному розв’язанні двовимірних і особливо тривимірних задач математичної фізики, в яких умови фізичної і геометричної апроксимації двовимірної і тривимірної області потребують використання досить дрібної розрахункової сітки з великим числом розрахункових вузлів за лінійним розміром.

Існуючі бібліотеки програм  на мовах високого рівня, розроблені на основі, так званих, прямих методів розв'язування СЛАР, типу методу Гауса і його модифікацій. Число арифметичних операцій множення для чисельного розв'язання СЛАР розмірністю за допомогою прямого методу - . Кубічна залежність числа арифметичних операцій від розміру матриці СЛАР приводить при до нереально великого часу розв’язання навіть на найсучасніших ЕОМ. Крім того, час вирішення непропорційно зростає при використанні прямих методів у випадку з причини недостатності обсягу оперативної пам'яті для зберігання даних завдання.

Ітераційні  методи розв'язання СЛАР набагато економніші, як за часом розв’язання, так і  щодо використання оперативної пам'яті. Так, якщо ітераційний метод швидко сходиться з числом ітерацій , то час розв’язання, пропорційне вже квадрату розміру матриці ~ , виявляється істотно менше, приблизно в раз для дійсної і раз для комплексної СЛАР. Крім того, потрібно зберігати в оперативній пам'яті, як правило, тільки одну матрицю, наприклад, матрицю переходу явного ітераційного методу. При використанні ітераційних методів, що швидко збігаються, можна розв’язати в реальному часі на сучасних ЕОМ СЛАР з комплексною матрицею розмірністю .

В даний час  відсутні бібліотеки підпрограм широкого призначення для чисельного розв’язання  великих і надвеликих СЛАР. Таким  чином, розробка ефективних ітераційних алгоритмів для широкого класу матриць СЛАР великої розмірності і бібліотек підпрограм на їх основі є актуальною задачею.

Найбільш алгоритмічно простими серед ітераційних методів  є стаціонарні ітераційні методи, такі як оптимальний метод простої ітерації та метод релаксації. У той же час можна досягти їх ефективної збіжності для досить широкого класу дійсних і комплексних матриць СЛАР. Для нестаціонарних ітераційних методів, таких як метод з Чебишевським набором параметрів, мінімальних нев'язок, зв'язаних градієнтів, збіжність доведена у вузькому класі матриць, наприклад, таких як дійсні симетричні позитивні матриці. І в цьому вузькому класі матриць збіжність оптимальних стаціонарних методів, що спираються на відомі спектральні матричні властивості, виявляється в деяких випадках навіть кращою. При цьому число арифметичних операцій стаціонарного алгоритму мінімальне. Ще однією перевагою оптимального методу простої ітерації є можливість природного розпаралелювання алгоритму при постановці його на сучасних багато паралельних ЕОМ, так як алгоритм по суті зводиться до одного множення матриці на вектор. Всі ці аргументи вказують на вибір стаціонарних ітераційних методів в якості алгоритмічної основи для бібліотеки підпрограм з розв’язання СЛАР з великими матрицями.

У курсовій роботі розглянуто метод релаксації для рішення СЛАР.

 

Розділ 1. Теоретична основа методу релаксації

1.  Постановка задачі

 

Нехай задано систему  алгебраїчних рівнянь виду

    (1.1)

де  — невідомі, які потрібно визначити;

— коефіцієнти системи

— вільні члени системи.

(1.1) можна записати  у вигляді

     (1.2)

де

А - матриця розмірності ,

x=(x₁,x₂,...,x_n)^T- вектор розв’язків,

f=(f₁,f₂,...,f_n)^T - вектор правих частин.

Чисельні методи розв’язання даної системи прийнято розділяти на два класи: прямі  методи та ітераційні.

Прямими методами називаються  методи, що дозволяють отримати розв'язок системи рівнянь (1.1) за кінцеве число арифметичних операцій.

До прямих методів належать метод Крамера, метод Гаусса, LU-метод, метод прогонки і ряд інших  методів. Основним недоліком прямих методів є те, що для знаходження  розв’язку необхідно виконати велику кількість операцій.

Суть ітераційних методів полягає в тому, що розв’язок системи (1.1) знаходиться як границя послідовних наближень x⁽ⁿ⁾ при n®¥, де n - номер ітерації. Застосування ітераційних методів вимагає завдання початкового значення невідомих х⁽⁰⁾ і точності обчислень e> 0. Обчислення проводяться до тих пір, поки не буде виконана оцінка

     (1.3)

Основна перевага ітераційних методів полягає  в тому, що точність шуканого рішення  задається.

Число ітерацій n=n(e), яке необхідно виконати для отримання заданої точності e, є основною оцінкою якості методу. З цього числа проводиться порівняння різних методів.

Головним недоліком  цих методів є те, що питання  збіжності ітераційного процесу  вимагає окремого дослідження. Доведені в даний час теореми про  збіжність ітераційних методів мають місце для систем, на матриці яких накладені обмеження.

Прикладом звичайних  ітераційних методів можуть слугувати  метод Якобі (метод простих ітерацій), метод Зейделя, метод верхніх  релаксацій.

До особливого класу ітераційних методів слід віднести варіаційні ітераційні методи: метод мінімальних нев'язок, метод швидкого спуску і т.д.

Ітераційні  методи також діляться на однокрокові, коли для визначення розв’язку на j+1 ітерації використовуються значення розв’язку, знайдені на j ітерації, і багатокрокові, коли для визначення рішення на j+1 ітерації використовується кілька попередніх ітерацій.

Зауважимо, що існують системи, для яких ітераційний процес сходиться, а вектор нев'язки, що виходить при  підстановці знайденого розв’язку  у вихідну систему

,    (1.4)

виходить більшим  за величиною, тобто знайдений розв’язок  не задовольняє вихідної системи. У  цьому випадку в якості критерію досягнення точності розв’язку може бути взята величина нев'язки, яка  оцінюється за однією з норм .

Завданням роботи є знаходження  з (1.1) за допомогою методу релаксації.

2. Теоретичні основи розв’язання  СЛАР

 

Одним з найбільш поширених однокрокових ітераційних  методів є метод верхніх релаксацій, який має наступний вигляд

,    (1.5)

де  > 0 - заданий числовий параметр. Цей параметр вибирається таким чином, щоб на кожному кроці ітераційного процесу зменшувалася величина, що характеризує похибку отриманого розв’язку системи відносно до потрібно.

Для отримання  розрахункових формул (1.5) перепишемо у вигляді

  ,  (1.6)

або в покомпонентному  записі отримаємо

    (1.7)

Наведемо кілька рядків по компонентного запису

,    (1.8)

    (1.9)

   (1.10)

   і т.д.

Практика застосування ітераційних методів показала, що ці методи призводять до точного розв’язку  для систем з матрицею А, що має  спеціальний вид. Наведемо ряд теорем про збіжність ітераційних методів. Доведення цих теорем наводяться в книзі [1].

Розглянемо  ітераційні методи з постійним ітераційним  параметром, записані у вигляді

     (1.11)

Теорема 1. Нехай  А - симетрична позитивно визначена матриця, t > 0 і нехай виконано нерівність В-0,5tА>0. Тоді ітераційний метод (1.11) сходиться.

Наслідок 1. Нехай А - симетрична позитивно визначена матриця з діагональним переважанням, тобто

     (1.12)

Тоді метод  Якобі сходиться.

Наслідок 2. Нехай А - симетрична позитивно визначена матриця. Тоді метод верхніх релаксацій сходиться  за умови 0<w<2. Зокрема, метод Зейделя сходиться (w= 1).

Теорема 2. Ітераційний метод (1.11) сходиться при будь-якому  початковому наближенні тоді і тільки тоді, коли всі власні значення матриці   за модулем менше одиниці.

Теорема 3. Нехай А і  В - симетричні позитивно визначені  матриці, для яких справедливі нерівності , де g₁,g₂ - позитивні постійні, g₁>g₂. При ітераційний метод (1.11) збігається і для похибки справедливі оцінки

, i=0,1,...,   (1.13)

де

      (1.14)

      (1.15)

      (1.16)

       (1.17)

 

Наслідок 1. Якщо А^Т=А>0, то для методу простої ітерації

       (1.18)

при

   (1.20)

справедлива оцінка

,     (1.21)

де

      (1.22)

      (1.23)

 

Наслідок 2. Для симетричної  матриці А і 

     (1.24)

справедлива рівність

  ,     (1.25)

де, . Часто зустрічаються завдання з погано обумовленою матрицею А, коли відношення велике. У цьому випадку число r₀ близько до одиниці, і метод простої ітерації сходиться повільно.

Оцінимо число  ітерацій n₀(e), яке потрібно для досягнення заданої точності e у випадку малих x, тобто для отримання оцінки

       (1.26)

З умови  отримуємо, що

,      (1.27)

і при малих x маємо

     (1.28)

Зауважимо, що як критерій збіжності ітераційного методу може використовуватися нев'язка, яка виходить при підстановці знайденого розв’язку в систему (1.2). Отже, можемо сказати, що метод релаксації є актуальним для використання при розв’язанні СЛАР.

3. Метод верхніх релаксацій

Серед явних  однокрокових ітераційних методів  найбільше поширення отримав метод верхніх релаксацій (1.5). Це пов'язано з тим, що метод верхніх релаксацій містить вільний параметр w, змінюючи який можна отримувати різну швидкість збіжності ітераційного процесу.

Найбільш ефективно  цей метод застосовується при  вирішенні безлічі близьких алгебраїчних систем лінійних рівнянь. На першому етапі проводиться розв’язання однієї з систем з різними значеннями ітераційного параметра w і з аналізу швидкості збіжності ітераційного процесу вибирається оптимальне значення цього параметра. Потім всі інші системи розв’язуються з обраним значенням w.

Ще одна перевага ітераційного методу верхніх релаксацій полягає в тому, що при його реалізації на ЕОМ алгоритм обчислень має  простий вигляд і дозволяє використовувати  всього один масив для невідомого вектора.

Основна обчислювальна  формула має вигляд

     (1.29)

У вираз (1.29) і входять однаковим чином, отже, при обчисленнях вони можуть записуватися в один і той же масив. При реалізації методу верхніх релаксацій використовується наступна форма запису алгоритму обчислень

     (1.30)

Дійсно, при  послідовному знаходженні елемента (i+1 ітерації) на кожному кроці будуть використовуватися знайдені раніше значення, які при k<j відповідають i+1 ітерації, а при k>j - i ітерації.

Сучасна обчислювальна техніка  дозволяє проводити дослідження  стійкості та збіжності ітераційного методу в залежності від параметрів задачі. Наприклад, можна проводити  дослідження впливу підвищення точності розв’язання задачі на число необхідних ітерацій, дослідження впливу початкового наближення, зміни коефіцієнтів матриці А і правих частин системи.

4. Похибка методу верхньої релаксації

Одне з основних питань застосування ітераційних методів пов'язаний з коректністю вибору точності методу e.

Аналізуючи обчислювальні  похибки виразу (1.30), отримаємо оцінку найменшого значення точності методу верхніх релаксацій.

Очевидно, що шукана похибка  обчислень буде визначатися похибкою задання коефіцієнтів вихідної системи і похибкою округлення.

Запишемо різницю двох ітераційних наближень розв’язку  і оцінимо її мінімальне значення

     (1.31)

Нехай коефіцієнти  і f_i задані з деякою відносною похибкою. Припустимо, що ітераційний метод сходиться, і нев'язка

      (1.32)

буває зі зростанням номера ітерації k, тобто . Оцінка абсолютної похибки правої частини виразу (1.5) може бути представлена в наступному вигляді

,      (1.33)

тут - модуль мінімального значення діагонального елемента . Звідси випливає, що похибка методу .

5. Метод блочної релаксації

Вихідна матриця  розбивається на блоки. Вектор правої частини і вектор невідомих розбиваються на блок-вектори відповідної розмірності. Наприклад, для розміру блоку  рівного двом, отримуємо:

     (1.34)

де

     (1.35)

    (1.36)

     (1.37)

Запишемо формулу  для блоків матриці  і блок-векторів і :

   (1.38)

Позначимо

   (1.39)

    (1.40)

Тоді, підставляючи (1.39) і (1.40) в (1.38) і помноживши зліва  на , для кожного блок-вектора одержуємо СЛАР:

      (1.41)

Рішення отриманих  систем (1.41) рекомендується виконувати з використанням факторизації матриці  , причому факторизацію слід виконувати 1 раз перед першою ітерацією.

 

Розділ 2. Практична реалізація методу релаксації

1.  Приклад розв’язання системи методом релаксації

Приклад 1: Розв’язати наступну систему методом релаксації

(2.1)

Обчислення  проводити з точністю до двох знаків після коми.

Розв’язання:

Наводимо систему (2.1) до виду, зручного для рішення методом релаксації

(2.2)

Задаємо початкове  наближення нульовими значеннями

(2.3)

Знаходимо значення невязок

 

Детальніше  див. таблицю 2.1

 

Таблиця 2.1 Процес знаходження  невідомих і нев’язок

	0	0,60	0	0,70	0	0,80
		0,16		0,16		-0,80
		0,76		0,86		0
		0,17	0,86	-0,86		0,09
		0,93		0		0,09
	0,93	-0,93		0,09		0,09
		0		0,09	0,18	0,18
		0,04		0,04		-0,18
		0,04	0,13	0,13		0
		0,03		-0,13		0,01
	0,07	0,07		0		0,01
		-0,07		0,01		0,01
		0		0,01	0,02	0,02
		0		0		-0,02
		0
		0	0,01	0,01		0
		0		-0,01		0
		0		0		0
	1,00		1,00		1,00

 

Отже,

І так далі. Підставляємо результати обчислені в таблиці. Підрахувавши всі прирощення , що містить значення коренів,

Для перевірки  підставляємо знайдені значення коренів  у вихідне рівняння; в цілому система вирішена точно.

2. Блок-схема алгоритму

Рисунок 2.1 Блок-схема методу релаксації

 

 

 

3. Опис програмної реалізації методу релаксації

 

Для запуску  програми, що реалізує метод релаксації потрібно запустити файл Relax.exe. Програма Relax реалізована мовою програмування паскаль.

Програма бере дані із файлу 'a.txt', який розміщається в каталозі, де розміщена сама програма. 'a.txt' містить в першому рядку  число n, що характеризує розмірність матриці. І в наступних n рядкам містить по n+1 елементу, які відповідають коефіцієнтам матриці і вільним членам.

У програмі задана точність у змінній sigma:=1e-7, тобто 1*10^-7. Програма також обмежена максимальною кількістю кроків 10000.

Перелічимо  основні функції і процедури даної програми і їх вхідні параметри:

1. norma – знаходить норму, тобто визначає похибку. У вхідних параметрах приймає вектор правої частини системи типу TVector і повертає число – похибку типу double.
2. copyVector – робить копію вектора, що передається в параметрах. Вхідні параметри даної процедури це вектора, який копіюємо і у який копіюємо , а також їх розмірність.
3. addVector – сумує два вектори, які приймає у вхідних параметрах.
4. subVector – процедура, що знаходить різницю двох векторів, які приймає у вхідних параметрах.
5. linearSolveFullRelax – головна функція програми знаходження розв’язку СЛАР методом релаксації.

До вихідних даних даної програми належить набір  значень, що є розв’язком СЛАР. Вихідні  дані записуються в файл у кореневій  директорії програми під назвою 'result.txt'.

Текстовий приклад  роботи програми.

Внесемо у файл 'a.txt' наступні дані:

3

10 -2 -2 6

-1 10 -2  7

-1 -1 10 8

Запускаємо  програму (рис. 2.2).

Рисунок 2.2 Демонстрація роботи програми Relax.exe

Перевіряємо файл 'result.txt' і бачимо , що він містить наступні дані (розв’язок СЛАР):

      1.000000

      1.000000

      1.000000

Зауваження. Для успішної компіляції програми в середовищі Turbo Pascal 7.0 необхідно установити режим компіляції  Numeric processing в 8087/80287 (рис. 2.3)

Рисунок 2.3 Режим компіляції Numeric processing 8087/80287.

 

Висновок

 

Можна стверджувати, що майже будь-яке завдання обчислювальної математики зводиться в кінцевому  підсумку до вирішення отриманої  системи лінійних або тензорних  алгебраїчних рівнянь (СЛАР).

Але такі системи рівнянь можуть бути, по-перше, дуже великого розміру, наприклад, NxN = 10000х10000, і навіть більше, по-друге, система рівнянь може виявитися недовизначеною, по-третє, вона може виявитися з лінійно залежними рівняннями, по-четверте, вона може виявитися перевизначеною і несумісною. Крім того, по-п'яте, обчислювальна техніка може мати далеко не рекордну швидкодію і обсяг оперативної пам'яті, і явно кінцеву розрядність двійкового представлення чисел та пов'язані з цим ненульові обчислювальні похибки. Тому ітераційні методи отримали велике застосування у розв’язанні СЛАР. Сучасна обчислювальна техніка дозволяє проводити дослідження стійкості та збіжності ітераційного методу в залежності від параметрів задачі.

Найбільш ефективно  метод релаксацій застосовується при розв’язанні багатьох близьких алгебраїчних систем лінійних рівнянь. На першому етапі проводиться розв’язання однієї з систем з різними значеннями ітераційного параметра w і з аналізу швидкості збіжності ітераційного процесу вибирається оптимальне значення цього параметра. Потім всі інші системи вирішуються з обраним значенням w.

Ще одна перевага ітераційного методу верхніх релаксацій полягає в тому, що при його реалізації на ЕОМ алгоритм обчислень має  простий вигляд і дозволяє використовувати всього один масив для невідомого вектора.

При виконанні  курсової роботи я навчився розв’язувати системи лінійних рівнянь методом релаксації (ослаблення) змінних, і закріпив набуті навички розробкою програми.

 

Список використаних джерел

 

1. Воеводин В.В. «Вычислительные основы линейной алгебры». Москва «Наука», 1977.
2. Фаддеев Д.К., Фаддеева В.Н. «Вычислительные методы линейной алгебры». Москва «Физматгиз», 1963.
3. Самарский А.А., Гулин А.В.» Численные методы». Москва «Наука», 1989.
4. Самарский А.А., Николаев Е.С. «Методы решения сеточных уравнений». Москва «Наука», 1978.
5. Самарский А.А. «Введение в численные методы». Москва «Наука», 1987.
6. Стренг Г. «Линейная алгебра и ее применение». Москва «Мир», 1980.
7. Карманов В.Г. «Математическое программирование». Москва «Наука», 1989.
8. Алексеев Е.Р. «Программирование на С++». Москва «НТ Пресс», 2007.
9. http://www.exponenta.ru/   -  сайт посвящен решению математических задач в прикладных программных пакетах.
10. Бахвалов Н.С., Жидков Н.П., Кобельков Г.М. Численные методы. - М.: Наука, 1987.- 600 с.

 

 

 

 

 

 

 

 

 

Додатки

 

Додаток А

Текст програми

uses crt;

const

   MAX_SIZE=20;                          {максимальний размір матриці}

 

type

   TMatrix= array [1..MAX_SIZE] of array [1..MAX_SIZE] of double;

   TVector= array [1..MAX_SIZE] of double;

 

var

   i, j: integer;

   a: TMatrix;                            {матрица коефіцієнтів}

   f: TVector;                            {вектор правої частини}

   res: TVector;                          {вектор розвязання системи}

   n: integer;                            {размірність системи}

   sigma: double;                         {точність обрахунків}

   max_step: longint;                     {максимальне число кроків}

   data_file: text;                       {файл даних}

   result_file: text;                     {файл результату}

 

{евклідова норма  вектору}

function norma(var r: TVector; n: integer): double;

var

   res: double;

 

begin

   for i:=1 to n do begin

      res:=res+sqr(r[i]);

   end;

 

   norma:=sqrt(res);

end;

 

 

{копіювання  вектора}

procedure copyVector(var v: TVector;

                     var r: TVector;

                     n: integer);

var

   i: integer;

begin

   for i:=1 to n do begin

      r[i]:=v[i];

   end;

end;

 

 

{добуток матриці  на вектор}

procedure multMatrixVector(var m: TMatrix;

                           var v: TVector;

                           var r: TVector;

                           n: integer);

var

   i, j: integer;

   res: double;

 

begin

   for i:=1 to n do begin

      res:=0;

      for j:=1 to n do begin

         res:=res+m[i,j]*v[j];

      end;

      r[i]:=res;

   end;

end;

 

 

{сума векторів}

procedure addVector(var v1: TVector;

                    var v2: TVector;

                    var r: TVector;

                    n: integer);

var

   i: integer;

begin

   for i:=1 to n do begin

      r[i]:=v1[i]+v2[i];

   end;

end;

 

 

{різниця векторів}

procedure subVector(var v1: TVector;

                    var v2: TVector;

                    var r: TVector;

                    n: integer);

var

   i: integer;

begin

   for i:=1 to n do begin

      r[i]:=v1[i]-v2[i];

   end;

end;

 

 

{розв’язання  системи}

procedure linearSolveFullRelax(var a: TMatrix;      {матрица коефіцієнтів}

                               var f: TVector;        {вектор  правої частини}

                               var res: TVector;     {результат}

                               n: integer;                {розмірність системи}

                               sigma: double);       {похибка}

var

   x, r, alpha, t: TVector;

   k, max_step: longint;

 

begin

   k:=0;

   max_step:=10000;

   copyVector(f,x,n);

 

   repeat

      multMatrixVector(a,x,t,n);

      subVector(f,t,r,n);

 

      for i:=1 to n do begin

         alpha[i]:=r[i]/a[i,i];

      end;

 

      addVector(x,alpha,t,n);

      copyVector(t,x,n);

 

      k:=k+1;

   until (k>=max_step) or (norma(r,n)<=sigma);

 

   writeln('Число  кроків: ',k);

 

   if(k=max_step) then begin

      writeln('Задана точність не досягнута');

   end;

 

   copyVector(x,res,n);

end;

 

 

begin

   sigma:=1e-7;

   max_step:=10000;

 

   clrscr;

   writeln('Розвязок  СЛАР');

   writeln('Метод  повної релаксації);

   writeln;

 

   assign(data_file, 'a.txt');

   reset(data_file);

   read(data_file, n);

 

   for i:=1 to n do begin

      for j:=1 to n do begin

         read(data_file, a[i,j]);

      end;

      read(data_file, f[i]);

   end;

 

   close(data_file);

 

   writeln('Розширена  матриця системи:');

   for i:=1 to n do begin

      for j:=1 to n do begin

         write(a[i,j]:14:6);

      end;

      writeln(f[i]:14:6);

   end;

   writeln;

 

   linearSolveFullRelax(a,f,res,n,sigma);

 

   writeln('Результат:');

   for i:=1 to n do begin

      writeln(res[i]:14:6);

   end;

   writeln;

 

   assign(result_file, 'result.txt');

   rewrite(result_file);

 

   for i:=1 to n do begin

      writeln(result_file, res[i]:14:6);

   end;

 

   close(result_file);

 

   readkey;

end.

 

 

 

 

 

 

 

 

 

 

 

Додаток Б

 

 

 

І-31     091851КП	Аркуш