Міністерство освіти України Державний університет "Львівська політехніка" МЕТОДИ РОЗВ'ЯЗУВАННЯ СИСТЕМ ЛІНІЙНИХ ДИФЕРЕНЦІАЛЬНИХ РІВНЯНЬ ІНСТРУКЦІЯ ДО ЛАБОРАТОРНОЇ РОБОТИ N 4 З КУРСУ "ЧИСЕЛЬНІ МЕТОДИ АНАЛІЗУ АВТОМАТИЧНИХ СИСТЕМ" Затверджено на засiданнi кафедри "Автоматика та телемеханiка" Протокол N 12 вiд 15.02.2000p. Львів 2000 Методи розв'язування систем лінійних диференціальних рівнянь: Iнструкцiя до лабораторної роботи N 4 з курсу "Чисельні методи аналізу автоматичних систем управління" для студентiв спецiальностi 6.0914 "Комп'ютеризовані системи, автоматика і управління" / Уклали Л.В.Мороз, А.Е.Лагун - Львiв:Державний університет "Львівська політехніка", 1998.- 11 c. Укладачi: Л.В.Мороз, канд.техн.наук, доцент А.Й.Наконечний, канд.техн.наук, доцент А.Е.Лагун, асистент Вiдповiдальний за випуск І.М.Ковела, канд.техн.наук, доц. Рецензенти: І.М.Ковела, канд.техн.наук, доц. В.В.Самотий, докт.техн.наук, доц. Мета роботи: вивчити основні чисельні методи розв'язування систем лінійних диференціальних рівнянь. ВСТУП Розв'язати систему диференціальних рівнянь першого порядку
(1) чисельним методом означає для заданої послідовності аргументів x
та числа y
, не знаходячи аналітичного виразу функцій-розв'язків y
=F(x), знайти такі значення y
, щоб y
=F(x
), де і=1,...,n. Отже, чисельні методи дозволяють отримати таблицю значень функцій y
=F(x
) для заданої послідовності аргументів (аргументи змінюються з кроком інтегрування h=x
-x
). Знаходження розв'язку системи (1) при початкових умовах y(x
)=y
називається задачею Коші. Ця задача має єдиний розв'язок, який залежить від початкових умов, якщо праві частини (1) неперервні і виконується умова Ліпшица:
(2) де L - константа Ліпшица. Зауважимо, що довільне диференціальне рівняння n-го порядку

(3) можна завжди звести до системи n диференціальних рівнянь першого порядку, шляхом заміни y
(x)=y
(x). Позначимо
(4) Складемо систему диференціальних рівнянь для функцій y
,y
,...,y
. Маємо
(5) Отже, для функцій y
,y
,...,y
одержимо таку систему диференціальних рівнянь:
(6) де y
(0)=y(0), y
(0)=y'(0),..., y
(0)=y
(0). Приклад. Звести диференціальне рівняння другого порядку до системи двох диференціальних рівнянь першого порядку
Введемо заміну y
=y, y
=y'. Тоді отримаємо
1. Методи розв'язування систем диференціальних рівнянь з сталим кроком інтегрування. 1.1. Метод Ейлера. Розглянемо задачу Коші для диференціального рівняння першого порядку:


. (7) Розкладемо розв'язок y(x) в ряд Тейлора в околі точки x
:
(8) Якщо x
змінюється з кроком h, а саме x
=x
+h, припустивши, що сітка x
,x
,...,x
є рівномірною, можна отримати наближене значення y(x) в точці x
: y
=
(9) За цією формулою можна знайти y
, маючи інформацію про попередню точку, але виникають труднощі при знаходженні похідних високих порядків
. Якщо обмежитися i=1, то отримаємо формулу методу Ейлера:
(10) При переході до системи диференціальних рівнянь скалярні величини замінюються на векторні:
(11) Похибка цього методу Ейлера має порядок h
. У модифікованому методі Ейлера використовується така формула:
(12) тобто здійснюємо приріст x з половинним кроком h/2 i, знайшовши нахил інтегральної кривої в цій точці
), за цим нахилом знаходимо приріст функції на цілому кроці h. Для вдосконаленого методу Ейлера знаходять нахил інтегральної кривої в точці x
і в точці x
+h, а потім за середнім нахилом на цілому кроці h визначають приріст функції:
(13) В практичних обчисленнях застосовують також метод Ейлера з ітераціями. В ньому згідно з грубим початковим наближенням
(14) організовують ітераційний процес
(15) де k=1,2,3,... . Ітерації виконують, доки два послідовних наближення не співпадуть із заданою похибкою
. Якщо після кількох ітерацій співпадіння не відбувається, то необхідно зменшити крок. У трьох останніх методах похибка наближено дорівнює h
. Для наближеної оцінки похибки методів Ейлера зручно скористатися так званим "подвійним перерахунком", а саме : знаходять значення
з кроком h, а потім визначають
завдяки подвійному перерахунку з кроком h/2. Різниця між
та
буде похибкою методу. 1.2. Метод Рунге-Кутта. Метод Рунге-Кутта дозволяє будувати схеми різних порядків точності і тому є найбільш вживаним. Формула методу Рунге-Кутта другого порядку точності співпадає з формулою модифікованого методу Ейлера (похибка h
). Алгоритм реалізації методу Рунге-Кутта четвертого порядку точності полягає в циклічному обчисленні y
на кожному кроці за формулами:
, (16) де
Похибка цього методу становить
. Її можна знайти за допомогою "подвійного перерахунку", описаного вище. Існує також метод Рунге-Кутта-Мерсона, похибка інтегрування якого має порядок h
. Він реалізується за наступним алгоритмом:
, (17) де
Цей метод дає наближену оцінку похибки на кожному кроці інтегрування за формулою:
. (18) Метод Рунге-Кутта-Фельберга має також четвертий порядок точності, як і два попередні, але дає більш точну оцінку похибки на кожному кроці інтегрування. Алгоритм реалізації цього методу такий:
(19) де
Похибка
(20) 2. Методи розв'язування систем диференціальних рівнянь зі змінним кроком інтегрування. При розв'язуванні систем диференціальних рівнянь із змінним кроком інтегрування формули методів пункту 1 залишаються тими ж самими, лише необхідно ввести автоматичну зміну кроку інтегрування в залежності від заданої похибки. Для методів Ейлера, а також методу Рунге-Кутта 4-го порядку це здійснюється наступним чином. Знайшовши похибку відповідного методу на кроці
, порівнюють її із заданою
. Якщо задана похибка є меншою, то крок зменшують у два рази і проводять обчислення з цим кроком. Потім знову визначають похибку і перевіряють умову
>
. Крок ділять до того часу, доки ця умова не виконається, після чого наступне обчислення проводиться з останнім кроком. Якщо знайдена перший раз похибка методу є меншою від заданої, то одержане значення виводять на друк, а крок подвоюють. Так роблять до того часу, доки не виконається умова
<
; тоді крок знову починають ділити. Для методів Рунге-Кутта-Мерсона і Рунге-Кутта-Фельберга оцінку похибки методів проводимо за залишковим членом R
. Якщо
то крок збільшують у 2 рази, якщо
, то крок залишають тим самим, що і при попередньому обчисленні, інакше при
крок зменшують в 2 рази. Отже, при змінному кроці інтегрування існує можливість розв'язування системи диференціальних рівнянь із заданою похибкою. ЗАВДАННЯ. Скласти програму для розв'язування диференціального рівняння другого порядку відповідним чисельним методом, попередньо звівши його до системи диференціальних рівнянь першого порядку, якщо задані початкові умови для функції y(0) та її похідної y'(0), проміжок інтегрування [a,b], на якому міняється аргумент x з кроком h (у випадку автоматичного вибору кроку інтегрування заданий крок h є початковим, а похибка -
). Варіанти завдань і метод розв'язку беруться за вказівкою викладача. У випадку сталого кроку інтегрування необхідно знайти похибку методу на кожному кроці інтегрування. Таблиця 1 №
Диференціальне рівняння Початкові умови Проміжок інтегрування Крок інтегрування Похибка  1 
y(0)=1 y’(0)=0 [0;0.5] 0.1 0.001  2 
y(0)=1 y’(0)=1 [0;0.5] 0.05 0.001  3 
y(0)=1 y’(0)=0 [0;0.5] 0.05 0.001  4 
y(0)=0 y’(0)=-4 [0;0.2] 0.02 0.001  5 
y(0)=1 y’(0)=-1 [0;1] 0.1 0.001  6 
y(0)=1.8 y’(0)=-0.5 [0;2] 0.2 0.001  7 
y(0)=2.2 y’(0)=0.8 [0;0.2] 0.02 0.001  8 
y(0)=0.8 y’(0)=2 [0;1] 0.1 0.001  9 
y(0)=1.43 y’(0)=-0.3 [0;1] 0.1 0.001  10 
y(0)=1 y’(0)=2 [0;1] 0.1 0.001  11 
y(0)=-0.9 y’(0)=3.2 [0;1] 0.1 0.001  12 
y(0)=1/9 y’(0)=-5/9 [0;1] 0.1 0.001  13 
y(0)=5/6 y’(0)=2/3 [1;2] 0.1 0.001  14 
y(0)=1.48 y’(0)=3.6 [0;0.5] 0.05 0.001  15 
y(0)=0 y’(0)=0 [0;0.2] 0.02 0.001  16 
y(0)=5.1 y’(0)=4.2 [0;1] 0.1 0.001   Продовження таблиці 1 №
Диференціальне рівняння Початкові умови Проміжок інтегрування Крок інтегрування Похибка  17 
y(0)=2.5 y’(0)=1.5 [0;1] 0.1 0.001  18 
y(1)=4 y’(1)=4 [1;1.5] 0.05 0.001  19 
y(1)=5 y’(1)=-1 [1;1.5] 0.05 0.001  20 
y(0)=1 y’(0)=2 [0;0.5] 0.05 0.001  21 
y(0)=2.6 y’(0)=3.2 [0;1] 0.1 0.001  22 
y(0)=2 y’(0)=3.5 [1;2] 0.1 0.001  23 
y(0)=1.38 y’(0)=-0.2 [1;2] 0.1 0.001  24 
y(0)=1.43 y’(0)=2.3 [1;2] 0.1 0.001  25 
y(0)=11/9 y’(0)=-9/8 [0;1] 0.1 0.001   ЛІТЕРАТУРА 1. Краскевич В.Е., Зеленский К.Х., Гречко В.И. Численные методы в инженерных исследованиях.- К.: Вища шк. Головное изд-во, 1986.- 263 c. 2. Прикладные программы и лабораторный практикум для персонального компьютера: Учебное пособие / Гаврилюк М.А., Галамай Т.Г., Коссак О.М.и др. - Киев: УМК ВО, 1988. - 205 c. 3. Дьяконов В.П. Справочник по алгоритмам и программам на языке бейсик для персональных ЭВМ: Справочник. - М.: Наука. Гл. ред. физ.-мат. лит., 1987. - 240 с. Лященко М. Я., Головань М. С. Чисельні методи. - К.: Либідь, 1996. 5. Микрокалькуляторы в курсе высшей математики: Практикум / Барановская Г.Г., Любченко И.Н.- К. : Вища шк. Головное изд-во, 1987 - 288 с. Алгоритмы и программы на Бейсике: Учеб. пособие для студентов пед. ин-тов по физ.-мат. спец. / Гринчишин Я.Т., Ефимов В.И., Ломакович А.Н. - М.: Просвещение, 1988. - 160 с. Самохин А. Б., Самохина А. С. Численные методы и программирование на Фортране для персонального компьютера. - М.: Радио и связь, 1996. Поджаренко В. О., Сокол В. М., Селезньова Р. В. Математичні обчислення в Maple for Windows. - Вінниця: Континент - ПРИМ, 1997, - 106 с.