Еліптичні рівняння

До еліптичних відносять рівняння виду:

     рівняння Пуассона          (11.8)

      рівняння Лапласа          (11.9)

Розглянемо класичну задачу Діріхле

    

у деякій області R і на границі цієї області, якою є крива G. Звичайно, рівняння (11.9) являє собою рівняння Лапласа, окремий випадок рівняння Пуассона.

Припустимо що область R має вигляд прямокутника ширини А і висоти В

Графічне представлення переходу до кінцево-різницевих формул

Розділимо спочатку ширину прямокутника А на n інтервалів, кожен розміром ; а висоту В – на m частин розмірами . Усередині області виходять при цьому перерізів сітки. Визначимо різницеве співвідношення для кожної внутрішньої точки і розв‘яжемо отриману систему рівнянь.

Відношення між похідними та значеннями функції у вузлах сітки

Похідна	Схема	Наближена формула	Порядок похибки
			h2 h h
			e2 e e
			h2 h4 h4
			e2 e4 e4

Нумерація точок перетинання сітки починається в горизонтальному напрямку зліво-направо, від нуля, крайня права точка буде при цьому n-й. Аналогічно у вертикальному напрямку нумеруємо точки знизу-вверх від нуля до m. Вузол з індексами (i, j) буде i-м ліворуч і j-м знизу.

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

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

          (11.11)

Якщо позначити , то диференціальне рівняння (11.9) зведеться до різницевого рівняння вигляду           (11.12)

для i=1, 2, 3, …, n-1 та j=1, 2, 3, …, m-1...

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

          (11.13)

Важливо те, що одержана система – це система лінійних алгебраїчних рівнянь. Усього є рівнянь щодо невідомих. Після того як 2(m+n) невідомих будуть виключені за допомогою граничної умови (11.11), залишається точно рівнянь щодо невідомих. Вираз (11.7) являє собою систему лінійних алгебраїчних (n-1)(m-1) рівнянь відносно (n+1)(m+1) невідомих. Після того, як 2(n+m) невідомі будуть виключені за допомогою граничних умов, залишається (n-1)(m-1) рівняння відносно (n-1)(m-1) невідомих.

У (11.11) і (11.12) кожна пара значень (і, j) визначає вузол, в якому рівняння (11.6) розв'язується відносно . Рівняння (11.12) можна легко запам'ятати, якщо знати різницеву схему п'яти вузлів (рис. 11.2, 11.4.д,е).

Розпишемо докладно деякі з рівнянь (11.12). Для зручності прийнявши , тобто , але спільність висновків не важко перевірити при кожному . Починаючи з і при незмінному пройдемо значення :

          (11.14)

Тепер збільшуючи до 2, пройдемо значення при :

          (11.15)

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

Система (11.14), (11.15) має дві важливі властивості, що допомагають вибрати для неї метод розв’язку. Властивості ці наступні:

1.     Переважаюча частина коефіцієнтів системи дорівнює нулю.

2.     У кожному рівнянні один з коефіцієнтів дорівнює + 4. Якщо в рівнянні є п'ять коефіцієнтів, відмінних від нуля, то сума інших чотирьох коефіцієнтів дорівнює – 4, якщо ж кількість ненульових коефіцієнтів менше п'яти, то сума інших дорівнює – 2 чи – 3.

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

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

          (11.16)

Для цієї системи рівнянь виконується умова зближення ітераційного методу Гауса - Зейделя. Щоб розв'язати цю систему, розроблено алгоритм (рис. 11.5), який використовує метод Гауса–Зейделя, який застосовують до еліптичних рівнянь різницевих, називають методом Лібмана (або методом послідовних зсувів).

    При порівнянні реальних систем ДР буває корисно збільшити чи зменшити чергове виправлення при розрахунку чергового наближення до кореня рівняння (аналогічний метод може виявитися корисним при розв’язанні систем лінійних алгебраїчних рівнянь). У дійсності це справедливо і для методу Лібмана. Щоб можна було збільшувати чи зменшувати чергові виправлення, замінимо в алгоритмі (рис. 11.5) визначення за формулою:

    

де х – значення u_i,j, обчислене з рівняння (11.12). Параметр ω називається прискорюючим множником. У загальному випадку вимагається, щоб значення ? лежало в межах від 1 до 2. Якщо ω = 1, то є звичайний метод Лібмана, при ω > 1 є прискорений метод Лібмана..

Метод Лібмана було розглянуто тільки для випадку рівняння Лапласа. Взагалі кажучи, будь-яке еліптичне рівняння без членів, що містять и_ху, зводиться до системи різницевих рівнянь, що задовільняє умовам збіжності.

Усе сказане дотепер відносилося до лінійних рівнянь. Питання про збіжність для нелінійних рівнянь розроблено дуже слабо. Потрібно відзначити, однак, що багаторазово робилися успішні спроби вирішення квазілінійних (квазілінійними називаються рівняння, лінійні по и_хх, u_yy і u_xy, але нелінійні по u_x, и_у, і, х, у) рівнянь за допомогою екстрапольованого методу Лібмана. Загальний підхід у цьому випадку зводиться до того, що задається навмання деяке значення і спостерігається швидкість збіжності інтераційного процесу. Потім, змінюючи ω і спостерігаючи викликані цим зміни швидкості збіжності, можна підібрати таке її значення, при якому ітераційний процес сходиться швидше всього. Вдавалося одержати збіжність навіть для рівнянь із членом типу и_ху, хоча в цьому випадку немає ніяких теоретичних основ очікування збіжності.

Є багато інших способів вирішення різницевих рівнянь (11.12). Найбільш часто використовуються лінійна релаксація, блокова релаксація і метод напрямків, що чергуються. Часто вони виявляються більш ефективними, ніж метод Лібмана. Нарешті, можна використовувати екстраполяційний перехід до границі.

Різницеві схеми для розв’язання рівняння теплопровідності, коливання струни.

Як приклад фізичної задачі, що приводить до рівняння цього типу, розглянемо процес теплопередачі по довгому стрижні, що лежить уздовж осі х від х = 0 до х = L. Припускаючи, що в точці х = 0 температура підтримується на рівні Те, а в тоці х = L температура підтримується на рівні Т_L. Припускаючи також, що в момент часу t = 0 розподіл температури уздовж стрижня задана функцією f(х). Тоді розподіл температури уздовж стрижня в усі наступні моменти часу отримується шляхом розв‘язку рівняння:      (12.1) де і – температура стрижня в даній точці в даний момент часу, а постійна , де с – теплоємність матеріалу стрижня, – щільність матеріалу стрижня і k – його теплопровідність. Для простоти припустимо а = 1, тоді рівняння зведеться до наступного вигляду:      , або .     (12.2)

Граничними умовами для цього рівняння є:           (12.3) і початковою умовою є:      .     (12.4)

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

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

При записі в різницевій формі граничні умови запишуться так:           (12.5)

Початкова умова (12.4) запишеться у вигляді      .     (12.6)

Щоб перетворити в різницеву форму (12.2), знову уявимо собі сітку, що охоплює область 0 < х < L і t > 0 з інтервалом розбиття k у напрямку х та інтервалом розбиття k у напрямку t. Використовуючи різницеві формули одержуємо:      ,     (12.7)      де     ,     (12.8)

та індекс i змінюється від 1 до n-1, а індекс j – від 1 до . Верхня межа помилки обмеження виражається співвідношенням      ,     (12.9)

де , за умови 1/6, тому що тоді це вираження перетворюється в нуль. При =1/6 верхня межа помилки обмеження виглядає так:      ,     (12.10) де .

Р ізницева схема представлена на рис. 12.2.

    

Різницева схема для параболічного рівняння (випадок 1)

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

          (12.11)

Це не завжди виявляється зручним.

Метод, що привів до рівняння (12.7), був заснований на використанні правих різниць (11.3, 11.5). Якщо використовувати ліві різниці (11.4, 11.6), то різницеве рівняння запишеться в наступному вигляді:           (12.12) де означає те ж, що й у (12.8).

Р ізницева схема цього методу виглядає так:

    

Різницева схема для параболічного рівняння (випадок 2)

Верхня межа помилки обмеження визначається формулою:      (12.13) де .

Розглянемо ще один спосіб побудови різницевого рівняння, еквівалентного (12.12). Згадаємо, що, згідно (11.5), у точці (i, j)

Аналогічно в точці (i, j+1)

Усереднюючи ці два наближення, одержимо

Якщо тепер для обчислення u_i скористатися правими різницями, то    

і різницеве рівняння, еквівалентне (12.12), запишеться у вигляді

          (12.14)

Цей спосіб побудови різницевого рівняння часто називають методом Кранка – Нікольсона.

Верхня межа помилки обмеження дається формулою           (12.15) де

Ц я верхня межа в 6 + 1 разів менше, ніж (12.13).

    

    Різницева схема для параболічного рівняння (випадок 3)

Оскільки межа помилки обмеження (12.15) менше, то по цій характеристиці різницеве рівняння (12.14) переважніше, ніж (12.8).

Вище розглянуто три різних способи різницевого представлення диференціального рівняння и_хх = и_t а саме формули (12.7), (12.12) і (12.14). Розглянемо тепер методи розв’язання кожного з цих різницевих рівнянь і зіставимо їхні переваги і недоліки. В усіх випадках вважається, що граничні і початкові умови задані у вигляді (12.5) і (12.6).

Перший спосіб, обумовлений формулою (12.7), являє собою явну систему рівнянь для u_i,j+1 у тому же самому сенсі, у якому (11.27) було явною системою для випадку гіперболічного рівняння. Маючи на підставі початкових і граничних умов перший ряд розв’язку, можна обчислити другий рядок, j = 1, безпосередньо з (12.3), при j = 0. Обчисливши другий рядок, можна в такий же спосіб обчислити третій рядок, поклавши цього разу j = 1. Таким же способом можна продовжувати обчислювати розв’язок настільки далеко по осі часу, яка у цьому є необхідність.

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

чи, те ж саме, при

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

І (12.12), і (12.14) являють собою неявні методи розв’язання, аналогічні тому, що був застосований для еліптичних рівнянь. Розглянемо спочатку (12.12). Можна виписати рівняння для першого рядка розв’язку j=1, використовуючи початкові і граничні умови (12.3) і (12.4):

Усього виходить n - 1 лінійних рівнянь з п - 1 невідомими для i = 1,2,...,n-1. Ці рівняння можна вирішити методом виключення; крім того, що переважна більшість коефіцієнтів дорівнює нулю, сам вид системи робить рішення методом виключення досить просто. Справа в тім, що єдині ненульові коефіцієнти розташовані в трьох діагональних рядах симетрично щодо головної діагоналі. Система рівнянь цього типу називається тридіагональною. Можна вирішувати цю систему і за допомогою ітераційного методу Гаусса – Зейделя, тому що діагональний коефіцієнт у кожному рівнянні дорівнює 1 + 2 , а сума недіагональних дорівнює –2 , що забезпечує збіжність. Для прискорення збіжності можна скористатися екстраполяцією.

Знайшовши розв’язання для першого рядка, тобто для t = k, j = 1, можна записати систему рівнянь для другого рядка і знову одержати тридіагональну систему з n – 1 рівнянь. Після вирішення цієї системи можна переходити до третього рядка і так далі. Для обчислення одного рядка рішення щораз приходиться вирішувати систему з п – 1 лінійних алгебраїчних рівнянь.

Оскільки трудомісткість цього способу набагато більше, ніж для явних формул (12.7), можна задати собі питання, чи має сенс узагалі, прибігати до такого способу? Відповідь на це питання полягає в тім, що (12.12) стійке і сходиться для всіх > 0. Тому, використовуючи метод (12.12), можна брати набагато більшу величину кроку за часом. Спробуємо визначити, чи можна дійсно домогтися економії машинного часу.

Відомо, що при використанні методу (12.12) витрати машинного часу на обчислення одного рядка рішення будуть удвічі більше, ніж при використанні методу (12.7). Якщо взяти = 1, то загальна витрата часу на розрахунок рішення до визначеного моменту буде вдвічі менше. Правда, потрібно врахувати, що помилка обмеження круто зростає з ростом кроку за часом, так що визначальним фактором може з'явитися саме максимально припустима помилка.

Остання різницева формула (12.14) є іншим неявним методом і для обчислення одного рядка також вимагає рішення тридіагональної системи з n - 1 рівнянь. Витрати машинного часу при використанні цієї формули будуть, отже, такими ж, як і для (12.12). Цей метод також стійкий для всіх .

Помітимо, що при подвоєнні у формулі (12.14) помилка обмеження подвоюється, у той час як у (12.12) помилка обмеження збільшувалась в 4 рази. Тому метод (12.14) (метод Кранка - Нікольсона) кращий, ніж (12.12).

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

На рис. 12.5 представлено чотири різницеві схеми які використовуються для рівнянь другого порядку в частинних похідних.

    

Різницеві схеми, для еліптичного (а), гіперболічного (б) і параболічного (в, г)

Приклад розв’язку відомої інженерної задачі моделювання розподілу температури в трубі квадратного перетину

Розглянемо довгу квадратну трубу з квадратним отвором, по якій тече гаряча рідина. Труба наполовину занурена в крижану ванну, так що температура нижньої половини поверхні труби дорівнює 0°С. Верхня площина труби знаходиться при постійній температурі 100°С. Крім того, ми припускаємо, що температура зовнішньої поверхні труби лінійно змінюється від 0° С до 100° С на ділянці між крижаною ванною і верхньою площиною труби. Рідина усередині труби має температуру 200°С. Зовнішній розмір труби дорівнює 10 дюймам, внутрішній розмір дорівнює 4 дюймам. Поперечний переріз труби зображений на рис. 12.6.

    

     Рисунок 12.6– Перетин труби в практичному прикладі

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

,

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

закінчитися, тобто тепловий режим труби став стаціонарним. Тоді і рівняння здобуває наступний вид:

Цьому рівнянню задовольняє розподіл температур усередині труби, а розподіл температур на границях труби задано граничними умовами. Таким чином, наша задача являє собою задачу Дирихле для рівняння Лапласа. Природно, мається на увазі двовимірна задача, тобто, наприклад, не враховується зниження температури рідини уздовж осі труби. Інакше кажучи, ми будемо обчислювати розподіл температури тільки для одного якогось перетину труби.

Припустимо, що по зовнішньому розмірі в напрямках х и у зроблена розбивка на 60 інтервалів, так що h = k = 0,16667 дюйма. Тоді граничні умови запишуться так:

Різницеві рівняння виглядають так:

при та

та

Ми можемо вирішити ці рівняння прискореним методом Лiбмана і знайти температуру у всіх внутрішніх крапках труби.

Формула прискореного методу Лiбмана переписана для цієї задачі в наступному виді:

Контрольні питання:

1. Різницеві схеми для розв'язування параболічного рівняння на ЕОМ.

2. Особливість алгоритму розв'язання параболічного рівняння на ЕОМ.

3. Що таке збижність методу?

4. Від чого залежить збижність методів розв'язання еліптичного, гіперболічного та параболічного рівнянь?