Чисельне вирішення диференціальних рівнянь
Багато завдань науки та техніки зводяться до вирішення звичайних диференціальних рівнянь (ОДП). ОДУ називаються такі рівняння, які містять одну або кілька похідних від функції, що шукається. У загальному вигляді ОДУ можна записати так:
Де x - незалежна змінна, - i-а похідна від шуканої функції. n – порядок рівняння. Загальне рішення ОДУ n–го порядку містить n довільних постійних, тобто. загальне рішення має вигляд.
Для виділення єдиного рішення необхідно встановити n додаткових умов. Залежно від способу завдання додаткових умов існують два різні типи завдань: завдання Коші та крайове завдання. Якщо додаткові умови задаються в одній точці, така задача називається завданням Коші. Додаткові умови завдання Коші називаються початковими умовами. Якщо ж додаткові умови задаються більш ніж у одній точці, тобто. при різних значеннях незалежної змінної, то таке завдання називається крайовим. Самі додаткові умови називаються крайовими чи граничними.
Зрозуміло, що з n=1 можна говорити лише завдання Коші.
Приклади постановки завдання Коші:
Приклади крайових завдань:
Вирішити такі завдання аналітично вдається лише деяких спеціальних типів рівнянь.
Численні методи вирішення задачі Коші для ОДУ першого порядку
Постановка задачі. Знайти рішення ОДУ першого порядку
На відрізку за умови
При знаходженні наближеного рішення вважатимемо, що обчислення проводяться з розрахунковим кроком , розрахунковими вузлами є точки проміжку [ x 0 , x n ].
Метою є побудова таблиці
|
x i |
x n |
|||
|
y i |
y n |
тобто. шукаються наближені значення y у вузлах сітки.
Інтегруючи рівняння на відрізку, отримаємо
Цілком природним (але не єдиним) шляхом отримання чисельного рішення є заміна в ньому інтеграла якоюсь квадратурною формулою чисельного інтегрування. Якщо скористатися найпростішою формулою лівих прямокутників першого порядку
,
то отримаємо явну формулу Ейлера:
Порядок розрахунків:
Знаючи, знаходимо, потім т.д.
Геометрична інтерпретація методу Ейлера:
Користуючись тим, що у точці x 0 відоме рішення y(x 0)= y 0 і значення його похідної можна записати рівняння дотичної до графіка шуканої функції в точці :. При досить малому кроці hордината цієї дотичної, отримана підстановкою в праву частину значення, повинна мало відрізнятися від ординати y(x 1) рішення y(x) Завдання Коші. Отже, точка перетину дотичної з прямою x = x 1 може бути приблизно прийнята за нову початкову точку. Через цю точку знову проведемо пряму, яка приблизно відображає поведінку дотичної в точці. Підставляючи сюди (тобто перетин з прямою x = x 2), отримаємо наближене значення y(x) у точці x 2: і т.д. У результаті для i-ї точки отримаємо формулу Ейлера.
Явний метод Ейлер має перший порядок точності або апроксимації.
Якщо використовувати формулу правих прямокутників: 
, то прийдемо до методу
Цей метод називають неявним методом Ейлераоскільки для обчислення невідомого значення за відомим значенням потрібно вирішувати рівняння, у загальному випадку нелінійне.
Неявний метод Ейлер має перший порядок точності або апроксимації.
У цьому методі обчислення складається з двох етапів:
Ця схема називається ще методом предиктор – коректор (передбачуване – виправляюче). У першому етапі наближене значення передбачається з невисокою точністю (h), але в другому етапі це передбачення виправляється, отже результуюче значення має другий порядок точності.
Методи Рунге - Кутта:ідея побудови явних методів Рунге-Кутти p-го порядку полягає в отриманні наближень до значень y(x i+1) за формулою виду
…………………………………………….
Тут a n , b nj , p n, - Деякі фіксовані числа (параметри).
При побудові методів Рунге-Кутти параметри функції ( a n , b nj , p n) підбирають таким чином, щоб отримати потрібний порядок апроксимації.
Схема Рунге – Кутта четвертого порядку точності:
приклад. Вирішити завдання Коші:
Розглянути три методи: явний метод Ейлера, модифікований метод Ейлера, метод Рунґе - Кутта.
Точне рішення:
Розрахункові формули за явним методом Ейлера для цього прикладу:
Розрахункові формули модифікованого методу Ейлера:
Розрахункові формули методу Рунге - Кутта:
y1 – метод Ейлера, y2 – модифікований метод Ейлера, y3 – метод Рунґе Кутта.
Видно, що найточнішим є метод Рунґе – Кутта.
Численні методи вирішення систем ОДУ першого порядку
Розглянуті методи можуть бути використані для вирішення систем диференціальних рівнянь першого порядку.
Покажемо це для випадку системи двох рівнянь першого порядку:
Явний метод Ейлера:
Модифікований метод Ейлера:
Схема Рунге - Кутта четвертого порядку точності:
До вирішення систем рівнянь ОДУ зводяться завдання Коші для рівнянь вищих порядків. Наприклад, розглянемо завдання Коші для рівняння другого порядку
Введемо другу невідому функцію. Тоді завдання Коші замінюється наступним:
Тобто. у термінах попереднього завдання: .
приклад. Знайти рішення задачі Коші:
На відрізку.
Точне рішення:
Дійсно:
Розв'яжемо завдання явним методом Ейлера, модифікованим методом Ейлера і Рунге - Кутта з кроком h = 0.2.
Введемо функцію.
Тоді отримаємо наступне завдання Коші для системи двох ОДУ першого порядку:
Явний метод Ейлера:
Модифікований метод Ейлера:
Метод Рунге - Кутта:
Схема Ейлера:
Модифікований метод Ейлера:
Схема Рунге - Кутта:
Max(y-y теор) = 4 * 10 -5
Метод кінцевих різниць вирішення крайових завдань для ОДУ
Постановка задачі: знайти рішення лінійного диференціального рівняння
задовольняє крайових умов:. (2)
Теорема.Нехай. Тоді існує єдине рішення поставленого завдання.
До цієї задачі зводиться, наприклад, завдання про визначення прогинів балки, яка на кінцях спирається шарнірно.
Основні етапи методу кінцевих різниць:
1) область безперервної зміни аргументу () замінюється дискретним безліччю точок, званих вузлами: .
2) Шукана функція безперервного аргументу x, приблизно замінюється функцією дискретного аргументу на заданій сітці, тобто. . Функція називається сітковою.
3) Вихідне диференціальне рівняння замінюється різницевим рівнянням щодо сіткової функції. Така заміна називається різницевою апроксимацією.
Таким чином, рішення диференціального рівняння зводиться до відшукання значень сіткової функції у вузлах сітки, які перебувають з розв'язання рівнянь алгебри.
Апроксимація похідних.
Для апроксимації (заміни) першої похідної можна скористатися формулами:
- права похідна різницева,
- ліва різницева похідна,
Центральна похідна різницева.
тобто, можливо безліч способів апроксимації похідної.
Всі ці визначення випливають із поняття похідної як межі: 
.
Спираючись на різницеву апроксимацію першої похідної можна побудувати різницеву апроксимацію другої похідної:
Аналогічно можна отримати апроксимацію похідних вищого порядку.
Визначення.Похибкою апроксимації n-ої похідної називається різниця: .
Для визначення порядку апроксимації використовується розкладання до ряду Тейлора.
Розглянемо праву різницеву апроксимацію першої похідної:
Тобто. права різницева похідна має перший по hпорядок апроксимації.
Аналогічно і для лівої похідної різницевої.
Центральна різницева похідна має другий порядок апроксимації.
Апроксимація другої похідної за формулою (3) також має другий порядок апроксимації.
Для того щоб апроксимувати диференціальне рівняння необхідно в ньому замінити всі похідні апроксимаціями. Розглянемо завдання (1), (2) і замінимо на (1) похідні:
В результаті отримаємо:
(4)
Порядок апроксимації вихідного завдання дорівнює 2, т.к. друга та перша похідні замінені з порядком 2, а решта – точно.
Отже, замість диференціальних рівнянь (1) (2) отримана система лінійних рівнянь для визначення у вузлах сітки.
Схему можна представити у вигляді:
тобто, отримали систему лінійних рівнянь із матрицею:
Ця матриця є трехдиагональной, тобто. всі елементи, які розташовані не на головній діагоналі та двох прилеглих до неї діагоналях дорівнюють нулю.
Вирішуючи отриману систему рівнянь, ми отримаємо вирішення вихідного завдання.
Диференціальні рівняння - це рівняння, які невідома функція входить під знаком похідної. Основне завдання теорії диференціальних рівнянь - вивчення функцій, що є рішеннями таких рівнянь.
Диференціальні рівняння можна розділити на звичайні диференціальні рівняння, у яких невідомі функції є функціями однієї змінної, і диференційні рівняння у приватних похідних, у яких невідомі функції є функціями двох і більшої кількості змінних.
Теорія диференціальних рівнянь у похідних більш складна і розглядається в більш повних або спеціальних курсах математики.
Вивчення диференціальних рівнянь почнемо з найпростішого рівняння - рівняння першого порядку.
Рівняння виду
F(x,y,y") = 0,(1)
де х - незалежна змінна; у - потрібна функція; у - її похідна, називається диференціальним рівнянням першого порядку.
Якщо рівняння (1) можна вирішити щодо у", то воно набуває вигляду
і називається рівнянням першого порядку, дозволеним щодо похідної.
У деяких випадках рівняння (2) зручно записати у вигляді f (х, у) dх - dy = 0, що є окремим випадком більш загального рівняння
P(x,y)dx+Q(x,y)dy=O,(3)
де Р(х,у) та Q(х,у) - відомі функції. Рівняння в симетричній формі (3) зручне тим, що змінні х і в ньому рівноправні, тобто кожну з них можна розглядати як іншу функцію.
Дамо два основні визначення загального та приватного рішення рівняння.
Загальним рішенням рівняння (2) в деякій області G площині Оху називається функція у=ц(х,С), яка залежить від х і довільної постійної З, якщо вона є рішенням рівняння (2) при будь-якому значенні постійної З, і якщо за будь-яких початкових умовах y x=x0 =y 0 таких, що (x 0 ;y 0)=G, існує єдине значення постійної З = З 0 таке, що функція у = ц (х, З 0) задовольняє цим початковим умовам у = ц (х 0, С).
Приватним рішенням рівняння (2) в області G називається функція у = ц (х, З 0), яка виходить із загального рішення у = ц (х, С) при певному значенні постійної З = З 0 .
Геометрично загальне рішення у=ц(х,С) є сімейством інтегральних кривих на площині Оху, що залежить від однієї довільної постійної З, а приватне рішення у=ц(х,С 0) - одну інтегральну криву цього сімейства, що проходить через задану точку (х 0; у 0).
Наближене розв'язання диференціальних рівнянь першого порядку методом Ейлера. Суть цього методу полягає в тому, що шукана інтегральна крива, що є графіком приватного рішення, приблизно замінюється ламаною. Нехай дані диференціальне рівняння
та початкові умови y | x = x0 = y 0 .
Знайдемо приблизно рішення рівняння на відрізку [х 0, b], що задовольняє заданим початковим умовам.
Розіб'ємо відрізок [х 0, b] точками х 0<х 1 ,<х 2 <...<х n =b на n равных частей. Пусть х 1 --х 0 =х 2 -- x 1 = ... =x n -- x n-1 = ?x. Обозначим через y i приближенные значения искомого решения в точках х i (i=1, 2, ..., n). Проведем через точки разбиения х i - прямые, параллельные оси Оу, и последовательно проделаем следующие однотипные операции.
Підставимо значення х 0 і у 0 праву частину рівняння y"=f(x,y) і обчислимо кутовий коефіцієнт y"=f(x 0 ,y 0) дотичної до інтегральної кривої в точці (х 0 ;у 0). Для знаходження наближеного значення у 1 шуканого рішення замінюємо на відрізку [х 0 x 1] інтегральну криву відрізком її дотичної в точці (х 0 у 0). При цьому отримуємо
y 1 - y 0 = f (x 0; y 0) (x 1 - x 0),
звідки, оскільки х 0 , х 1 , у 0 відомі, знаходимо
y1 = y0 + f (x0; y0) (x1 - x0).
Підставляючи значення х 1 і y 1 в праву частину рівняння y"=f(x,y), обчислюємо кутовий коефіцієнт y"=f(x 1 ,y 1) дотичної до інтегральної кривої в точці (х 1 ;y 1). Далі, замінюючи на відрізку інтегральну криву відрізком дотичної, знаходимо наближене значення рішення у 2 у точці х 2:
y 2 = y 1 +f(x 1; y 1)(x 2 - x 1)
У цьому рівні відомими є х 1 , у 1 , х 2 , а й у 2 виражається них.
Аналогічно знаходимо
y 3 = y 2 +f(x 2 ;y 2) ?x, …, y n = y n-1 +f(x n-1 ;y n-1) ?x
Таким чином, приблизно побудована шукана інтегральна крива у вигляді ламаної та отримані наближені значення y i шуканого рішення в точках х i . При цьому значення у i обчислюються за формулою
y i = y i-1 + f (x i-1; y i-1) x (i = 1,2, …, n).
Формула і є основною розрахунковою формулою методу Ейлера. Її точність тим вища, що менше різниця?x.
Метод Ейлера ставитися до чисельних методів, що дає рішення у вигляді таблиці наближених значень шуканої функції у(х). Він є порівняно грубим і застосовується переважно для орієнтовних розрахунків. Однак ідеї, покладені в основу методу Ейлера, є вихідними для інших методів.
Ступінь точності методу Ейлера, власне кажучи, невелика. Існують набагато точніші методи наближеного розв'язання диференціальних рівнянь.
Для вирішення диференціальних рівнянь необхідно знати значення залежної змінної та її похідних при деяких значеннях незалежної змінної. Якщо додаткові умови задаються за одного значення невідомої, тобто. незалежною змінною., то таке завдання називається завданням Коші. Якщо початкові умови задаються при двох або більше значеннях незалежної змінної, завдання називається крайової. При вирішенні диференціальних рівнянь різних видів, функція, значення якої потрібно визначити, обчислюється у вигляді таблиці.
Завдання Коші - однокрокові: методи Ейлера, методи Рунге-Кутта; - Багатокрокові: метод Майна, Метод Адамса. Кроєве завдання – метод зведення кроєвої задачі до задачі Коші; -Метод кінцевих різниць.
При розв'язанні задачі Коші мають бути задані дифр. ур. порядку n чи система дифр. ур. першого порядку з n рівнянь та n додаткових умов її вирішення. Додаткові умови повинні бути задані при тому самому значенні незалежної змінної. При вирішенні кроєвої задачі мають бути задані ур. n-ого порядку чи система з n рівнянь і n додаткових умов за двох чи більше значення незалежної змінної. При розв'язанні задачі Коші функція визначається дискретно у вигляді таблиці з деяким заданим кроком . При визначенні кожного чергового значення можна використовувати інформацію про одну попередню точку. У цьому випадку методи називають однокроковим, або можна використовувати інформацію про кілька попередніх точок - багатокрокові методи.
Задано: g(x,y)y+h(x,y)=0, y=-h(x,y)/g(x,y)= f(x,y), x 0 , y( x 0) = y 0. Відомо: f(x,y), x 0 , y 0 . Визначити дискретне рішення: x i, y i, i = 0,1, ..., n. Метод Ейлера заснований на розкладанні функції до ряду Тейлора околиці точки x 0 . Околиця описується кроком h. y(x 0 +h)y(x 0)+hy(x 0)+…+ (1). У методі Ейлера враховуються лише два складові ряду Тейлора. Введемо позначення. Формула Ейлера набуде вигляду: y i + 1 = y i + i y , y y = hy (x i) = hf (x i, y i), y i + 1 = y i + hf (x i, y i) (2), i = 0,1,2…, x i+1 =x i +h
Формула (2) є формулою найпростішого методу Ейлера.
Для отримання чисельного рішення використовується ф-ла дотичної, що проходить через урав. дотичної: y=y(x 0)+y(x 0)(x-x 0), x=x 1 ,
y 1 =y(x 0)+f(x 0 ,y 0) (x-x 0), т.к.
x-x 0 = h, то y 1 = y 0 + hf (x 0, y 0), f (x 0, y 0) = tg £.
Задано: y=f(x,y), y(x0)=y0. Відомо: f(x,y), x 0 , y 0 . Визначити: залежність y від x у вигляді табличної дискретної функції: x i, y i, i = 0,1, ..., n.
1) обчислимо тангенс кута нахилу у початковій точці
tg £=y(x n ,y n)=f(x n ,y n)
2) Обчислимо значення y n+1 на
наприкінці кроку за формулою Ейлера
y n+1 =y n +f(x n ,y n) 3) Обчислимо тангенс кута нахилу
дотичної в n+1 точці: tg £=y(x n+1 , y n+1)=f(x n+1 , y n+1) 4) Обчислимо середнє арифметичне кутів
нахилу: tg £=½. 5) Використовуючи тангенс кута нахилу, перерахуємо значення функції в n+1 точці: y n+1 =y n +htg £= y n +½h=y n +½h – формула модифікованого методу Ейлера. Можна показати, що отримана ф-ла відповідає розкладу ф-ії до ряду Тейлора, включаючи доданки (до h 2). Модифікований метод Ейлнра на відміну простого є методом другого порядку точності, т.к. похибка пропорційна h2.
Кафедра фізхімії ПФУ (РГУ)
ЧИСЛІВІ МЕТОДИ І ПРОГРАМУВАННЯ
Матеріали до лекційного курсу
Лектор – ст. викл. Щербаков І.М.
Постановка задачі
При вирішенні наукових та інженерно-технічних завдань часто буває необхідно математично описати будь-яку динамічну систему. Найкраще це робити у вигляді диференціальних рівнянь ( ДК) або системи диференціальних рівнянь. Найбільш часто вони таке завдання виникає при вирішенні проблем, пов'язаних з моделюванням кінетики хімічних реакцій та різних явищ перенесення (тепла, маси, імпульсу) – теплообміну, перемішування, сушіння, адсорбції, при описі руху макро- та мікрочастинок.
Звичайним диференціальним рівнянням(ОДУ) n-го порядку називається наступне рівняння, яке містить одну або кілька похідних від функції y(x):
Тут y(n)позначає похідну порядку n деякої функції y(x), x це незалежна змінна.
У ряді випадків диференціальне рівняння можна перетворити на вид, у якому старша похідна виражена у явному вигляді. Така форма запису називається рівнянням, дозволеним щодо старшої похідної(при цьому у правій частині рівняння старша похідна відсутня):
Саме така форма запису прийнята як стандартноюпід час розгляду чисельних методів рішення ОДУ.
Лінійним диференціальним рівняннямназивається рівняння, лінійне щодо функції y(x) та всіх її похідних.
Наприклад, нижче наведені лінійні ОДУ першого та другого порядків
Рішенням звичайного диференціального рівнянняназивається така функція y(x), яка за будь-яких х задовольняє цьому рівнянню в певному кінцевому або нескінченному інтервалі. Процес розв'язання диференціального рівняння називають інтегруванням диференціального рівняння.
Загальне рішення ОДУ n-го порядку містить n довільних констант C1, C2, …, Cn
Це очевидно випливає з того, що невизначений інтеграл дорівнює першорядному підінтегральному виразу плюс константа інтегрування.
Оскільки рішення ДУ n -го порядку необхідно провести n інтегрувань, то загальному рішенні з'являється n констант інтегрування.
Приватне рішенняОДУ виходить із загального, якщо константам інтегрування надати деякі значення, визначивши деякі додаткові умови, кількість яких дозволяє обчислити всі невизначені константи інтегрування.
Точне (аналітичне) рішення (загальне чи приватне) диференціального рівняння передбачає отримання шуканого рішення (функції y(x)) як висловлювання від елементарних функцій. Це можливо далеко не завжди навіть для рівнянь першого порядку.
Чисельне рішення ДУ (приватне) полягає у обчисленні функції y(x) та її похідних у деяких заданих точках , що лежать на певному відрізку. Тобто, фактично, рішення ДУ n-го порядку виду виходить у вигляді наступної таблиці чисел (стовпець значень старшої похідної обчислюється підстановкою значень рівняння):
Наприклад, для диференціального рівняння першого порядку таблиця рішення буде два стовпці – x і y .
Безліч значень абсцис у яких визначається значення функції, називають сіткою, де визначена функція y(x) . Самі координати при цьому називають вузлами сітки. Найчастіше, для зручності, використовуються рівномірні сітки, у яких різниця між сусідніми вузлами постійна і називається кроком сіткиабо кроком інтегруваннядиференціального рівняння
Або , i= 1, …, N
Для визначення приватного рішеннянеобхідно встановити додаткові умови, які дозволять обчислити константи інтегрування. Причому таких умов має бути рівно n. Для рівнянь першого порядку – одне, другого – 2 тощо. Залежно від способу їхнього завдання при вирішенні диференціальних рівнянь існують три типи завдань:
· Завдання Коші (початкове завдання): Необхідно знайти таке приватне рішеннядиференціального рівняння, яке задовольняє певним початковими умовами, заданими в одній точці:
тобто, задано певне значення незалежної змінної (х 0) і значення функції і всіх її похідних аж до порядку (n-1) у цій точці. Ця точка (х 0) називається початковий. Наприклад, якщо вирішується ДК 1-го порядку, то початкові умови виражаються у вигляді кількох чисел (x 0 , y 0)
Такі завдання зустрічається при вирішенні ОДУ, які описують, наприклад, кінетику хімічних реакцій В цьому випадку відомі концентрації речовин у початковий момент часу ( t = 0) і необхідно знайти концентрації речовин через деякий проміжок часу ( t). Як приклад можна також навести завдання про теплоперенесення чи масопереносі (дифузії), рівняння руху матеріальної точки під впливом сил тощо.
· Крайове завдання . У цьому випадку відомі значення функції та (або) її похідних у більш ніж одній точці, наприклад, початковий і кінцевий момент часу, і необхідно знайти приватне рішення диференціального рівняння між цими точками. Самі додаткові умови у цьому випадку називаються крайовими (граничними) умовами. Природно, що крайове завдання може вирішуватися для ОДУ не нижче 2-го порядку. Нижче наведено приклад ОДУ другого порядку з граничними умовами (задані значення функції у двох різних точках):
· Завдання Штурма-Ліувіля (завдання на власні значення). Завдання цього схожі на крайову завдання. При їх вирішенні необхідно знайти, за яких значень будь-якого параметра рішення ДКзадовольняє крайовим умовам (власні значення) та функції, які є рішенням дистанційного керування при кожному значенні параметра (власні функції). Наприклад, багато завдань квантової механіки є завданнями на власні значення.
Чисельні методи вирішення задачі Коші ОДУ першого порядку
Розглянемо деякі чисельні методи розв'язання завдання Коші(Початкового завдання) звичайних диференціальних рівнянь першого порядку. Запишемо дане рівняння у загальному вигляді, дозволеному щодо похідної (права частина рівняння не залежить від першої похідної):
(6.2)
Необхідно знайти значення функції y у заданих точках сітки, якщо відомі початкові значення, де є значення функції y(x) у початковій точці x0.
Перетворимо рівняння множенням на d x
І проінтегруємо ліву та праву частини між i-им та i+ 1-им вузлами сітки.
(6.3)
Ми отримали вираз для побудови рішення в i+1 вузлі інтегрування через значення x і y в i-му вузлі сітки. Складність, однак, полягає в тому, що інтеграл у правій частині є інтеграл від неявно заданої функції, знаходження якого в аналітичному вигляді у загальному випадку неможливе. Численні методи рішення ОДУ різним способом апроксимують (наближають) значення цього інтеграла для побудови формул чисельного інтегрування ОДУ.
З безлічі розроблених на вирішення ОДУ першого порядку методів розглянемо методи , і . Вони досить прості і дають початкове уявлення про підходи до вирішення цієї задачі в рамках чисельного рішення.
Метод Ейлера
Історично першим і найпростішим способом чисельного розв'язання задачі Коші для ОДУ першого ладу є метод Ейлера. У його основі лежить апроксимація похідної ставленням кінцевих прирощень залежною ( y) та незалежної ( x) змінних між вузлами рівномірної сітки:
де y i+1 це потрібне значення функції в точці x i+1 .
Якщо тепер перетворити це рівняння і врахувати рівномірність сітки інтегрування, то вийде ітераційна формула, за якою можна обчислити y i+1якщо відомо y i в точці х i :
Порівнюючи формулу Ейлера із загальним виразом, отриманим раніше, видно, що для наближеного обчислення інтеграла в методі Ейлера використовується найпростіша формула інтегрування - формула прямокутників по лівому краю відрізка.
Графічна інтерпретація методу Ейлера також не становить труднощів (див. малюнок нижче). Дійсно, виходячи з виду розв'язуваного рівняння () слід, що значення є значення похідної функції y(x) у точці x=x i - , і, таким чином, дорівнює тангенсу кута нахилу дотичної, проведеної до графіка функції y(x) у точці x = x i.
З прямокутного трикутника на малюнку можна знайти
звідки і виходить формула Ейлера. Таким чином, суть методу Ейлера полягає в заміні функції y(x) на відрізку інтегрування прямою лінією, що стосується графіка в точці x=x i . Якщо потрібна функція сильно відрізняється від лінійної на відрізку інтегрування, то похибка обчислення буде значною. Помилка методу Ейлера прямо пропорційна кроку інтегрування:
Помилка~ h
Процес обчислень будується в такий спосіб. За заданих початкових умов x 0і y 0можна обчислити
Таким чином, будується таблиця значень функції y(x) з певним кроком ( h) за xна відрізку. Помилка у визначенні значення y(x i)при цьому буде тим менше, чим менше вибрано довжину кроку h(що визначається точністю формули інтегрування).
При великих методах Ейлера дуже неточний. Він дає дедалі точніше наближення зі зменшенням кроку інтегрування. Якщо відрізок занадто великий, то кожна ділянка розбивається на N відрізків інтегрування і до кожного з них застосовується формула Ейлера з кроком, тобто крок інтегрування h береться менше за крок сітки, на якій визначається рішення.
Приклад:
Використовуючи метод Ейлера, побудувати наближене рішення для наступного завдання Коші:
На сітці з кроком 0,1 в інтервалі (6.5)
Рішення:
Дане рівняння вже записано у стандартному вигляді, дозволеному щодо похідної шуканої функції.
Тому для вирішуваного рівняння маємо
Приймемо крок інтегрування рівним кроці сітки h = 0,1. При цьому для кожного вузла сітки буде обчислено лише одне значення (N=1). Для перших чотирьох вузлів сітки обчислення будуть наступними:
Повні результати (з точністю до п'ятого знака після коми) наведені в третій колонці - h = 0,1 (N = 1). У другій колонці таблиці для порівняння наведено значення, обчислені за аналітичним рішенням даного рівняння 
.
У другій частині таблиці наведено відносну похибку отриманих рішень. Видно, що за h =0,1 похибка дуже велика, досягаючи 100% для першого вузла x =0,1.
Таблиця 1 Рішення рівняння методом Ейлера (для колонок вказано крок інтегрування та кількість відрізків інтегрування N між вузлами сітки)
| x | Точне Рішення | 0,1 | 0,05 | 0,025 | 0,00625 | 0,0015625 | 0,0007813 | 0,0001953 |
|---|---|---|---|---|---|---|---|---|
| 1 | 2 | 4 | 16 | 64 | 128 | 512 | ||
| 0 | 0,000000 | 0,000000 | 0,000000 | 0,000000 | 0,000000 | 0,000000 | 0,000000 | 0,000000 |
| 0,1 | 0,004837 | 0,000000 | 0,002500 | 0,003688 | 0,004554 | 0,004767 | 0,004802 | 0,004829 |
| 0,2 | 0,018731 | 0,010000 | 0,014506 | 0,016652 | 0,018217 | 0,018603 | 0,018667 | 0,018715 |
| 0,3 | 0,040818 | 0,029000 | 0,035092 | 0,037998 | 0,040121 | 0,040644 | 0,040731 | 0,040797 |
| 0,4 | 0,070320 | 0,056100 | 0,063420 | 0,066920 | 0,069479 | 0,070110 | 0,070215 | 0,070294 |
| 0,5 | 0,106531 | 0,090490 | 0,098737 | 0,102688 | 0,105580 | 0,106294 | 0,106412 | 0,106501 |
| 0,6 | 0,148812 | 0,131441 | 0,140360 | 0,144642 | 0,147779 | 0,148554 | 0,148683 | 0,148779 |
| 0,7 | 0,196585 | 0,178297 | 0,187675 | 0,192186 | 0,195496 | 0,196314 | 0,196449 | 0,196551 |
| 0,8 | 0,249329 | 0,230467 | 0,240127 | 0,244783 | 0,248202 | 0,249048 | 0,249188 | 0,249294 |
| 0,9 | 0,306570 | 0,287420 | 0,297214 | 0,301945 | 0,305423 | 0,306284 | 0,306427 | 0,306534 |
| 1 | 0,367879 | 0,348678 | 0,358486 | 0,363232 | 0,366727 | 0,367592 | 0,367736 | 0,367844 |
Відносні похибки обчислених значень функції при різних h |
||||||||
| x | h | 0,1 | 0,05 | 0,025 | 0,00625 | 0,0015625 | 0,0007813 | 0,0001953 |
| N | 1 | 2 | 4 | 16 | 64 | 128 | 512 | |
| 0,1 | 100,00% | 48,32% | 23,76% | 5,87% | 1,46% | 0,73% | 0,18% | |
| 0,2 | 46,61% | 22,55% | 11,10% | 2,74% | 0,68% | 0,34% | 0,09% | |
| 0,3 | 28,95% | 14,03% | 6,91% | 1,71% | 0,43% | 0,21% | 0,05% | |
| 0,4 | 20,22% | 9,81% | 4,83% | 1,20% | 0,30% | 0,15% | 0,04% | |
| 0,5 | 15,06% | 7,32% | 3,61% | 0,89% | 0,22% | 0,11% | 0,03% | |
| 0,6 | 11,67% | 5,68% | 2,80% | 0,69% | 0,17% | 0,09% | 0,02% | |
| 0,7 | 9,30% | 4,53% | 2,24% | 0,55% | 0,14% | 0,07% | 0,02% | |
| 0,8 | 7,57% | 3,69% | 1,82% | 0,45% | 0,11% | 0,06% | 0,01% | |
| 0,9 | 6,25% | 3,05% | 1,51% | 0,37% | 0,09% | 0,05% | 0,01% | |
| 1 | 5,22% | 2,55% | 1,26% | 0,31% | 0,08% | 0,04% | 0,01% | |
Зменшити крок інтегрування вдвічі, h = 0.05, у цьому випадку для кожного вузла сітки обчислення буде проводитися за два кроки (N = 2). Так, для першого вузла x = 0,1 отримаємо:
(6.6)
Ця формула виявляється неявною щодо y i+1 (це значення є і в лівій і в правій частині виразу), тобто є рівнянням щодо y i+1 , вирішувати яке можна, наприклад, чисельно, застосовуючи будь-який ітераційний метод (у такому його можна розглядати як ітераційну формулу методу простої ітерації). Однак, можна вчинити інакше і приблизнообчислити значення функції у вузлі i+1за допомогою звичайної формули:
,
яке потім використовувати при обчисленні (6.6).
Таким чином виходить метод Гюначи метод Ейлера з перерахунком. Для кожного вузла інтегрування проводиться наступний ланцюжок обчислень
(6.7)
Завдяки точнішій формулі інтегрування, похибка методу Ґюна пропорційна вже квадрату кроку інтегрування.
Помилка~ h 2
Підхід, використаний у методі Гюна, використовується для побудови так званих методів прогнозу та корекції, які будуть розглянуті пізніше.
Приклад:
Проведемо обчислення рівняння () з допомогою методу Гюна.
При кроці інтегрування h = 0,1 у першому вузлі сітки x 1 отримаємо:
Що набагато точніше значення, набутого методом Ейлера при тому ж кроці інтегрування. У таблиці 2 наведено нижче порівняльні результати обчислень при h = 0,1 методів Ейлера і Гюна.
Таблиця 2 Рішення рівняння методами Ейлера та Гюна
| x | Точне | Метод Гюна | Метод Ейлера | ||
|---|---|---|---|---|---|
| y | отн. похибка | y | отн. похибка | ||
| 0 | 0,000000 | 0,00000 | 0,00000 | ||
| 0,1 | 0,004837 | 0,00500 | 3,36% | 0,00000 | 100,00% |
| 0,2 | 0,018731 | 0,01903 | 1,57% | 0,01000 | 46,61% |
| 0,3 | 0,040818 | 0,04122 | 0,98% | 0,02900 | 28,95% |
| 0,4 | 0,070320 | 0,07080 | 0,69% | 0,05610 | 20,22% |
| 0,5 | 0,106531 | 0,10708 | 0,51% | 0,09049 | 15,06% |
| 0,6 | 0,148812 | 0,14940 | 0,40% | 0,13144 | 11,67% |
| 0,7 | 0,196585 | 0,19721 | 0,32% | 0,17830 | 9,30% |
| 0,8 | 0,249329 | 0,24998 | 0,26% | 0,23047 | 7,57% |
| 0,9 | 0,306570 | 0,30723 | 0,21% | 0,28742 | 6,25% |
| 1 | 0,367879 | 0,36854 | 0,18% | 0,34868 | 5,22% |
Зазначимо істотне збільшення точності обчислень методу Гюна проти методом Эйлера. Так, для вузла x =0,1 відносне відхилення значення функції, визначеного методом Гюна, виявляється у 30 (!) разів менше. Така ж точність обчислень за формулою Ейлера досягається при числі відрізків інтегрування N приблизно 30. Отже, при використанні методу Гюна при однаковій точності обчислень знадобиться приблизно в 15 разів менше ЕОМ часу, ніж при використанні методу Ейлера.
Перевірка стійкості рішення
Рішення ОДУ у певній точці x i називається стійким, якщо знайдене у цій точці значення функції y iмало змінюється у разі зменшення кроку інтегрування. Для перевірки стійкості, таким чином, треба провести два розрахунки значення ( y i) – з кроком інтегрування h та при зменшеній (наприклад, двоє) величині кроку
Як критерій стійкості можна використовувати трохи відносної зміни отриманого рішення при зменшенні кроку інтегрування (ε – наперед задана мала величина)
Така перевірка може здійснюватись і для всіх рішень на всьому інтервалі значень x. Якщо умова не виконується, крок знову ділиться навпіл і знаходиться нове рішення і т.д. до одержання стійкого рішення.
Методи Рунге-Кутти
Подальше поліпшення точності рішення ОДУ першого порядку можливе за рахунок збільшення точності наближеного обчислення інтеграла у виразі.
Ми вже бачили, яка перевага дає перехід від інтегрування за формулою прямокутників до використання формули трапецій при апроксимації цього інтеграла.
Скориставшись формулою Симпсона, що добре зарекомендувала себе, можна отримати ще більш точну формулу для вирішення завдання Коші для ОДУ першого порядку - широко використовуваного в обчислювальній практиці методу Рунге-Кутти.
Перевагою багатокрокових методів Адамса під час вирішення ОДУ у тому, що у кожному вузлі розраховується лише одне значення правої частини ОДУ - функції F(x,y ). До недоліків можна віднести неможливість старту багатокрокового методу з єдиної початкової точки, так як для обчислень за k -кроковою формулою необхідне знання значення функції в k вузлах. Тому доводиться (k-1) рішення у перших вузлах x 1 , x 2 , …, x k-1 отримувати за допомогою якогось однокрокового методу, наприклад методу
При вирішенні наукових та інженерно-технічних завдань часто буває необхідно математично описати будь-яку динамічну систему. Найкраще це робити у вигляді диференціальних рівнянь ( ДК) або системи диференціальних рівнянь. Найбільш часто вони таке завдання виникає при вирішенні проблем, пов'язаних з моделюванням кінетики хімічних реакцій та різних явищ перенесення (тепла, маси, імпульсу) – теплообміну, перемішування, сушіння, адсорбції, при описі руху макро- та мікрочастинок.
У ряді випадків диференціальне рівняння можна перетворити на вид, у якому старша похідна виражена у явному вигляді. Така форма запису називається рівнянням, дозволеним щодо старшої похідної (при цьому у правій частині рівняння старша похідна відсутня):
Рішенням звичайного диференціального рівняння називається така функція y(x), яка за будь-яких х задовольняє цього рівняння в певному кінцевому або нескінченному інтервалі. Процес розв'язання диференціального рівняння називають інтегруванням диференціального рівняння.
Історично першим і найпростішим способом чисельного розв'язання задачі Коші для ОДУ першого порядку є метод Ейлера. В його основі лежить апроксимація похідної ставленням кінцевих прирощень залежної (y) та незалежної (x) змінних між вузлами рівномірної сітки:
де y i+1 це потрібне значення функції в точці x i+1 .
Точність методу Ейлера можна підвищити, якщо скористатися для апроксимації інтеграла точнішою формулою інтегрування – формулою трапецій.
Ця формула виявляється неявною щодо y i+1 (це значення є і в лівій і в правій частині виразу), тобто є рівнянням щодо y i+1 , вирішувати яке можна, наприклад, чисельно, застосовуючи будь-який ітераційний метод (у такому його можна розглядати як ітераційну формулу методу простої ітерації).
Склад курсової: Курсова робота складається з трьох частин. У першій частині короткий опис методів. У другій частині постановка та вирішення задачі. У третій частині – програмна реалізація мовою ЕОМ
Мета курсової роботи: вивчити два методи розв'язання диференціальних рівнянь-метод Ейлера-Коші та вдосконалений метод Ейлера.
Чисельне диференціювання
Диференціальним називається рівняння, що містить один або кілька похідних. Залежно кількості не залежних змінних, диференціальні рівняння діляться на дві категорії.
Звичайні диференціальні рівняння (ОДП)
Диференціальні рівняння у приватних похідних.
p align="justify"> Звичайними диференціальними рівняннями називаються такі рівняння, які містять одну або кілька похідних від шуканої функції . Їх можна записати вигляді
незалежна змінна
Найвищий порядок, що входить до рівняння (1), називається порядком диференціального рівняння.
Найпростішим (лінійним) ОДУ є рівняння (1) порядку дозволене щодо похідної
Рішенням диференціального рівняння (1) називається будь-яка функція, яка після її підстановки до рівняння перетворює його на тотожність.
Основне завдання, пов'язане з лінійною ОДУ, відомо як завдання Каші:
Знайти рішення рівняння (2) у вигляді функції, що задовольняє початковій умові (3)
Геометрично це означає, що потрібно знайти інтегральну криву, яка проходить через точку ) під час виконання рівності (2).
Чисельний з погляду завдання Каші означає: потрібно побудувати таблицю значень функції, що задовольняє рівняння (2) і початкова умова (3) на відрізку з деяким кроком . Зазвичай вважається, що початкова умова задано в лівому кінці відрізка.
Найпростішим із чисельних методів розв'язання диференціального рівняння є метод Ейлера. В його основі лежить ідея графічної побудови рішення диференціального рівняння, проте цей метод дає одночасно і спосіб знаходження функції в чисельній формі або таблиці.
Нехай дано рівняння (2) з початковою умовою тобто поставлено завдання Каші. Вирішимо спочатку наступне завдання. Знайти найпростішим способом наближене значення рішення в деякій точці де досить малий крок. Рівняння (2) спільно з початковою умовою (3) задають напрямок дотичної шуканої інтегральної кривої в точці з координатами
Рівняння дотичної має вигляд
Рухаючись уздовж цієї дотичної, отримаємо наближене значення рішення у точці :
Маючи в своєму розпорядженні наближене рішення в точці можна повторити описану раніше процедуру: побудувати пряму проходить через цю точку з кутовим коефіцієнтом , і по ній знайти наближене значення рішення в точці
. Зауважимо, що ця пряма не є дотичною до реальної інтегральної кривої, оскільки точка нам не доступна, проте якщо досить мало одержувані наближені будуть близькі до точних значень рішення.
Продовжуючи цю ідею, побудуємо систему і віддалених точок
Отримання таблиці значень шуканої функції
за методом Ейлера полягає у циклічному застосування формули
Малюнок 1. Графічна інтерпретація методу Ейлера
p align="justify"> Методи чисельного інтегрування диференціальних рівнянь, в яких рішення виходять від одного вузла до іншого, називаються покроковими. Метод Ейлера – найпростіший представник покрокових методів. Особливістю будь-якого покрокового методу і те, що з другого кроку вихідне значення у формулі (5) саме є наближеним, тобто похибка кожному наступному кроці систематично зростає. Найбільш використовуваним методом оцінки точності покрокових методів наближеного чисельного рішення ОДУ є спосіб подвійного проходження заданого відрізка з кроком та з кроком
1.1 Удосконалений метод Ейлера
Основна ідея цього методу: обчислюване за формулою (5) чергове значення буде точніше, якщо значення похідної, тобто кутовий коефіцієнт прямої заміщає інтегральну криву на відрізку обчислюватиметься не по лівому краю (тобто в точці), а по центру відрізка . Але так як значення похідної між точками не обчислюється, то перейдемо до здвоєних ділянок центром, в яких є точка, при цьому рівняння прямої набуває вигляду:
А формула (5) набуває вигляду
Формула (7) застосована тільки для , отже, значення нею отримати не можна, тому знаходять за методом Ейлера, при цьому для отримання більш точного результату надходять так: з початку за формулою (5) знаходять значення
(8)
У точці а потім знаходиться за формулою (7) з кроком
(9)
Після того як знайдено подальші обчислення при 
проводиться за формулою (7)
