Диференційних рівнянь

Лекція 8

Багатокрокові методи розв’язання

диференційних рівнянь

Методи Ейлера, Гюна, Тейлора, Рунге-Кута називаються однокроковими методами, тому що в них використовується лише інформація про значення диференційного рівняння в одній попередній точці для того, щоб обрахувати значення в наступній, тобто лише значення в початковій точці (t ₀; y ₀) використовується для обчислення значення диференційного рівняння в точці (t ₁; y ₁), і загалом, щоб обрахувати у_{к +1} необхідне лише у_к. Але після знаходження значень диференційного рівняння в декількох точках можна їх використовувати для отримання більш точних результатів у наступній вузловій точці. Наприклад, для чотирикрокового методу Адамса-Бешфорса-Маултона необхідні значення в точках у_{к -3}, у_{к -2}, у_{к -1}, та у_к для обчислення у_{к +1}. У цьому методі значення диференційного рівняння в чотирьох початкових точках, (t ₀, y ₀), (t ₁, y ₁), (t ₂, y ₂) та (t ₃, y ₃), розраховуються спочатку для подальшої генерації точок {(t _k, y _k): k 4} з використанням багатокрокового методу.

Корисною властивістю багатокрокового методу є можливість визначити локальну помилку відсікання (ЛПВ) та ввести корегуючий елемент, який підвищує точність відповіді на кожному кроці. Також можна визначити, чи буде довжина кроку достатньою, щоб отримати значення у_к+1 із заданим рівнем точності, або задати більший крок, який виключить непотрібні обрахування і забезпечить необхідну точність розрахунків.

 

Метод Адамса-Бешфорса-Маултона. Метод прогнозу-корекції Адамса-Бешфорса-Маултона – це багатокроковий метод, виведений із фундаментальної теореми аналізу:

. (8.1)

Прогноз використовує наближення поліномом Лагранжа для функції f (t, y (t)), що побудована по точках (t_{k -3}, y_{k -3}), (t_{k -2}, y_{k -2}), (t_{k -1}, y_{k -1}) i (t_k, y_k). Для функції f (t, y (t)) з формули (8.1), яка інтегрується на інтервалі [ t_k, t_{k +1}] прогноз Адамса-Бешфорса, має вигляд

. (8.2)

Як будь-який прогноз, він потребує корекції. Коректор отримується аналогічно. А саме, як тільки значення p_{k +1} обраховано, його використовують для побудови наступного поліному Лагранжа для функції f (t, y (t)), який будується по точкам (t_{k -2}; y_{k -2}), (t_{k -1}; y_{k -1}) (t_k; y_k) і новій точці (t_{k +1}; y_{k +1}), в якій замість y_{k +1} використовується прогноз p_{k +1}, тобто функція в точці t_{k +1}; обраховується як f (t_{k +1}, p_{k +1}). Після розрахунку поліному шляхом інтегрування на інтервалі [ t_k; t_{k +1}] отримаємо коректор Адамса-Маултона:

. (8.3)

 

Оцінка помилки і корекція. Залишковий елемент формули чисельного інтегрування використовується, щоб отримати і прогноз, і коректор порядку О(h ⁵). Локальна помилка відсікання для формул (8.2) та (8.3) має вигляд:

, (локальна похибка відсікання для прогнозу) (8.4)

. (локальна похибка відсікання для коректора) (8.5)

Припустимо, що h мале і y ⁽⁵⁾ (t) – майже усталена на інтервалі. Тоді можна виключити члени, що містять похідну 5-го порядку у формулах (8.4) та (8.5). В результаті отримаємо:

. (8.6)

 

Переваги методу прогнозу-корекції Адамса-Бешфорса-Маултона у тому, що формула (8.6) дає наближену оцінку помилки, яку отримуємо при обчисленні значень p_{k +1} i y_{k +1}, не використовуючи у ⁽⁵⁾(t).

Коректор (8.3) використовує наближення при обчисленні у_{к +1}. Оскільки у_{к +1} також є оцінкою для у (t_{k +1}), її можна використовувати у коректорі (8.3) для генерування нового наближення для f_{k +1}, яке, в свою чергу, буде генерувати нове значення для у_{к +1}. Окрім того, якщо продовжити цю ітерацію у коректорі, вона буде збігатись до фіксованої точки у формулі (8.3) швидше, ніж при класичному розв‘язку диференційного рівняння. Якщо необхідно підвищити точність, то процедура корекції більш ефективна, ніж зменшення довжини кроку.

 

Формулу (8.6.) можна використовувати для покращення значень прогнозу. Якщо припустити, що різниця між значеннями прогнозу і корекції на кожному кроці змінюється повільно, то можна підставити і замість і в (8.6) і отримати наступний управляючий параметр:

(8.7)

Його значення використовується замість на кроці корекції, і формула (8.3) приймає вигляд:

 

(8.8)

 

Формулу (8.6) також можна використовувати для визначення адаптивної довжини кроку. Можливо зменшити довжину кроку, наприклад, до h /2 і збільшити до 2 h. Припустимо ε = 5x10^-6 – критерій відносної похибки, і нехай kε = 10^-5

 

Якщо , тоді h = h /2 (8.9)

Якщо , тоді h =2 h (8.10)

 

Якщо прогноз і корекція значень не дають відповідностей розв‘язків у п’ять значущих цифр, то за формулою (8.9) зменшують довжину кроку. Якщо ж розв‘язки збігаються до семи або більше значущих цифр, то за формулою (8.10) довжина кроку збільшується.

Зменшення довжини кроку потребує чотирьох нових початкових значень. Для заміщення значень, яких не вистачає і які ділять інтервали [ t_{k -2}, t_{k -1}] і [ t_{k -1}, t_k ] навпіл, використаємо інтерполяцію функції f (t,y (t)) поліномом четвертого ступеню. Нові чотири вузлові точки, t_{k -3/2}, t_{k -1}, t_{k -1/2} і t_k, використовуються в наступних розрахунках.

Інтерполяційна формула, необхідна для отримання нових початкових значень (рис. 8.1) для кроку h/2 має вигляд:

(8.11)

Рис. 8.1. Зменшення довжини кроку до h /2.

 

Набагато простіше збільшити довжину кроку, тоді для розрахунків за одно-кроковим методом потрібно розрахувати сім попередніх точок (рис. 8.2), щоб подвоїти крок.

Рис. 8.2. Збільшення довжини кроку до 2 ∙ h.

 

 

Метод Мілна-Сімпсона. Інший розповсюджений метод прогнозу-корекції це метод Мілна-Сімпсона. Його прогноз базується на інтегруванні функції f (t, y (t)) на інтервалі [ t_{k -3}, f_{k +1}]:

(8.12)

Прогноз використовує наближення поліномом Лагранжа для f(t, y(t)), побудований по точкам (t_{k -3}, f_{k -3}), (t_{k -2}, f_{k -2}), (t_{k -1}, f_{k -1}) та (t_k, f_k). Інтегруємо його по інтервалу [ t_{k -3}, f_{k +1}] і отримуємо прогноз Мілна:

(8.13)

 

Коректор отримуємо аналогічно. На цьому етапі використаємо значення p_{k +1}. Будуємо інший поліном Лагранжа для функції f (t, y (t)) по точках (t_{k -1}, f_{k- 1}), (t_k, f_k) і у новій точці (t_{k +1}, f_{k +1}) = (t_{k +1, f} (t_{k +1}, p_{k +1})). Інтегруємо поліномом на відрізку [ t _{k -1}, t_{k +1}] і в результаті отримаємо відому формулу Сімпсона:

 

(3.14)

 

Оцінка помилка та корекція. Залишковий член формули чисельного інтегрування використовується для отримання як прогнозу, так і коректору порядку О(h ⁵). Остаточна загальна помилка (О.З.П.) для формул (8.13) та (8.14) дорівнює

(О.З.П. для прогнозу) (8.15)

(О.З.П. для коректора) (8.16)

Припустимо, що крок h настільки малий, що y ⁽⁵⁾(t) майже стала на інтервалі . Тоді із формул (8.15) і (8.16) можемо виключити елементи, які містять похідну п’ятого порядкую, і в результаті отримати формулу:

(8.17)

 

Формула (3.17) дає оцінку прогнозу, яка ґрунтується на обчисленні значень і і не використовує . Її можна застосовувати для покращення значень прогнозу. Якщо припустити, що різниця між значеннями прогнозу і корекції на кожному кроці змінюється повільно, то можна підставити і замість і в формулу (8.17) і отримати наступний управляючий параметр:

(8.18)

Його значення використовується замість на кроці корекції, і формула (8.14) приймає вигляд:

(8.19)

Таким чином, в покращеному модифікованому методі Мілна-Сімпсона послідовно використовуються наступні формули:

(прогноз)

(управляючий параметр)

(8.20)

(коректор)

 

Метод Хеммінга. Ще один важливий метод для розв’язання задачі Коші y’=f(t,y) з початковою умовою y(a)=y₀ на інтервалі [a, b] – метод Хеммінга. Прогноз, як і в попередніх методах, базується на наближенні поліномом Лагранжа для f(t, y(t)). Він побудований по точкам (t_{k -3}, f_{k -3}), (t_{k -2}, f_{k -2}), (t_{k -1}, f_{k -1}) та (t_k, f_k). При інтегруванні отримуємо прогноз Хеммінга:

(8.21)

Управляючий параметр:

(8.22)

 

Коректор отримуємо, будуючи другий поліном Лагранжа для функції f (t, y (t)) по точках (t_{k -1}, f_{k- 1}), (t_k, f_k) і у новій точці (t_{k +1}, f_{k +1}) = (t_{k +1, f} (t_{k +1}, p_{k +1})). В результаті одержуємо формулу:

(8.23)

 

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

 

(8.24)

 

Методи Адамса-Бешфорса-Мултона і Хемінга дозволяють отримати результати з високим рівнем точності за невелику кількість кроків, тобто швидко збігаються до розв‘язку.

 

Наведемо приклад розв‘язання задачі Коші для рівняння з початковими умовами y(0) = 1 на інтервалі [0;10]. Всі три розглянуті вище багатокрокові методи мають порядок похибки розв‘язку О(h⁵). Побудуємо графіки розв‘язків рівняння при мінімальній (меншій за визначений за критерієм стійкості методу) кількості кроків і при кількості кроків, визначеній за формулами (8.24).

Для методу Адамса-Бешфорса-Маултона ; тому N = 65 (рис. 8.3). *// У прикладі

Для методу Мілна-Сімпсона ; тому N = 110 (рис. 8.4).

Для методу Хемінга ; тому N = 70 (рис. 8.5).

Рис. 8.3. Розв‘язок задачі за методом Адамса-Бешфорса-Маултона.

Рис. 8.4. Розв‘язок задачі за методом Мілна-Сімпсона.

Рис. 8.5. Розв‘язок задачі за методом Хемінга.