69316

ЧИСЕЛЬНЕ ІНТЕГРУВАННЯ ЖОРСТКИХ СИСТЕМ ДИФЕРЕНЦІЙНИХ РІВНЯНЬ. ЧИСЕЛЬНІ МЕТОДИ РОЗВ’ЯЗУВАННЯ КРАЄВИХ ЗАДАЧ

Лекция

Математика и математический анализ

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

Украинкский

2014-10-03

1.14 MB

8 чел.

Лекція 18. ЧИСЕЛЬНЕ ІНТЕГРУВАННЯ ЖОРСТКИХ СИСТЕМ ДИФЕРЕНЦІЙНИХ РІВНЯНЬ. ЧИСЕЛЬНІ МЕТОДИ РОЗВ’ЯЗУВАННЯ КРАЄВИХ ЗАДАЧ.

Жорстка задача Коші. Коефіцієнт жорсткості

При побудові і дослідженні математичних моделей об'єктів для підвищення їх точності й адекватності необхідно враховувати велику кількість факторів і явищ, що неминуче приводить до явища жорсткості і описуючих його жорстких рівнянь. Сутність цього явища визначається необхідністю врахування в поведінці (реакції) проектованої системи складових з великими значеннями похідних, які швидко зменшуються, і повільних складових з малими значеннями похідних. Добре відомо, що необґрунтоване неврахування різного роду «малих величин» при математичному моделюванні реальних об'єктів і процесів може істотно спотворити правдиву картину явища. Тому в математичних моделях, що описують ще не вивчені об'єкти і процеси, дослідники змушені враховувати велику кількість, на перший погляд, другорядних факторів. У результаті, як правило, підвищується порядок системи рівнянь моделі і їх жорсткість.

Жорсткість є властивістю математичної задачі, а не застосовуваного чисельного методу.

Розглянемо спочатку лінійну систему диференціальних рівнянь (10.1) з постійною, тобто не залежною від t, матрицею A:

. (10.1)

Нехай {1,2,…, m}, k = k ik — множина власних чисел матриці A. Система диференціальних рівнянь (10.1) з постійною матрицею mхm називається жорсткою, якщо система асимптотично стійка по Ляпунову, тобто :

а відношення найбільшої за модулем дійсної частини власного числа до найменшої за модулем дійсної частини власного числа, тобто:

К =  (10.2)

є великим.

У випадку лінійної системи із змінними коефіцієнтами, тобто при A = A(t) в рівнянні (10.1), власні значення k і, відповідно, число жорсткості K, є функціями t, і визначення жорсткої системи може бути переформульоване в такий спосіб: система

 (10.3)

називається жорсткою на інтервалі (0,T), якщо для всіх t(0,T) виконується умова:

і число  велике.

В основі визначення жорсткості нелінійної системи

 (10.4)

лежить останнє визначення лінійної системи із змінними коефіцієнтами, роль матриці A(t) якої надається матриці Якобі , тобто матриці J, елементи якої рівні частинним похідним:

,

де fj, xj, — компоненти вектор–функцій F, x відповідно. Таким чином, жорсткість визначається властивостями матриці Якобі для правої частини розв'язуваної системи рівнянь (10.4).

Якщо об'єкт дослідження фізично стійкий то:

(10.5)

і модулі розв’язку  будуть згасати при зростанні t із швидкістю, пропорційною 1/Re( - λk), що часто називають локальною «часовою постійною» системи. Мірою жорсткості задачі Коші є інтервал, в якому знаходяться її локальні «часові постійні». Є повна аналогія, наприклад, з числом обумовленості матриці А, яка характеризує погано обумовлені лінійні задачі:

(10.6)

Звичайно кажуть про жорсткість задачі Коші (10.1), якщо К(J)>>1. Однак, у кожному конкретному випадку різні величини К(J) можуть вважатися великими.

З іншого боку, максимальний крок обчислень явними методами обмежений нерівністю:

. (10.7)

Для використання виразів (10.7) і (10.6) необхідні максимальні значення модулів власних λi або максимальні і мінімальні значення їх дійсних частин Re(λi).

Величину λi можна оцінити зверху величиною норми матриці, а їх дійсну частину — величиною сліду матриці Якобі (10.2). Малі за модулем λi , щоб не знаходити потім оберненої матриці J - 1(x), можна також оцінювати інтервалом інтегрування T задачі Коші (10.1), обраного звичайно, виходячи із знання фізичних особливостей досліджуваного об'єкта чи процесу.

Жорсткі системи диференціальних рівнянь з'являються при математичному моделюванні там, де розкид часових характеристик закладено в самій фізичній природі задачі. Прикладами можуть бути опис складних систем шляхом формального об'єднання різнотипних підсистем (наприклад, інерційний об'єкт із малоінерційним контролером), коли вектор розв’язку містить «швидкі» і «повільні» компоненти. Однак, подібні явища можуть проявлятися і як результат об'єднання і взаємодії, у загальному випадку, нежорстких вихідних підсистем, при цьому поведінка системи (її рух) істотно відрізняється від поведінку (руху) кожної підсистеми окремо.

Вперше властивість жорсткості диференціальних рівнянь було відзначено в зв'язку з обчислювальними труднощами їх чисельного інтегрування традиційними явними методами. Розв’язок з малим кроком, обумовленим співвідношенням (10.7), приводить до того, що по завершенні «швидких» складових розв’язку, коли значення x(t) починають мало змінюватися від кроку до кроку, будь–яка спроба збільшення h викликає експоненційне зростання похибки розв’язку через втрату його стійкості. Тому були розроблені спеціальні неявні методи чисельного розв’язку жорстких рівнянь, зацікавість до яких і область застосування яких безупинно зростають. Основна задача таких методів: забезпечувати можливість зміни кроку обчислень у широких межах і адаптація його значень до інформації, що міститься на різних часових ділянках розв’язку.

Жорстку задачу іноді називають «задачею з дуже рознесеними часовими сталими» чи «задачею з великою константою Ліпшиця», під якою розуміють норму матриці Якобі для правої частини рівняння (10.1):

. (10.8)

Задачу можна вважати жорсткою, якщо коефіцієнт жорсткості (10.6) задовольняє нерівність K(J) > 10. Однак, у практичних задачах, властивих моделюванню процесів керування, теорії електричних і електронних ланцюгів, хімічній кінетиці та ін., коефіцієнт жорсткості часто досягає значень 106 - 108.

10.2. Неявні методи Ейлера і Рунге–Кутта

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

Наприклад, для неявного методу Ейлера замість виразу (8.10) буде справедливим наступне:

 (10.9)

де f(un + 1, tn) — оцінка правої частини розв'язуваного диференціального рівняння (10.1) у момент часу tn + 1 , для якого власне і шукається значення un + 1. Формула (10.9) для згаданого випадку має наочну графічну інтерпретацію (рис.10.1), аналогічну до тієї, яка була приведена на рис.10.1 для явного методу.

Рис. 10.1. Графічна інтерпретація неявного методу Ейлера

На рис. 10.1 показано процес виконання чергового кроку наближення, коли з початкової точки з координатами (tn, un) проводиться пряма з нахилом, обумовленим нахилом дотичної до точного розв’язку в точці (tn + 1, xn + 1), у результаті чого знаходиться наступна точка наближення (tn + 1, un + 1). Слід зазначити, що для опуклої форми кривої точного розв’язку наближення неявним методом Ейлера будуть відставати від точних значень функцій, а наближення явним методом, що використовують нахил дотичної у вихідній точці, — перевищувати точне значення. Тому природним є бажання усереднити наближення для отримання більш точного результату:

 (10.10)

Формула (10.10) відома як неявний метод трапецій, і вона деякою мірою відповідає формулі (8.35) явного методу Хейна, в якому точне значення похідної в кінці інтервалу обчислень замінено на наближення через параметр k2.

Аналогії можна продовжити, і, користаючись методикою Рунге–Кутта для апроксимації скороченого ряду Тейлора для f(xn + 1, tn + 1 - h), отримаємо неявні методи Рунге–Кутта:

1)другого порядку:

(10.11)

де k1* = f(tn + 1,un + 1), k2* = f(tn + 1 + h,un + 1 - hk1*).

2)третього порядку:

 (10.12)

де k1* = f(xn + 1,tn + 1), k2* = f(tn + 1 - h/3,un + 1 - h/3k1*), k3* = f(tn + 1 - 2/3h,un + 1 - 2/3hk2*),

3) четвертого порядку:

(10.13)

де k1* = f(xn + 1,tn + 1), k2* = f(tn + 1 - h/2,un + 1 - h/2k1*), k3* = f(tn + 1 - h/2,un + 1 - h/2k2*), k4* = f(tn + 1,un + 1 - hk3*),

Порівнюючи вирази для неявних методів (10.4)–(10.13) з формулами явних методів (8.24)–(8.30), неважко побачити їх велика подібність, за винятком того, що точки початку і кінця інтервалу помінялися місцями, і відлік часу ведеться в оберненому напрямку (назад). Користаючись принципом дуальності, можна за аналогією з формулами (8.23)–(8.30) записувати формули неявних вкладених методів. Як і раніше, глобальна похибка неявного методу k–го порядку оцінюється величиною O(hk), а локальна — величиною O(hk + 1).

Застосування неявних методів Ейлера та Рунге–Кутта перетворює задачу інтегрування жорстких диференціальних рівнянь у задачу розв’язку нелінійних алгебраїчних рівнянь (10.9), (10.11)–(10.13) відносно шуканого наближення функції un + 1 у кожній часовій точці. При цьому складність розв'язуваного нелінійного алгебраїчного рівняння зростає з ростом числа оцінок ki у різних формулах.

Приклад 10.1.

Розв’язати неявним методом Ейлера наступне рівняння: , t0 = 0, x0 = 1, розв’язане раніше в прикладі 8.1 явним методом Ейлера. Потрібно, як і раніше, з кроком h = 0,1 знайти значення х(0,2).

Скориставшись формулою (10.9) , запишем:

а потім послідовно формуємо і розв’язуємо нелінійні алгебраїчні рівняння:

,

,

10.3. Неявні лінійні багатокрокові методи

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

10.3.1. Неявні багатокрокові методи Адамса–Мултона

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

1)неявний метод другого порядку:

який співпадає з неявним методом трапецій (10.10);

2)неявний метод третього порядку: 

 (10.14)

де fn + i, i = - 1,0,1 оцінки правої частини розв'язуваного диференціального рівняння в моменти часу tn + 1, tn, tn - 1.

3) неявний метод четвертого порядку:

 (10.15)

який (при тій же точності) значно простіше для реалізації в порівнянні з неявним методом Рунге–Кутта четвертого порядку (10.13), оскільки тільки один член праворуч fn + 1 залежить від шуканого значення un + 1.

10.3.2. Неявні багатокрокові методи Ракитського Ю.І.

На відміну від методів Рунге–Кутта не вимагають багаторазово оцінювати значення похідної:

1.метод другого порядку:

(10.16)

2.метод четвертого порядку:

 (10.17)

10.3.3. Неявні багатокрокові методи «диференціювання назад» Куртиса–Гіра

На відміну від неявних методів Адамса–Мултона використовують сукупність попередніх значень функції і тільки одну оцінку похідної в поточний момент часу tn + 1, тому що базуються на застосуванні інтерполяційного поліному Лежандра. При постійному кроці обчислень h отримуємо:

1. метод другого порядку:

 (10.18)

2. метод третього порядку:

 (10.19)

3. метод четвертого порядку:

 (10.20)

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

10.4. Багатокрокові неявні методи змінного порядку і змінного кроку

Узагальнюючи формули (10.18)–(10.20), можна записати загальний вираз для неявного методу «диференціювання назад», відомого як формула Гіра із змінним кроком h та змінним порядком k:

 (10.21)

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

 (10.22)

Коефіцієнти αi у формулах (10.21) і (10.22), як буде показано нижче, визначаються

з урахуванням зміни кроку в процесі обчислень:

(10.23)

Оскільки вираз (10.21) являє собою нелінійне алгебраїчне рівняння (невідома величина un + 1 знаходиться в його лівій і правій частинах одночасно), то для поліпшення збіжності методу Ньютона–Рафсона визначається початкове наближення у вигляді прогнозованого значення:

, (10.24)

де коефіцієнти

(10.25)

Відома формула диференціювання назад Брайтона (ФДН) відрізняється від формули Гіра (10.21) тим, що замість повних значень функції використовуються тільки їх скінчені різниці:

 (10.26)

Тоді замість виразів (10.22) і (10.24) отримуємо еквівалентні співвідношення:

(10.27)

у яких коефіцієнти взаємозалежні:

(10.28)

Якщо піти далі і від перших скінчених різниць значень функції перейти до скінчених різниць високих порядків, то отримаємо нові еквівалентні співвідношення:

, (10.29)

 (10.30)

у яких коефіцієнти і скінчені різниці високих порядків обчислюються за наступними рекурентними формулами:

 (10.31)

10.5. Визначення коефіцієнтів неявних формул наближення

Виберемо для ілюстрації методики визначення значень коефіцієнтів формулу (10.22):

.

На її основі побудуємо систему лінійних рівнянь для обчислення значень αi, по черзі застосовуючи цю формулу до послідовності спеціально підібраних поліномів:

,

Формула (10.22) точно визначає значення похідних від цих поліномів

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

 (10.32)

Підставляючи обчислені значення (10.32) в розгорнуту формулу (10.22)

будуємо шукану систему лінійних рівнянь відносно невідомих коефіцієнтів:

Застосовуючи матричну форму запису з визначником Ван–дер–Монда, остаточно отримуємо:

, (10.33)

звідки

  (10.34)

Аналогічно будується система лінійних рівнянь для визначення коефіцієнтів формули прогнозу (10.24):

 (10.35)

з якої отримують значення

 (10.36)

які співпадають з (10.25).

Приклад 10.2.

Знайти значення коефіцієнтів формули наближення неявного багатокрокового методу Гіра, який відповідно до виразу (10.22) має вигляд:

Оскільки порядок формули p = 2, то шукана система лінійних рівнянь відносно невідомих коефіцієнтів при постійному кроці обчислень формується у вигляді:

Вона легко приводиться до наступної лінійної системи:

з якої отримуємо α2 = - 1/2, α1 = 2, α0 = - 3/2.

Отже, для розглянутого метода при постійному кроці:

або ,

що збігається з виразом (10.18), відомим як неявний метод Куртиса–Гіра–Шихмана другого порядку.

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

Приклад 10.3.

Знайти значення коефіцієнтів формули наближення неявного багатокрокового методу, який має вигляд:

.

Порядок формули теж p = 2, і шукана система лінійних рівнянь для невідомих коефіцієнтів при постійному кроці обчислень набуває вигляду:

З першого рівняння знаходиться коефіцієнт α1 = 1. Друге і третє рівняння перетворяться в систему рівнянь:

з розв’язку якої знаходимо b1 = b2 = 1, тому шукана формула наближення дорівнює:

.

10.6. Стійкість неявних методів

Визначимо області стійкості неявних методів і порівняємо їх з умовами стійкості явних методів, розглянутих у попередніх главах 8 і 9. Почнемо з неявних методів Ейлера і трапецій, скориставшись методикою аналізу стійкості різницевих рівнянь, яка вже була використана в главі 9 і базується на визначенні коренів характеристичного рівняння (9.41).

10.6.1. Неявний метод Ейлера

На підставі виразу (10.9) з урахуванням тестового рівняння x’ = λx формуємо різницеве рівняння методу un + 1 = xn + hλun + 1 звідки (1 - hλ)un + 1 - un = 0. Відповідне йому характеристичне рівняння набуває вигляду (1 - hλ)z - 1 = 0 корінь якого для Reλ < 0 і довільного значення h завжди задовольняє умову стійкості:

(10.37)

10.6.2. Неявний метод трапецій

На основі виразу (10.10) з урахуванням тестового рівняння x’ = λx формуємо різницеве рівняння методу:

 чи

Відповідне характеристичне рівняння отримуємо у вигляді:

з якого, як і раніше для неявного методу Ейлера, випливає, що для Reλ < 0 і довільного значення h умова стійкості завжди задовольняється:

 (10.38)

10.6.3. Неявний метод Куртиса–Гіра–Шихмана другого порядку

На підставі виразу (10.18) 

з урахуванням тестового рівняння x’ = λx формуємо різницеве рівняння методу:

характеристичне рівняння для якого:

Оцінимо корені цього рівняння:

, (10.39)

тому що при < 0 знаменник кореня завжди більше одиниці за значенням (для ||≤1/2), або за модулем (для || > 1/2).

Запропоновано називати А — стійкими розглянуті тут та інші неявні методи, в яких область стійкості при розв’язку тестового модельного рівняння: x’ = λx містить всю ліву напівплощину площини , тобто при Reλ<0 розв’язок буде асимптотично стійким при будь–якому позитивному h. Було також показано, що ніякий явний лінійний багатокроковий метод не може бути А-стійким, і що не існує А-стійкого неявного лінійного багатокрокового методу з порядком k>2.

Дослідження стійкості методу Куртиса–Гіра–Шихмана третього порядку зручніше проводити графічно за допомогою побудови годографа стійкості аналогічно визначенню областей стійкості багатокрокових методів у главі 9 (див рис. 9.2).

Записуючи характеристичне рівняння (10.19) у вигляді: (11 - 6) - 18z - 1 + 9z - 2 - 2z - 3 = 0

і знаходимо з нього μ = = (11 - 18z - 1 + 9z - 2 - 2z - 3)/6.

Підставляючи z = q = e,отримаємо рівняння межі області стійкості:
μ = hλ = (11 - 18e + 9e2 - 2e3)/6

Рис. 10.2. Область стійкості трьохкрокового методу Гіра

На відміну від попереднього прикладу, існують точки межі області стійкості, розташовані в лівій напівплощині. Тому чотирикроковий метод не є А — стійким

Для неявних методів з порядком k>2 вимоги А — стійкості послаблюються. Вводиться поняття «жорстко стійких» методів (рис. 10.3), для яких їх область стійкості при розв’язку «тестового» рівняння містить не всю ліву напівплощину , а тільки дві області R1 і R2, де R1 визначається умовою:

- aRe(λh)b, ||Im()||≤d (10.40)

а R2  умовою: 

Re(λh) < - a, a > 0, d > 0, b 0 (10.41)

Рис. 10.3. До визначення жорстко–стійких методів

Крім того, для hλR1, забезпечується задана точність. Суть вимоги R1R2 полягає в тому, що при великих h|λ| забезпечується загасання «швидких» складових розв’язку, а при малих h|λ|, якщо R1, зберігається необхідна точність «повільних» складових розв’язку.

До жорстко–стійких методів відносяться методи Гіра диференціювання назад, неявні методи Ракитського Ю.В. і неявні методи Рунге–Кутта при k>2.

Приведемо для прикладу умови стійкості неявних методів Рунге–Кутта:

1. другого порядку:

 (10.42)

2. третього порядку:

 (10.43)

3. четвертого порядку:

 (10.44)

Рис. 10.4. Області стійкості неявних методів Рунге–Кутта

Порівнюючи умови (10.42)–(10.44) з умовами стійкості явних методів тих же порядків (8.59), знаходимо, що умови стійкості неявних методів Рунге–Кутта можуть бути побудовані з областей стійкості явних методів (рис.8.7) заміною на — і зміною знаку нерівності. Однак при цьому змінюється на протилежне правило оцінки властивостей будь–які точки площини: тепер точки всередині годографа характеризують нестійкі розв’язки, а поза годографом стійкі. Області стійкості багатокрокових неявних методів змінного порядку Гіра, отримані на підставі виразу (10.21), приведені на рис.10.5.

Як і у випадку рис. 10.4, точки всередині годографа характеризують нестійкі розв’язки, а поза годографом стійкі. З приведеного рисунку видно, що методи першого і другого порядку (k  2) є А — стійкими, тобто області їх абсолютної стійкості включають всю ліву напівплощину Re(h ) 0 , де власні значення матриці Якобі для правої частини розв'язуваного диференціального рівняння. Для k = 3,…,6 вони, з одного боку, А(α) стійкі, тому що області їх абсолютної стійкості включають необмежений клин arg( - h ) < α, а з іншого боку — жорстко стійкі в сформульованому вище контексті, тому що області їх стійкості обмежені частиною лівої напівплощини , при цьому в околі початку координат забезпечується не тільки стійкість, але і точність.

Рис. 10.5. Області стійкості k–крокових методів Гіра а — годографи методів 1–го — 4–го порядків; б — годографи методів 4–го — 6–го порядків

Жорстку стійкість зручно характеризувати максимальним кутом дотичної α до кривої годографа в лівій напівплощині (для k2). У табл. 10.1 приведені значення кута α для формул Гіра різних порядків k  6 при відносній зміні кроку обчислень q = hj / hj - 1.

Таблиця 10.1 Значення кута α для годографів методів Гіра різних порядків

Зміна кроку q

Порядок методу k

1

2

3

4

5

6

7

1

90o

90o

86o

74o

52o

18o

0o

1.01

90o

90o

85o

72o

48o

10o

0o

1.02

90o

90o

85o

71o

45o

1o

0o

1.05

90o

90o

83o

65o

32o

0o

0o

1.1

90o

89o

81o

55o

1o

0o

0o

1.2

90o

87o

70o

25o

0o

0o

0o

1.5

90o

75o

25o

0o

0o

0o

0o

2.0

90o

50o

0o

0o

0o

0o

0o

З табл. 10.1 випливає, що в процесі обчислень часовий крок потрібно змінювати поступово (а не подвоювати або зменшувати вдвічі, як це прийнято), якщо ми хочемо використовувати при обчисленнях неявні методи високих порядків. Підвищити порядок А(α) стійких методів можна, якщо у формулу наближення ввести другу похідну, як це зроблено в методах Енрайта:

(10.45)

при цьому жорстка стійкість зберігається аж до методу з порядком k = 12 (рис.10.6). Друга похідна функції, яка потрапила до формули (10.45), обчислюється через матрицю Якобі для правої частини розв'язуваного рівняння так:

Рис. 10.6. Області стійкості методів Енрайта, які використовують другу похідну

10.7. Оцінка локальної похибки і організація обчислень

Як було відзначено в главі 9, головний член локальної похибки наближення багатокрокового лінійного методу k–го прядку визначається першим неврахованим членом відповідного еквівалентного ряду Тейлора, який дорівнює:

= [Ck + 1 x (k + 1) ( ) h k + 1].

Цей вираз, що містить похідні розв’язку високих порядків, можна безпосередньо використовувати у випадках, коли застосовуються формули методів Гіра змінного порядку і змінного кроку (10.29) і (10.30), виражені через кінцеві різниці високих порядків, при цьому:

un + 1 =  (un + 1 - un). (10.46)

В цих випадках принципово можлива одночасна паралельна оцінка похибки обчислень для кожного із стійких методів з порядками k = 1,…,6, що полегшує вибір найбільш ефективного з них на кожному кроці обчислень.

Якщо ж застосовуються формули методів Гіра (10.22) або (10.27), виражені через значення функції або її перші скінчені різниці, то доречніше використовувати оцінку похибки через «вкладені методи» (див. параграф 8.4). При цьому найбільш зручним є використання формул прогнозу (10.24), (10.27,б) і (10.30) , з одного боку, і формул наближень (10.22), (10.27,а) і (10.30), з іншого. В результаті похибка виражається через різницю прогнозованого і обчисленого значень:

 (10.47)

Оцінка похибки через розв’язок з різними кроками, описана в главі 8, для методу змінного кроку принципово неможлива.

Нагадаємо, що прогнозоване значення u(0)n + 1 використовується також при застосуванні методу Ньютона для обчислення наближення u(k)n + 1 при розв’язку нелінійних алгебраїчних рівнянь (10.22), (10.27,а) і (10.30). З метою економії часу часто застосовується модифікована процедура Ньютона, для якої обчислення матриці Якобі і розв’язок отриманої лінійної системи методом LU–розкладання виконується лише один раз на кожному часовому кроці. Це, звичайно, веде до зростання локальної похибки і, як наслідок, до зростання числа виконуваних часових кроків на заданому інтервалі часу, однак, загальний час розрахунку може виявитися меншим, тому що істотно зменшуються витрати на кожен такий крок. У табл.. 10.2 ця ситуація ілюструється на декількох практичних прикладах, для яких обчислення прискорюються в 1.5–2 рази при похибці розв’язку  = 10 - 6.

Таблиця 10.2. Порівняльні характеристики розрахунків

Схема

Число рівнянь

Багато ітерацій

Одна ітерація

N

NS

NIT

RJC

NS

NIT

RJC

Генератор

9

1241

3500

128

1616

1812

196

Лінія затримки

22

800

2108

80

824

899

96

Підсилювач

32

203

514

21

210

232

22

Тригер

32

591

1447

73

613

682

69

Схема керування

61

107

450

8

271

307

36

Імпульсна схема

102

790

2183

70

811

895

84

Тут: NS  кількість часових кроків; NITзагальна кількість ітерацій; RJCкількість відкинутих кроків.

Звичайно, така організація обчислень можлива за умов, які практично збігаються з загальними умовами можливості розв'язання диференціального рівняння (див. параграф 8.1) і полягають в тому, що права частина диференціального рівняння, принаймні, повинна бути двічі диференційована, і її значення повинні бути обмежені. До того ж для всіх t повинна бути монотонною наступна функція:

. (10.48)

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

10.8. Автоматичний вибір часового кроку і порядку методу при обчисленнях

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

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

10.8.1. Алгоритм ALLTED (Україна)

Застосовуваний у випадках, коли використовуються формули методів Гіра змінного порядку і змінного кроку (10.29) і (10.30), виражені через скінені різниці високих порядків. Він передбачає наступні кроки:

1. На кожному часовому кроці за формулами (10.46) прогнозується похибка розв’язку для методів 1, 2,... , k + 1 порядків :

(10.49)

2. Вибирається найменша з можливих похибка:

(10.50)

яка визначає вибір порядку методу.

3. Визначається нове значення часового кроку через питому задану користувачем похибку розв’язку (8.56) з врахуванням обраного прогнозованого мінімального значення локальної похибки (10.50) і залежності похибки розв’язку від порядку методу ( 8.9):

, (10.51)

де — задана користувачем відносна похибка обчислень; T — часовий інтервал розв’язку; hj — попередній крок.

10.8.2 Алгоритм DIFSUB (США)

Запропонований Гіром, і застосовуваний у випадках, коли використовуються формули методів Гіра змінного порядку і змінного кроку (10.22) і (10.27), а похибка визначаається через різницю прогнозованого і обчисленого значень (10.47). Передбачає наступні кроки:

1. Прогнозуються похибки для трьох сусідніх методів:

Після k кроків порядок методу змінюється у відповідності із значенням, що забезпечує мінімальне значення EiD

2. Визначається новий часовий крок за формулою (10.51).

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

Нагадаємо, що коефіцієнти формул наближення методів змінного порядку і змінного кроку, обумовлені формулами (10.23), (10.25), (10.28) і (10.31), обчислюються на кожному часовому кроці, що збільшує обсяг обчислень. Тому була запропонована спрощена схема організації обчислень із змінним кроком, коли використовуються формули наближення з постійними коефіцієнтами, але перерозраховуються шляхом інтерполяції збережені раніше обчислені значення з урахуванням зміненого значення часового кроку. Це досягається збереженням попередніх значень у вигляді вектора обернених скінчених різниць високих порядків:

 (10.52)

або вектора Нордсика:

, (10.53)

У першому випадку старі і нові значення пов'язані через коефіцієнт q, обумовлений виразом (10.51):

або в матрично–векторній формі

. (10.54)

В другому випадку для вектора Нордсика перерозрахунок збережених значень для нового значення кроку спрощується:

, (10.55)

тому що матриця перерозрахунку значень Q є, на відміну від матриці T(q), чисто діагональною Q = diag(qi), ненульові елементи в ній визначаються коефіцієнтом q в i–му степені.

10.9. Спільні явно–неявні процедури розв’язку

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

Доцільно побудувати такі спільні явно–неявні процедури розв’язку, в яких швидко змінювані ділянки розв’язку, для яких значення часових кроків обмежене, обраховувалися б за допомогою явних методів, а малоінформативні маломінливі ділянки — за допомогою неявних методів, які дозволяють довільно змінювати (і збільшувати) часовий крок. Основним питанням для таких процедур є питання автоматичного розпізнання швидких і повільних ділянок розв’язку, і автоматичного вибору відповідного методу розв’язку.

Рис. 10.7. Типова перехідна характеристика з виділеними ділянками кращого використання явних і неявних методів розв’язку

Одним з можливих розв’язків задачі розпізнавання локальної «жорсткості» розв’язку в околі поточної часової точки tn[t0,T] є перевірка виконання умови:

 (10.56)

де

R(Jn) — спектральний радіус матриці Якобі Jn, обумовлений її найбільшим власним значенням λmax, яке знаходиться процедурою ітерацій довільного вектора Δu0 = f(tn,un) матрицею Jn( див. главу 4).

Умова переходу від неявного методу до явного визначається співвідношенням:

 (10.57)

Умова зворотного переходу від явного методу до неявного набуває вигляду:

, (10.58)

де γ співвідношення обчислювальних витрат на виконання кроку неявним і явним методами (γ = 2).

У табл. 10.3 цей підхід ілюструється на декількох практичних прикладах, для яких обчислення прискорюються на 20–30 % при похибці розв’язку  = 10 - 6.

Таблиця 10.3. Результати спільного застосування явного і неявного методів

Схема

N

T, з

NS

NS*

Виграш, %

Формувач адреси

306

20

200

177/37

30,2

Схема синхронізації

388

50

333

281/26

21,5

Формувач імпульсів

548

80

546

500/38

27,4

Підсилювач

783

50

464

395/34

26,1

Контролер

10770

100

824

697/26

22,6

Тут: N — кількість компонент sd; T — інтервал інтегрування, сек; NS — кількість часових кроків базового методу; NS* — співвідношення кількості кроків неявного і явного методів.

При експериментальних обчисленнях, результати яких включені в табл.10.3, отримані при спільному застосуванні неявного методу Гіра з скінченими різницями високих порядків (10.29) і явного методу PECE:

(10.59)

10.10. Явні нелінійні методи розв’язку з А–стійкістю

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

10.10.1. Явний нелінійний метод першого порядку

Відмовимося від експонентної апроксимації розв’язку диференціального рівняння, обумовленого власними значеннями матриці Якобі від його правої частини (10.4), і перейдемо до часової апроксимації на інтервалі [tn,tn - 1] типу,

.

Тоді при k = 1 маємо гіперболічну апроксимуючу функцію:

,

що відповідає наступній нелінійній формулі наближення, яка запропонована Федулою:

(10.60)

Покажемо, користуючись тестовим рівнянням x’ = λx, що при Reλ < 0 формула (10.60) має А — стійкість, яка збігається зі стійкістю неявного методу Ейлера. Дійсно, на підставі співвідношення (10.60) можна записати:

. (10.61)

З виразу (10.61) випливає, що при Reλ < 0 при будь–якому довільному кроці h

. (10.62)

10.10.2. Явний нелінійний метод другого порядку.

Українським математиком Глинським була запропонована наступна явна формула другого порядку для розв’язку жорсткого диференціального рівняння:

(10.63)

Астійкість цієї формули можна довести, розглядаючи дві сусідні ітерації, для яких при використанні тестового модельного рівняння  повинна виконуватися умова:

. (10.64)

Використовуючи модельне рівняння, перетворимо вираз (10.63):

. (10.65)

Проаналізуємо знаменник формули (10.65) і переконаємося в наступному:

1. якщо 1/3<Kn() < 1, то Kn + 1() < 1 і умова (10.64) виконується;

2. якщо Kn() = 1/3, то Kn + 1() = 1 і умова (10.64) все ж виконується;

3. якщо 0<Kn() < 1/3, то при /2(3 - 1/Kn()) = 1/2 досягається мінімальне значення знаменника виразу (10.65), рівне 3/4, у результаті чого Kn + 1max = 4/3 і умова (10.64), як і раніше, виконується.

Приклад 10.2.

Розв’язуючи систему рівнянь

методом Глинського і Адамса–Башфорта другого порядку, порівняємо поведінку розв’язку при обчисленнях з часовим кроком h = 0,05 Неважко бачити, що власні значення матриці Якобі дорівнюють λ1 = - 0,1, λ2 = - 200. Тому припустимий крок для методу Адамса–Башфорта відповідно до умови (9.47) обмежений максимальною величиною h<1/λ2 = 0,005.

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

Таблиця 10.4. Результати оцінки похибки явних лінійного і нелінійного методів

h

t

Помилка методу Глінського

Помилка методу Адамса–Башфорта

0,05

0,15

- 0,22 E - 3

- 0,2 E - 3

0,05

1,15

- 0,22 E - 3

E + 20

0,05

2,15

- 0,21 E - 3

E + 45

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

ВИСНОВКИ:

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

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

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

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

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

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

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

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

ВПРАВИ:

1. Вирішить неявним методом трапецій (10.10) наступне рівняння:яке було розв’язане , , яке було розв’язане раніше в прикладі 10.1 неявним методом Ейлера. Необхідно, як і раніше, з кроком h = 0,1 знайти значення х(0,2).

2. Найдіть умови стійкості формули наближення методу середньої точки:  для випадку 

3. Визначить, чи є серед наведених різницевих рівнянь такі, що забеспечують стійкі рішення:

а)

б)

в)

4. Покажить, що всі рішення різницевого рівняння стійкі для n →∞ , якщо - 1 < λ < 1. Для будьяких інших значень λ в комплексній площині існує, хоч би, одне нестійке рішення.

5. Побудувати годограф стійкості в комплексній площині для неявного методу Гіра , користуючись пакетом Mathematica і характеристичним рівнянням (10.39): .

6. Знайдіть значення коефіцієнтів  и  , яких бракує в формулі 3 - го порядку: , користуючись методикою параграфа 10.5.

7. Отримайте значення коефіцієнтів методів Гіра третього, четвертого и п’ятого порядків для постійного часового крока.

Відповіді:

1. х(0,2) = 0,83438 (точне значення 0,83578)

2. <1

3. Ні , тому що значення коренів для випадку а (3, - 1), для випадку б (2±i) і для
випадку в ( - 1, - 1).

5.

Рис. 10.8. Область стійкості двокрокового методу Гіра

6. ,

7. Формула Гіра 3–го порядку:
.
Формула Г
іра 4–го порядку:
.
Формула Г
іра 5–го порядку:
.

КРАЙОВІ ЗАДАЧІ ДЛЯ ЗВИЧАЙНИХ ДИФЕРЕНЦІАЛЬНИХ РІВНЯНЬ

Звичайні диференціальні рівняння (ЗДР) задаються з деякими допоміжними умовами, що використовуються для визначення констант інтегрування при розв’язанні цих рівнянь. Для рівняння n-го порядку потрібно n таких умов. Якщо всі умови задані для того ж самого значення незалежної змінної (зокрема, для однієї часової точки), то реалізується задача про початкові умови Коші, докладно розглянута а главах 8-10. Навпаки, якщо умови задаються для різних значень незалежної змінної, то реалізується крайова задача, чи задача про межові умови.

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

,  (11.1)

із межовими умовами

 (11.2)

і називаються двоточковими крайовими задачами.

Перш ніж застосовувати будь-який чисельний метод, варто перевірити умови, що гарантують існування розв’язку задачі (11.1), (11.2). Наступна теорема дає загальні умови, що забезпечують існування й одиничність розв’язку. [25].

Теорема 1. Припустимо, що f(x,y,z) неперервна в області

R = {(x,y,z):axb, - ∞ < y < ∞, - ∞ < z < ∞} і, що

 і

неперервні на R. Якщо існує постійна M > 0, для якої задовольняють умовам

fy(x,y,z) > 0 для всіх (x,y,z)R і fz(x,y,z)M для всіх (x,y,z)R (11.3)

то крайова задача (11.1), (11.2) має єдиний розв’язок y(x)для axb.

Найбільш часто зустрічаються і найкраще вивчені двоточкові лінійні крайові задачі наступного вигляду:

,  , (11.4)

, , (11.5)

де , .

Умови, яким повинні задовольнять функції p(x), g(x) і f(x), щоб задача (11.4), (11.5) мали єдиний розв’язок, виходять з теореми як наслідок.

Наслідок. Припустимо, що p(x) і g(x) неперервна на R і g(x) < 0, то задача (11.4), (11.5) має єдиний розв’язок на axb.

Крайові умови (11.5) визначають третю крайову задачу для рівняння (11.4). Якщо
β0 = β1 = 0, то умови (11.5) визначають першу, а при α0 = α1 = 0 другу крайові задачі.

Точне (аналітичне) розв’язання крайових задач викликає більші труднощі, чим розв’язок задач Коші. Звідси — велике розмаїття наближених методів розв’язання таких задач. Ці методи можна розділити на дві групи: наближено-аналітичні, що дають наближене розв’язання крайової задачі на відрізку [a,b] у вигляді конкретної аналітичної функції, і чисельні методи, що визначають розв’язок у вигляді табличної функції, заданої на сітці [a,b]. Нижче будуть розглянуті ідеї і реалізації методів розв’язання обох типів.

11.1. Розв’язання лінійної крайової задачі композицією двох задач Коші

Будемо шукати розв’язок задачі (11.4), (11.5) у вигляді:

y(x) = Au(x) + v(x), (11.6)

де A — деяка константа, u(x) — задовольняє однорідному рівнянню

L[u] = 0 (11.7)

а v(x) — задовольняє неоднорідному рівнянню

L[v] = f(x). (11.8)

У цьому випадку функція y(x), обумовлена формулою, (11.6) буде розв’язком рівняння (11.3) при будь-якому A, в силу його лінійності. Справді

L[y] = L[Au + v] = AL[u] + L[v] = f(x).

Якщо припустити, щоб розв’язок (11.6) задовольняє першій крайовій умові (11.5) при будь-якому A, то отримаємо:

la[y] = la[Au + v] = Ala[u] + la[v].

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

la[y] = α0u(a) + β0u’(a) = 0 (11.9)

і la[y] = α0v(a) + β0v’(a) = γ0 (11.10)

Рівність (11.9) буде виконана, якщо прийняти, наприклад,

u(a) = β0, u(a) = - α0 (11.11)

Щоб задовольнити рівності (11.10), можна прийняти:

, якщо   (11.12)

чи , якщо . (11.13)

Врахуємо, що одночасно α0 і β0 у нуль не обертаються в силу умови (11.5).

Таким чином, для розв’язання крайової задачі (11.4), (11.5) необхідно розв’язати наступні задачі Коші:

,  (11.14)

і  (11.15)

з початковими умовами (11.12) чи (11.13), використовуючи будь-який чисельний метод розв’язку задачі Коші для рівнянь другого порядку. Наближений розв’язок цих задач отримуємо на відрізку [a,b], у результаті чого стають відомими значення u(b), u(a), v(b), v(a) . Це дозволяє підібрати константу A так, щоб функція (11.6) задовольняла не тільки рівнянню (11.12) і першій крайовій умові, але і другій крайовій умові (11.5). Маємо

lb[y] = lb[Au + v] = Alb + lb[v] = γ1,

звідки

, якщо . (11.16)

Якщо α1u(b) + β1u(b) = 0, то це означає, що однорідна крайова задача

L[u] = 0, la[u] = 0, lb[u] = 0

має нетривіальний розв’язок u(x), що слугує ознакою вирожденості вихідної задачі (11.4) , (11.5).

Приклад 11.1.

Розглянемо розв’язок наступної крайової задачі

,

y( -1 ) = y(1) = 0 (11.17)

Порівнюючи умови (11.17) із загальними крайовими умовами бачимо, що

a = - 1, b = 1, α0 = 1, β0 = 0, α1 = 1, β1 = 0, γ0 = 0, γ1 = 0.

Шукаємо розв’язок у вигляді (11.6). Тоді для функції u(x), відповідно до (11.7), маємо наступну задачу Коші

Розв’язуючи будь-яким чисельним методом останню задачу, отримуємо , використовувавши Mathematica-5:

a = -1; b = 1,α0 = 1;β0 = 0; α1 = 1; β1 = 0; γ0 = 0; γ1 = 0;
S1 = NDSolve[{u”[x]-(1+x2)u[x] =  = 0, u[-1] =  = 0, u’[-1] =  = -1}, u[x],
{x, -1, 1}; U[x_] = u[x]/.S1;

отримуємо розв’язок U(x). 

Задача Коші (11.15) з початковими умовами (11.12) для розглянутого випадку має вигляд:

Розв’язуючи аналогічно цю задачу, отримаємо функцію V(x):

S2 = NDSolve[{v”[x]-(1+x2)v[x] =  = x2+2x, v[-1] =  = 0, v’ [-1] =  = 0}, v[x],
{x, -1, 1} V[x_] = v[x]/.S2

Для визначення константи A по формулі (11.16) знайдемо значення похідних U(b), V(b) і обчислимо значення A 

dU[x_] = D[U[x],x]; dV[x_] = D[V[x],x]
A = (γ1-α1V[b]-β1dV[b])/(α1U[b]+β1dU[b])
Out = {-0.250928}

відповідно до рівняння отримаємо розв’язок у вигляді :

Y(x) = - 0,250928U(x) + V(x),

Графік якого приведений на рис. 11.1

Рис. 11.1. Розв’язок крайової задачі (11.17)

11.2. Метод прицілювання для зведення крайової задачі до задач Коші

Викладений вище композиційний метод редукції до задач Коші має наступні істотні недоліки:

1  він не дозволяє використовувати сучасні методи розв’язку задачі Коші змінного порядку і змінного кроку, а вимагає визначення розв’язків u(x) і v(x) з однаковим постійним кроком, щоб скористатися в кожній часовій точці їх композицією (11.6);

2  він призначений тільки для одномірної лінійної задачі, тому що у випадку системи рівнянь буде потрібне не обчислення константи А (11.16), а матриці А по матрично-векторному виразу типу (11.16), що являє собою далеко не просту задачу;

3  він не підходить для розв’язку нелінійних крайових задач.

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

y = f(x,y,y)

із заданими крайовими умовами

y(a) = γ0, y(b) = γ1, x[a,b]

знаходять у такий спосіб. Підбирається початкове значення першої похідної y(a), при якому виконується друга крайова умова y(b) = γ1 , шляхом ітераційного розв’язання задачі Коші у вигляді:

y = (x,y,y) (11.18)

y(a) = γ0 і y(a) = .

Спочатку вибирається довільне значення 0 і розв’язується задача Коші (11.18). Бажано вибирати 0 так, щоб y(b) < γ1 (рис.11.2). Потім вибирається 1 і розв’язок задачі Коші повторюється. Тепер бажано вибрати його так, щоб y(b) > γ1 (рис.11.2).

Рис. 11.2. Ілюстрація методу прицілювання

Після цього інтерполяцією уточнюється значення i для наступних задач Коші з початковими умовами:

2 = 1 - (g(1) - ) ,

3 = 2 - (g(2 - ) і т.д. (11.19)

Такий метод є універсальним і поширюється на випадки нелінійних диференціальних рівнянь n-го порядку. Слід зазначити, що через довільний вибір 0 задача (11.18) може виявитися жорсткою (глава 10), хоча вихідна задача (11.1),(11.2) буде добре обумовленою.

Приклад 11.2.

Розв’язати крайову задачу для нелінійного диференціального рівняння:

y = , y(0) = 1, y(1) = 1,5 

Припустимо, що були знайдені наступні розв’язки задачі Коші:

y(0) = 0  y(1) = 1,543081

y(0) = - 0,05  y(1) = 1,425660.

Тоді після використання процедури інтерполяції (11.19) знаходимо шукане початкове значення y(0) = - 0,01833.

Для прискорення обчислень можна організувати паралельний пошук, виділивши на інтервалі [a,b] n відрізків [xi - 1, xi], так що x0 = a і xn = b. Вводимо змінну t, t[0,1], через яку описуємо вектор змінних y(t) розмірності 2n:

y(t) = [y1(t),…,yp(t),(t),…,(t)],

yi(t) = y(xi - 1 + t(xi - xi - 1)) для кожного відрізка [xi ,xi - 1].

Тоді вихідне диференціальне рівняння (11.1) заміняється системою диференціальних рівнянь першого порядку:

y = (y,t),

а крайові умови для i = 1,2,…,n - 1 задаються умовами “стикування” окремих відрізків:

yi + 1(0) = yi(1),   (11.20)

додатково до крайових умов (11.2).

При застосуванні методу прицілювання початкові значення для складових вектора y(0) вибирають у такий спосіб:

 та

Для наступної k-тої ітерації прицілювання з урахуванням крайових умов (11.20) вибирають .

11.3. Метод кінцевих різниць

Розглянемо розв’язання крайової задачі (11.4), (11.5). Ідея методу полягає в тому, що похідні в диференціальному рівнянні (11.4) і крайових умовах (11.5) заміняються їх кінцево-різницевими апроксимаціями. Введемо на відрізку [a,b] сітку з кроком:

.

Позначимо через yi = y(xi) значення точного розв’язку задачі (11.1) у i-му вузлі сітки, а через ui  значення наближеного розв’язку в цій точці. Заміняючи в кожному внутрішньому вузлі сітки похідні різницями, отримаємо різницеві рівняння :

(11.21)

Симетричним різницевим апроксимаціям похідних першого і другого порядку властива глобальна похибка другого порядку відносно h, тобто О(h2), а отже, локальна похибку третього порядку відносно h, чи О(h3). Це легко доводиться розкладанням у ряд Тейлора виразів:


з різниці яких отримуємо шуканий результат:

.

Тому різницеве рівняння (11.21) апроксимує вихідне диференціальне рівняння (11.4) на його розв’язку також із другим порядком відносно h. Знайдемо нев'язку різницевого рівняння;

. (11.22)

Апроксимуємо крайові умови наступними різницевими рівняннями:

, (11.23)

Знайдемо похибку апроксимації крайових умов. Нев'язки рівнянь (11.23) мають вигляд:

.

Асиметрична апроксимація першої похідної на відміну від симетричної має глобальну похибку першого порядку по h, тобто О(h), і отже, локальну похибку другого порядку, чи О(h2). Це безпосередньо випливає з розкладання в ряд Тейлора виразу

,

із якого отримуємо

Отже, крайові умови (11.23) апроксимовані з першим порядком по h. Порядок їх апроксимації можна підвищити до другого, наприклад, наступними кінцево-різницевими співвідношеннями:

(11.24)

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

,

з яких при множенні першого розкладання на 4 і вирахуванні його з другого рівняння отримуємо:

.

Розглянемо ще одну можливість врахування крайових умов типу (11.5) на прикладі умови:

.

Для цього вводиться додаткова точка xn + 1 = b + h за межами інтервалу [a,b], для якої з глобальною похибкою O(h2) справедливе співвідношення із симетричною апроксимацією:

, (11.25)

Цю точку можна виключити, скориставшись апроксимацією диференціального рівняння y(x) = f(x) в кінцевій точці інтервалу xn = x0 + nh:

,

звідки отримуємо рівняння

, (11.26)

яким можна замінити останнє рівняння в системі алгебраїчних рівнянь, одержаній для рівняння y(x) = f(x) кусочно-різницевою апроксимацією похідних. Аналогічно можна вчинити з першою умовою із (11.5) і першим апроксимуючим рівнянням для x1 . Варто підкреслити, що врахування крайових умов різних типів впливає тільки на перше й останнє рівняння цієї системи.

Приведемо подібні члени в рівнянні (11.21), отримаємо тоді стандартне трьохточкове різницеве рівняння:

. (11.27)

Додаючи до системи рівнянь (11.25) різницеві рівняння (11.23) чи (11.24), отримаємо в кожному випадку системи рівнянь, що містять n + 1 рівняння з n + 1 невідомих ui. Порівняємо ці два варіанти апроксимації крайової задачі. У першому з них СЛАР, утворена рівняннями (11.21) і (11.23), має трьохдіагональну матрицю коефіцієнтів, і до неї можна застосувати метод прогонки (глава 3). В другому випадку відповідна методу прогонки трьохдиагональна матриця ще повинна бути утворена. Для цього потрібно з першого рівняння (11.27) при i = 1

і першого рівняння (11.24) виключити u2. Виконавши це виключення, отримаємо:

Отримано рівняння, що містить дві невідомі u0 і u1. Замінимо ним перше рівняння (11.24). Проробимо аналогічні перетворення з другим рівнянням (11.24) і останнім рівнянням (11.25) при i = n - 1:

.

Виключимо з них un - 2, знаходимо:

Одержали рівняння, що містить дві невідомі un і un - 1. Замінимо ним друге рівняння (11.27). Додаючи два останніх рівняння до рівнянь (11.27), отримаємо систему рівнянь з трьохдіагональною матрицею, що апроксимує вихідну крайову задачу (11.4), (11.5) з порядком O(h2), яку також можна розв’язати методом прогонки. Вище в главі 3 було показано, що стійкість методу прогонки можна гарантувати, коли матриця коефіцієнтів має діагональну домінантність. Виберемо крок h таким, щоб забезпечити діагональну домінантність. Для цього досить зажадати виконання наступних умов для рівнянь (11.27):

і .

Підсилюючи останні нерівності, приходимо до наступних обмежень:

. 11.28)

Для крайових умов (11.23) повинні виконуватися умови:

. (11.29)

Наявність обмежень (11.28) і (11.29) свідчить про умовну стійкість розглянутого методу апроксимації.

Приклад 11.3.

Розглянемо розв’язок різницевим методом задачі

(11.30)

Відповідна різницева задача має вигляд:

Виберемо крок, задовольнивши умові (11.26), що забезпечує стійкість різницевої схеми:

.

Функція Ln(x) монотонна, тому вона досягає найбільшого по модулю значення на одному з кінців інтервалу. Визначимо на якому: |Ln(0,5)| = 0,693147, |Ln(1,5)| = 0,405465.

Таким чином, h0,693147. Умови (11.29) виконані. Виберемо крок рівним h = 0,1, щоб розв’язувати з придатною точністю. Формули для прогонки в розглянутому випадку мають вигляд:

Виконавши прямий і зворотний прогони, отримаємо чисельні розв’язки задачі, графік якого приведений на рис.11.3

Рис. 11.3. Графіки розв’язку різницевим методом крайової задачі (11.30)

На рис. 11.3. приведені графіки чисельного розв’язку задачі (11.30) різницевим методом і стандартним оператором пакета Mathematica -5.

Приклад 11.4.

Проілюструємо ефективність застосування екстраполяції Річардсона при розв’язку крайових задач на прикладі рівняння

y’’ + y = x

з крайовими умовами y(0) = , y(π/2) = π/2 - 1.

Диференціальне рівняння апроксимується лінійною системою рівнянь:

Ay = h2(х - у) - r, (11.31)

де

A = ; r = .

Результати розв’язку системи (11.31) при розбивці інтервалу інтегрування [1, π/2 - 1] на n = 5 і n = 10 відрізків систематизовані в табл. 11.1. Там же приведені покращення, отримані екстраполяцією Річардсона (8.23) з поправками ∆/(qp1 - 1) = /3, тому що вибрано p1 = 2 і q = 2. Оскільки відомий точний розв’язок крайової задачі y(x) = cosx - sinx + x, то в останньому стовпці таблиці приводиться оцінка реалізованої похибки ε.

Таблиця 11.1. Застосування екстраполяції Річардсона

m = 5

 y(x,0.1)

m = 10

y(x,0.05)

Річардсон

106

0

1,000000

1,000000

0

1,000000

0

0,1

0,988402

0,2

0,956572

0,956291

- 94

0,956197

- 2

0,3

0,908337

0,4

0,849741

0,849597

- 48

0,849549

- 1

0,5

0,785398

0,6

0,721056

0,721199

+ 48

0,721247

1

0,7

0,662460

0,8

0,614224

0,614505

+ 94

0,614599

1

0,9

0,582395

1,0

0,570796

0,570796

0

0,570796

0

11.4. Власні значення однорідної крайової задачі

Однорідна крайова задача

y’’ + y = 0 , y(0) = 0 , y(1) = 0 (11.32)

завжди має тривіальний розв’язок y(х) = 0. Однак для практичного застосування цікаві нетривіальні розв’язки, існування і вигляд яких визначається параметром λ. Ці виняткові значення параметра λ називають власними значеннями задачі , а відповідні їм нетривіальні розв’язки — власними функціями задачі. Наприклад, задачі (11.32) відповідає розв’язок

y(x) = acos(x ) + bsin(x ),

аналіз якого з урахуванням крайових умов приводить до наступного:

y(0) = 0  a = 0

y(1) = 0  = n, n = 0,  1,  2,…

Отже, власними значеннями задачі (11.32) будуть числа

= n22 (11.33)

а власними функціями

y(x) = bsin(nx). (11.34)

Приклад 11.5.

Знайдемо власні значення крайової задачі (11.32) чисельно. Для цього на інтервалі [0,1] виділимо три відрізки (n = 3) і з кроком h = 1/3 проведемо різницеву апроксимацію другої похідної

.

Тоді рівняння (11.32) заміщається системою рівнянь:

, (11.35)

в якій використані координати проміжних точок інтервалу (рис. 11.4).

Рис. 11.4. Розбиття відрізка до задачі 11.5

Якщо є нетривіальний розв’язок задачі (y1≠0, y2≠0), то визначник системи рівнянь (11.35) повинний дорівнювати нулю (11.35):

звідки h2 - 2 = 1.

При h = 1/3 отримуємо значення 1 = 9 (точне значення 2 = 9,8696) і 2 = 27 (точне значення 42 = 39,48).

Похибка визначення власних значень визначається виразом :

(h) = + c2h2 + c3h3 + C4h4 + …

Можна істотно підвищити точність визначення власних значень задачі, застосувавши екстраполяцію Річардсона (8.23). Для цього проведемо ще одне додаткове розбиття інтервалу і оцінимо визначник системи рівнянь, побудованої для h = 1/4. Отримані результати при оцінці значення λ1 зведені в таблицю 11.2 , з якої випливає, що отримана оцінка λ1 = 9,8517 досить добре наближається до точного розв’язку (λ1 = 9,8696) .

Таблиця 11.2. Уточнення власних значень задачі

H

λ

∆/(7/9)

Річардсон

1/3

9

0,4791

9,8517

- 0,0176

¼

9,3726

З табл.11.2 випливає, що, завдяки застосуванню екстраполяції Річардсона, отримано оцінку λ1 з малою похибкою при виборі невеликої кількості точок при розбивці інтервалу, що дуже важливо для спрощення обчислення необхідних визначників.

11.5. Метод колокацій

Розв’язок крайової задачі (11.4), (11.5) шукається у вигляді функції

(11.36)

де φi(x), i =  — лінійно незалежні, двічі диференційовані базисні функції, визначені на відрізку [a,b]. Функція φ0(x) повинна задовольняти заданим крайовим умовам (11.5):

(11.37)

а функції φi(x), i =  — відповідним однорідним крайовим умовам, тобто

(11.38)

В силу лінійності крайових умов функція u(x), обумовлена виразом (11.36) задовольняє крайовим умовам (11.24) при будь-яких значеннях коефіцієнтів ci. Наприклад, у точці x = a маємо

.

Аналогічно для x = b отримаємо

.

Ідея методу колокацій полягає в тому, що, вибравши n різних точок

на відрізку [a,b], названих вузлами колокацій, підбирають значення ci так, щоб отримана при цьому функція u(x) (11.36) задовольняла рівнянню (11.4) у кожному з вузлів інтерполяції. Тобто для кожного вузла колокацій xj повинне виконуватися рівність:

, (11.39)

де .

Покладемо

, (11.40)

тоді (11.39) здобуває стандартний вигляд системи лінійних алгебраїчних рівнянь:

(11.41)

відносно коефіцієнтів ci. Розв’язавши цю систему будь-яким методом і підставивши отримані значення коефіцієнтів у вираз (11.36), отримаємо наближений розв’язок u(x).

Точність розв’язку крайової задачі методом колокацій сильно залежить від вдалого вибору базисних функцій φi(x). У конкретних задачах вибір цих функцій повинний здійснюватися на основі змістовних представлень про їх розв’язки чи на основі емпіричних даних. При їх відсутності можна використовувати запропонований в [4] метод. В якості φ0 візьмемо лінійну функцію

, (11.42)

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


. (11.43)

Функції φi(x) можна взяти у вигляді:

. (11.44)

Очевидно, що при будь-яких Si функція (11.43) задовольняє першій умові (11.38). Значення Si, при якому буде задовольнятися друга умова (11.38), дорівнює

(11.45)

Приклад 11.6.

Розглянемо розв’язок задачі

L[y(x)] = , , (11.46)

,  (11.47)

Порівнюючи умови (11.47) із загальними крайовими умовами ,бачимо, що:

a = 1,b = 2,α0 = 0, β0 = 0, α1 = 3, β1 = 1, γ1 = 1, γ1 = 0,5

Перейдемо до вибору базисних функцій. Знайдемо коефіцієнти функції (11.42) ,

розв’язавши систему рівнянь (11.43) для умов (11.47), отримаємо:

φ0(x) = 0,625x + 1.625.

Сформуємо базисні функції φi(x) у вигляді (11.43), (11.44)

φ1(x) = 1,25(x - 1) + (x - 1)2,φ2(x) = - 1,2(x - 1)2 + (x - 1)3 

Визначимо диференціальний оператор L[u(x)], що відповідає лівій частинні рівняння (11.4) для наближеного розв’язку:

.

Виразимо нев'язку диференціального рівняння для наближеному розв’язку u(x):

= 6-3x2+(-2c1+8.4c2)x4+(1.625+2.25-8.2c2)x5+
+(0.+0.c1+0.c2)x6+(-c1+4.2c2)x7-2c2x8

Виберемо дві точки колокацій з координатами x1 = 1,33, x2 = 1,667.

Поставимо вимогу, щоб нев'язки в точках колокацій були рівні нулю й отримаємо рівняння для визначення c1, c2:

5.7334-4.32347c1+3.48198c2 = 0
13.0213-22.2528
c1-9.71251c2 = 0

Розв’яжемо отримані рівняння і сформуємо наближений розв’язок:

= -0.597(-2.47684+x)*
*(3.27523-3.14045x+x2)

Знайдемо точний розв’язок крайової задачі (11.46), (11.47) у пакеті Mathematica-5:

Dsolve[{-x5y[x]+-x6y’[x]-x4y”[x] = = f[x], y[a] = = 1, y’[b]+3y[b] = = 0.5}, y[x],x]
Out = {{y[x]→1/(x
2}}

Приведемо графіки наближеного і точного розв’язку, що дорівнює y = x - 2.

Рис. 11.5. Графіки точного розв’язку крайової задачі (11.46), (11.47), який дорівнює y = x - 2, і наближеного розв’язку u[x], отриманого методом колокацій

Використовуючи три базисних функції, отримано досить точний наближений розв’язок крайової задачі (11.46), (11.47). Максимальна помилка досягається при x = 2 і дорівнює
ε = U[2] - 0,25 = - 0.33

11.6. Метод Гальоркіна

Як і у випадку методу колокацій, будемо шукати наближений розв’язок крайової задачі (11.4) і (11.5) у вигляді (11.36):

, (11.36)

де φi(x), i = 0,1,2,…  лінійно незалежні, двічі диференційовані базисні функції, визначені на відрізку [a,b]. Функція φ0(x) повинна задовольняти заданим крайовим умовам (11.37), а функції φi(x), i = 0,1,2,… — відповідним однорідним крайовим умовам (11.38).

Додатково зажадаємо, щоб система базисних функцій φi(x), i = 0,1,2,… була

ортогональної на відрізку [a,b] тобто

при  і . (11.48)

і повною, тобто не існує ніякої іншої відмінної від нуля ортогональної функції до всіх функцій φi(x), i = 0,1,2,…

Підставимо наближений розв’язок (11.36) у рівняння (11.3) і знайдемо нев'язку:

. (11.49)

Виберемо коефіцієнти ci таким чином, щоб значення інтеграла від квадрата нев'язки

було найменшим.

Доведено [3, т.2] , що це досягається лише в тому випадку, коли нев'язка R(x,c1,c2,…,cn) ортогональна до всіх базисних функцій φi(x).

Запишемо умову ортогональності:

або, у більш докладній формі:

. (11.50)

Одержали систему лінійних алгебраїчних рівнянь відносно коефіцієнтів ci.

Приклад 11.6.

Розглянемо методом Гальоркіна розв’язок задачі:

, (11.51)

y( - π) = y(π) = 2. (11.52)

Виберемо як систему базисних функцій φi(x), i =  наступні тригонометричні функції:

(11.53)

Ці функції лінійно незалежні на відрізку [ - π,π], причому функція φ0(x) задовольняє крайовій умові (11.52), а інші функції — нульовим крайовим умовам. Будемо шукати розв’язок у вигляді:

Підставимо наближений розв’язок в рівняння (11.51) і знайдемо нев'язку (11.49):

Отримаємо систему лінійних алгебраїчних рівнянь (11.50) відносно коефіцієнтів ci

Звідки випливає с3 = 0, с1 = 0, с4 = - 2/7, с2 = - 15/14 і наближений розв’язок має вигляд

Графік розв’язку приведено на рис.11.6

Рис. 11.6. Графік наближеного розв’язку U[x] крайової задачі (11.51), (11.52) , отриманого методом Гальоркіна

11.7. Метод найменших квадратів

Будемо шукати наближений розв’язок крайової задачі (11.4) і (11.5), як і в попередньому параграфі, у вигляді (11.36):

, (11.54),

де φi(x), i = 0,1,2,… — лінійно незалежні, двічі диференційовані базисні функції, визначені на відрізку [a,b]. Функція φ0(x) повинна задовольняти заданим крайовим умовам (11.37), а функції φi(x), i = 0,1,2,… — відповідним однорідним крайовим умовам (11.38).

Підставимо наближений розв’язок (11.54) у рівняння (11.4) і знайдемо нев'язку:

, (11.55)

яка повинна бути при 0≤xb якомога меншою по абсолютній величині. Тому зажадаємо, щоб

виконувалася умова:

. (11.56)

Для мінімуму інтеграла I необхідне виконання умов:

.

У результаті формується система лінійних рівнянь відносно коефіцієнтів c1,c2,…,cn, з якої вони і визначаються.

Приклад 11.7.

Розглянемо розв’язок задачі:

, , (11.57)

y(1) = 1, 3y(2) + y’(2) + y’(2) = 0,5 (11.58)

розв’язаної вище в прикладі 11.5 методом колокацій, використовуючи ті ж самі базисні функції, і ту ж форму представлення розв’язку:

Підставимо u[x] в (11.36) і знайдемо нев'язку:

R(x) = 6-3x2+(-2c1+8.4c2)x4+(1.625+2.25-8.2c2)x5+(0.+0.c1+0.c2)x6+
+(-c1+4.2c2)x
7-2c2x8

Отримаємо рівняння для визначення коефіцієнтів c1, c2:

Розв’яжемо, отриману систему рівнянь і знайдемо:

с1 = 0,789449, с2 = - 0,369488.

Сформуємо наближений розв’язок:

u(x) = - 0,369488 ( - 2,71191 + x) (4,20565 - 3,6247x + x2)

Приведемо графіки розв’язку, отриманого методом найменших квадратів u[x] і точне значення x - 2.

Рис. 11.7. Графіки точного розв’язку крайової задачі (11.57), (11.58), рівного y = x - 2, і
наближеного розв’язку
w[x], отриманого методом найменших квадратів у Mathematica-5

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

ВИСНОВКИ:

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

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

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

4 Метод кінцевих різниць дозволяє також розв’язувати задачу про власні значення і власні функції крайової задачі, які визначають нетривіальні розв’язки однорідної крайової задачі.

5 Застосування методу кінцевих різниць до нелінійних крайових задач можливе, але вимагає крім алгебраїзації похідних також лінеаризації нелінійностей, що досягається при обчисленнях організацією “цикл в циклі”, що передбачає в кожній точці розбивки інтервалу виконання ітерацій Ньютона-Рафсона для уточнення значень використовуваних змінних.

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

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

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

Вправи:

1 Переконайтеся, що функція y[x] = x2-0,2525826491x-2,528442297xln(x) є розв’язком крайової задачі  з межовими умовами y(0,5) = 1 і y(3) = -1.

2 Розв’язати крайову задачу  межовими умовами y(0) = 1,25, y(4) = -0,95 на інтервалі 0<x<4 методом зведення до задач Коші. Порівняйте отриманий розв’язок з чисельним розв’язком, визначеним у Mathematica -5. Побудуйте графіки розв’язків.

3 Розв’яжіть крайову задачу  з межовими умовами y(0) = 1, y(1) = 0,1 на інтервалі 0<x<1 різницевим методом. За допомогою Mathematica-5 знайти точний розв’язок і визначити похибку наближеного розв’язку.

4 Розв’язати крайову задачу  з межовими умовами y(0) = 1, y(1) = -0,5 на інтервалі 0<x<1 різницевим методом. За допомогою Mathematica-5 знайти точний розв’язок і визначити похибку наближеного розв’язку.

5 Розв’язати крайову задачу  з межовими умовами y(0) = 0, y(1) = 0 на інтервалі 0<x<1 методом колокацій. За допомогою Mathematica-5 знайти точний розв’язок і порівняти обидва розв’язки.

6 Розв’язати крайову задачу  з межовими умовами y(0) = 1, y(1) = 0 на інтервалі 0<x<1 методом Гальоркіна. За допомогою Mathematica-5 знайти наближений розв’язок і порівняти обидва розв’язки.

7 Розв’язати крайову задачу з межовими умовами y(1) = 0, y(2) = 1на інтервалі 1<x<2 методом найменших квадратів. За допомогою Mathematica-5 знайти наближений розв’язок і порівняти обидва розв’язки.

8 Розв’язати крайову задачу  з межовими умовами y(1) = 0, y(2) = 1 на інтервалі 1<x<2 методом колокацій і Гальоркіна. Порівняти обидва розв’язки.

9 Розв’язати крайову задачу  з межовими умовами y(1) = 0, y(2) = 1 на інтервалі 1<x<2 методом найменших квадратів і колокацій. Порівняти обидва розв’язки.

10  Знайти мінімальне власне значення min задачі , y(-1) = y(1) = 0, застосовуючи метод кінцевих різниць при h = 2/3 і h = 2/5 і апроксимацію типу  для лівої частини рівняння .
Відповідь:
 3,64

PAGE  217

EMBED Equation.3  


 

А также другие работы, которые могут Вас заинтересовать

66615. Разработка программного обеспечения автоматизированной системы проектирования операций сверления отверстий 198 KB
  Ввод исходных данных о материале детали и режущей части инструмента; геометрии, погрешности заточки и размерах инструмента; об условиях обработки (осевая подача, точность оборудования, размеры и погрешность заготовки и др.).
66616. Характеристика выпрямителя и сглаживающего фильтра блока питания ПК 5.21 MB
  Классификация источников питания СВТ В зависимости от характера преобразования энергии в источнике питания выполняемого при получении на его выходе требуемого напряжения источники питания подразделяются на: Первичные источники питания; Вторичные источники питания.