MACROBUTTON MTEditEquationSection2 Equation Chapter 1 Section 1 SEQ MTEqn \ r \ * MERGEFORMAT SEQ MTSec \ r 1 \ * MERGEFORMAT SEQ MTChap \ r 1 \ * MERGEFORMAT Реферат
У курсовій роботі розглядається метод сіток рішення параболічних рівнянь. Теоретична частина включає опис загальних принципів методу, його застосування до розв'язання параболічних рівнянь, дослідження розв'язності одержуваної системи різницевих рівнянь. У практичній частині розробляється програма для чисельного вирішення поставленого завдання. У додатку представлений текст програми і результати виконання тестових розрахунків.
Обсяг курсової роботи: 33 с.Ілюстрацій: 5.
Графіків: 1.
Джерел: 4.
Ключові слова: параболічне рівняння, рівняння теплопровідності, метод сіток, крайова задача, кінцеві різниці.
Зміст
Введення
1. Теоретична частина
1.1 Метод сіток рішення рівнянь параболічного типу
1.2 Метод прогонки рішення різницевої задачі для рівнянь параболічного типу
1.3 Оцінка похибки і збіжність методу сіток
1.4 Доказ стійкості різницевої схеми
2. Реалізація методу
2.1 Розробка програмного модуля
2.2 Опис логіки програмного модуля
2.3 Приклад роботи програми
Висновок
Список джерел
Додаток
Введення
До диференціальним рівнянням з приватними похідними приходимо при вирішенні найрізноманітніших завдань. Наприклад, за допомогою диференціальних рівнянь з частинними похідними можна вирішувати задачі теплопровідності, дифузії, багатьох фізичних і хімічних процесів.
Як правило, знайти точне рішення цих рівнянь не вдається, тому найбільш широке застосування отримали наближені методи їх вирішення. У даній роботі обмежимося розглядом диференціальних рівнянь з частинними похідними другого порядку, а точніше диференціальними рівняннями з приватними похідними другого порядку параболічного типу, коли ці рівняння є лінійними, а шукана функція залежить від двох змінних. У загальному випадку таке рівняння записується таким чином:
Зауважимо, що чисельними методами доводиться вирішувати і нелінійні рівняння, але знаходити їх рішення багато важче, ніж рішення лінійних рівнянь.
введемо в розгляд величину
1. Теоретична частина
1.1 Метод сіток рішення рівнянь параболічного типу
Для вирішення диференціальних рівнянь параболічного типу існує декілька методів їх чисельного розв'язання на ЕОМ, однак особливе положення займає метод сіток, так як він забезпечує найкращі співвідношення швидкості, точності отриманого рішення і простоти реалізації обчислювального алгоритму. Метод сіток ще називають методом кінцевих різниць. Нехай дано диференціальне рівняння
Потрібно знайти функцію
Сітка може складатися з клітин різної конфігурації: квадратних, прямокутних, трикутних і інших. Після побудови сітки вихідне диференціальне рівняння замінюється різницевим рівнянням у всіх внутрішніх вузлах сітки. Потім на підставі граничних умов встановлюються значення шуканого рішення в граничних вузлах. Приєднуючи граничні умови сіткової завдання до різницевих рівнянь, записаних для внутрішніх вузлів, отримуємо систему рівнянь, з якої визначаємо значення шуканого рішення у всіх вузлах сітки.
Заміна диференціального рівняння різницевим може бути здійснена різними способами. Один із способів апроксимації полягає в тому, що похідні, що входять в диференціальне рівняння, замінюються лінійними комбінаціями значень функції
Розглянемо неоднорідне рівняння теплопровідності, що є окремим випадком рівнянь параболічного типу:
Будемо шукати рішення цього рівняння в області
Зауважимо, що цю півсмуги завжди можна привести до півсмузі, коли
де
Для вирішення завдання область
Вузли сітки, що лежать на прямих
Для довільної
Чи можемо отримати три види різницевих рівнянь:
Різницеві рівняння (1.5) апроксимують рівняння (1.2) з похибкою
У різницевої схемою (1.5) задіяні 4 вузли. Конфігурація схеми (1.5) має вигляд:
У схемі (1.6) також беруть участь 4 вузли, і ця схема має вигляд:
У схемі (1.7) беруть участь 5 вузлів, і ця схема має вигляд:
Перша і третя схеми - явні, друга схема неявна. У випадку явних схем значення функції у вузлі чергового шару можна знайти, знаючи значення у вузлах попередніх шарів. У разі неявних схем для знаходження значень розв'язку у вузлах чергового шару доводиться розв'язувати систему рівнянь.
Для вузлів початкового (нульового) шару
Для граничних вузлів, що лежать на прямих
Рівняння (1.9) апроксимують граничні умови (1.4) з похибкою
Приєднуючи до системи різницевих рівнянь, записаних для внутрішніх вузлів, початкові та граничні умови (1.8) і (1.9) для різницевої задачі отримаємо повні різницеві схеми трьох видів. Для проведення обчислень найпростішою схемою надається перша: достатньо на підставі початкової умови знайти значення функції у вузлах шару
Третя схема також вельми проста для проведення обчислень, але при її використанні необхідно окрім значень рішення у вузлах шару
З точки зору точкової апроксимації третя схема найточніша.
Введемо в розгляд параметр
У будь-якому випадку згідно з методом сіток будемо мати стільки рівнянь, скільки є невідомих (значення шуканої функції у вузлах). Число невідомих дорівнює числу всіх вузлів сітки. Вирішуючи систему рівнянь, отримуємо рішення поставленої задачі.
Розв'язність цієї системи для явних схем питань не викликає, тому що всі дії виконуються в явно певній послідовності. У разі неявних схем розв'язність системи слід досліджувати в кожному конкретному випадку. Важливим питанням є питання про те, на скільки знайдені рішення добре (адекватно) відображають точні рішення, і чи можна необмежено згущуючи сітку (зменшуючи крок по осях) отримати наближені рішення, як завгодно близькі до точних рішень? Це питання про збіжність методу сіток.
На практиці слід застосовувати сходяться різницеві схеми, причому тільки ті з них, які є стійкими, тобто при використанні яких невеликі помилки в початкових або проміжних результатах не призводять до великих відхилень від точного рішення. Завжди слід використовувати стійкі різницеві схеми, проводячи відповідні дослідження на стійкість.
Перша з побудованих вище різницевих схем у разі першої крайової задачі буде стійкою при
Явні схеми прості для організації обчислювального процесу, але мають один дуже вагомий недолік: для їх стійкості доводиться накладати сильні обмеження на сітку. Неявні схеми вільні від цього недоліку, але є інша проблема - треба розв'язувати системи рівнянь великої розмірності, що на практиці при знаходженні вирішення складних рівнянь в дальній області з високим ступенем точності може потребувати великих обсягів пам'яті ЕОМ і часу на очікування кінцевого результату. На щастя, прогрес не стоїть на місці і вже зараз потужності сучасних ЕОМ цілком достатньо для вирішення поставлених перед ними завдань.
1.2 Метод прогонки рішення різницевої задачі для рівнянь параболічного типу
Розглянемо окремий випадок задачі, поставленої в попередньому розділі. В області
знайти рішення рівняння
з граничними умовами
і початковою умовою
Розглянемо стійку обчислювальну схему, для якої величина
Запишемо різницеве рівняння, що апроксимує диференціальне рівняння (1.10) у всіх внутрішніх вузлах шару
Ці формули має похибку
Перепишемо (1.13) у вигляді:
Дана обчислювальна схема має наступну конфігурацію:
Система (1.14) - (1.16) являє собою різницеву завдання, відповідну крайовій задачі (1.10) - (1.12).
За величину
(1.14) - (1.16) є система лінійних алгебраїчних рівнянь з 3-діагональною матрицею, тому її резонно вирішувати методом прогонки, так як він у декілька разів перевершує за швидкістю метод Гаусса.
Тут
Підставивши рівняння (1.18) в (1.14) отримаємо:
Порівнявши (1.17) і (1.19) знайдемо, що:
Покладемо в (1.14)
Зауважимо, що в другій формулі (1.21) величина
За допомогою формул (1.21) і (1.20) проводимо прогонку в прямому напрямку. У результаті знаходимо величини
Потім здійснюємо зворотний хід. При цьому скористаємося другий з формул (1.15) і формулою (1.17). Отримаємо наступний ланцюжок формул:
Таким чином, відправляючись від початкового шару
Отже, ми побудували неявну схему вирішення диференціальних рівнянь параболічного типу методом сіток.
1.3 Оцінка похибки і збіжність методу сіток
При вирішенні задачі методом сіток ми допускаємо похибка, що складається з похибки методу та обчислювальної похибки.
Похибка методу - це та похибка, яка виникає в результаті заміни диференціального рівняння різницевим, а також похибка, що виникає за рахунок зносу граничних умов з
Обчислювальна похибка - це похибка, що виникає при рішенні системи різницевих рівнянь, за рахунок практично неминучих машинних заокруглень.
Існують спеціальні оцінки похибки для вирішення задач методом сіток. Однак ці оцінки містять максимуми модулів похідних шуканого рішення, тому користуватися ними вкрай незручно, однак ці теоретичні оцінки хороші тим, що з них видно: якщо необмежено подрібнювати сітку, то послідовність рішень буде сходитися рівномірно до точного розв'язання. Тут ми зіткнулися з проблемою збіжності методу сіток. При використанні методу сіток ми повинні бути впевнені, що, необмежено згущуючи сітку, можемо отримати рішення, як завгодно близьке до точного.
Отже, на прикладі розв'язку крайової задачі для диференціального рівняння параболічного типу розглянемо основні принципи методу сіток. Відзначимо, що якщо при вирішенні різницевої задачі невеликі помилки в початкових і крайових умовах (або в проміжних результатах) не можуть привести до великих відхилень шуканого рішення, то говорять, що задача поставлена коректно в сенсі стійкості за вхідними даними. Різницеву схему називають стійкою, якщо обчислювальна похибка необмежено не зростає. В іншому випадку схема називається нестійкою.
1.4 Доказ стійкості різницевої схеми
Нехай
і граничним умовам
Тут
Розглянемо похибка
Похибка
(В силу лінійності рівняння (1.14)), а також наступними граничними та початковими умовами:
Приватне рішення рівняння (1.23) будемо шукати у вигляді
Тут числа
При цілому
Підставимо рівняння (1.26) в рівняння (1.24). При цьому отримаємо:
або
Вираз у квадратних дужках одно
Підставляючи цей вираз у попереднє рівняння замість виразу у квадратних дужках і проводячи скорочення на
звідки знаходимо
Таким чином, згідно з рівнянням (1.26), отримуємо лінійно-незалежні рішення рівняння (1.23) у вигляді
Зауважимо, що це приватне рішення задовольняє однорідним крайовим умовам (1.24). Лінійна комбінація цих приватних рішень також є рішенням рівняння (1.23):
причому
У результаті отримуємо систему рівнянь
містить
Для стійкості досліджуваної різницевої схеми необхідно, щоб за будь-яких значеннях коефіцієнтів
Аналізуючи (1.28) бачимо, що ця нерівність виконується для будь-яких значень параметра
2. Реалізація методу
2.1 Розробка програмного модуля
Поставлено мету: розробити програмний продукт для знаходження наближеного рішення параболічного рівняння:
в області
задовольняє умовам
Розіб'ємо область
де
Замінивши у кожному вузлі похідні кінцево-різницевими відносинами по неявній схемі, отримаємо систему види:
Перетворивши її, отримаємо:
де
У граничних вузлах
У початковий момент
Ця різницева схема стійка при будь-якому
де
Аналогічно
Підставивши значення (1.35) в (1.32) отримаємо:
Звідки
З порівняння (1.35) і (1.37) видно, що
Для
Звідки
або
Звідки, використовуючи (1.35), отримаємо:
Використовуючи даний метод, ми всі обчислення проведемо в наступному порядку для всіх
1) Знаючи значення функції
2) Знайдемо
3) Знайдемо
4) Знайдемо значення шуканої функції на
2.2 Опис логіки програмного модуля
Лістинг програми наведено в додатку 1. Нижче будуть описані функції програмного модуля та його призначення.
Функція main () є базовою. Вона реалізує алгоритм методу сіток, описаного в попередніх розділах роботи.
Функція f (x, y) представляє собою вільну функцію двох змінних диференціального рівняння (1.29). Як аргумент на неї передаються два дійсних числа з плаваючою точкою типу float. На виході функція повертає значення функції
Функції mu_1 (t) і mu_2 (t) представляють собою крайові умови. У них передається по одному аргументу (t) дійсного типу (float).
Функція phi () є відповідальною за початковий умови.
У функції main () визначено такі константи:
Варіюючи
Програма забезпечена трьома механізмами виведення результатів роботи: на екран у вигляді таблиці, в текстовий файл, а також у файл списку математичного пакета Waterloo Maple. Це дозволяє наочно уявити отримане рішення.
Програма написана на мові програмування високого рівня Borland C + + 3.1 у вигляді додатку MS-DOS. Забезпечується повна сумісність програми з усіма широко відомими операційними системами корпорації Майкрософт: MS-DOS 5.x, 6.xx, 7.xx, 8.xx, Windows 9x/Me/2000/NT/XP.
2.3 Приклад роботи програми
В якості прикладу розглянемо чисельне рішення наступного диференціального рівняння параболічного типу:
в області
задовольняє умовам
Задавши прямокутну сітку з кроком осі
2.10 1.91 1.76 1.63 1.53 1.44 1.37 1.31 1.26 1.22 1.18
2.11 1.75 1.23 1.20 1.15 1.10 1.07 1.04 1.04 1.07 1.21
2.12 1.61 0.95 0.96 0.93 0.91 0.90 0.90 0.94 1.03 1.24
2.13 1.51 0.79 0.81 0.81 0.80 0.81 0.83 0.89 1.03 1.27
2.14 1.45 0.69 0.73 0.74 0.74 0.76 0.80 0.88 1.04 1.31
2.15 1.41 0.64 0.69 0.70 0.71 0.74 0.79 0.89 1.05 1.34
У таблиці вісь x розташована горизонтально, а вісь t розташована вертикально і спрямована вниз.
На виконання програми на середньостатистичному персональному комп'ютері витрачається час, що дорівнює декількох мілісекунд, що говорить про високу швидкість алгоритму.
Докладно вихідний файл output.txt, що містить таблицю значень функції
Висновок
У роботі було розглянуто метод сіток рішення параболічних рівнянь в приватних похідних. Розкрито основні поняття методу, апроксимація рівняння та граничних умов, досліджено розв'язність і збіжність одержуваної системи різницевих рівнянь.
На підставі вивченого теоретичного матеріалу була розроблена програмна реалізація методу сіток, проаналізовано її збіжність і швидкодію, проведений тестовий розрахунок, побудований графіки отриманого чисельного рішення.
Список джерел
1. Березін І.С., Жидков Н.П. Методи обчислень. Т.2. - М.: Физматгиз, 1962.
2. Тихонов О.М., Самарський О.А. Рівняння математичної фізики. - М.: Наука, 1972.
3. Пірумов Н.Р. Чисельні методи. - М.: Видавництво МАІ, 1998.
4. Каліткін М.М. Чисельні методи. - М.: Наука, 1976.
Додаток
Текст програми
/ / - / /
# Include <stdio.h>
# Include <conio.h>
# Include <math.h>
void main (void);
float f (float x, float t);
float mu_1 (float t);
float mu_2 (float t);
float phi (float x);
/ / - / /
void main (void)
{
clrscr ();
FILE * myfile;
FILE * plotter;
float a [120] [120];
float b [120] [120];
float u [120] [120];
float T = 0.05;
float l = 1;
float h = 0.1;
float tau = 0.01;
int n, i, j, k;
float s = pow (h, 2) / tau;
n = ceil (l / h);
for (i = 0; i <= 119; i + +)
{
for (j = 1; j <= 119; j + +)
{
u [i] [j] = 0;
a [i] [j] = 0;
b [i] [j] = 0;
}
}
for (i = 0; i <= n; i + +)
{
u [i] [0] = phi (i * h);
}
for (j = 0; j <= floor (T / tau); j + +)
{
u [0] [j] = mu_1 (tau * j);
u [n] [j] = mu_2 (tau * j);
}
for (j = 0; j <= floor (T / tau); j + +)
{
a [1] [j + 1] = 1 / (2 + s);
for (i = 2; i <= n - 1; i + +)
{
a [i] [j + 1] = 1 / (2 + s - a [i - 1] [j + 1]);
}
b [1] [j + 1] = mu_1 ((j + 1) * tau) + s * u [1] [j] + pow (h, 2) * f (h, (j + 1) * tau) ;
for (i = 2; i <= n - 1; i + +)
{
b [i] [j + 1] = a [i - 1] [j + 1] + s * u [i] [j] + pow (h, 2) * f (i * h, (j + 1) * tau);
}
u [n] [j + 1] = mu_2 ((j + 1) * tau);
for (k = 1; k <= n - 1; k + +)
{
u [n - k] [j + 1] = a [n - k] [j + 1] * (b [n - k] [j + 1] + u [n - k + 1] [j + 1] );
}
}
myfile = fopen («output.txt», «w +»);
plotter = fopen («3dplot.txt», «w +»);
fprintf (myfile, «Таблиця значень функції u = u (x, t) в області D = {0 <= X <=% g, 0 <= T <=% g}: \ n», l, T);
printf («Значення функції u (x, t) в області D = {0 <= X <=% g, 0 <= T <=% g}: \ n \ n», l, T);
for (j = 0; j <= floor (T / tau); j + +)
{
for (i = 0; i <= n; i + +)
{
printf ("% .2 f», u [i] [j]);
fprintf (myfile, «u (% g) (% g) =% g; \ n», i * h, j * tau, u [i] [j]);
if (i <n & & j <floor (T / tau))
{
fprintf (plotter, «[[% g,% g,% g], [% g,% g,% g], [% g,% g,% g], [% g,% g,% g]] », i * h, j * tau, u [i] [j], (i + 1) * h, j * tau, u [i + 1] [j], i * h, (j + 1) * tau, u [i] [j + 1], (i + 1) * h, (j + 1) * tau, u [i + 1] [j + 1]);
if (i> = n - 1 & & j> = floor (T / tau) - 1)
{
}
else
{
fprintf (plotter ,»,»);
}
}
}
printf ("\ n");
}
fclose (myfile);
fclose (plotter);
printf ("\ nОсь x розташована горизонтально; вісь t розташована вертикально і спрямована вниз»);
printf ("Крок по осі x дорівнює% g; крок по осі t дорівнює% g. \ n», h, tau);
printf ("\ nДля виходу натисніть ENTER ...»);
while (getch ()! = 13);
}
/ / - / /
float f (float x, float t)
{
return x * t;
}
/ / - / /
float mu_1 (float t)
{
return 2.1 + t;
}
/ / - / /
float mu_2 (float t)
{
return 3.2 * (t + 1 / 2.71828);
}
/ / - / /
float phi (float x)
{
return (1.1 * pow (x, 2) + 2.1) * exp (-x);
}
Існують спеціальні оцінки похибки для вирішення задач методом сіток. Однак ці оцінки містять максимуми модулів похідних шуканого рішення, тому користуватися ними вкрай незручно, однак ці теоретичні оцінки хороші тим, що з них видно: якщо необмежено подрібнювати сітку, то послідовність рішень буде сходитися рівномірно до точного розв'язання. Тут ми зіткнулися з проблемою збіжності методу сіток. При використанні методу сіток ми повинні бути впевнені, що, необмежено згущуючи сітку, можемо отримати рішення, як завгодно близьке до точного.
Отже, на прикладі розв'язку крайової задачі для диференціального рівняння параболічного типу розглянемо основні принципи методу сіток. Відзначимо, що якщо при вирішенні різницевої задачі невеликі помилки в початкових і крайових умовах (або в проміжних результатах) не можуть привести до великих відхилень шуканого рішення, то говорять, що задача поставлена коректно в сенсі стійкості за вхідними даними. Різницеву схему називають стійкою, якщо обчислювальна похибка необмежено не зростає. В іншому випадку схема називається нестійкою.
1.4 Доказ стійкості різницевої схеми
Нехай
і граничним умовам
Тут
Розглянемо похибка
Похибка
(В силу лінійності рівняння (1.14)), а також наступними граничними та початковими умовами:
Приватне рішення рівняння (1.23) будемо шукати у вигляді
Тут числа
При цілому
Підставимо рівняння (1.26) в рівняння (1.24). При цьому отримаємо:
або
Вираз у квадратних дужках одно
Підставляючи цей вираз у попереднє рівняння замість виразу у квадратних дужках і проводячи скорочення на
звідки знаходимо
Таким чином, згідно з рівнянням (1.26), отримуємо лінійно-незалежні рішення рівняння (1.23) у вигляді
Зауважимо, що це приватне рішення задовольняє однорідним крайовим умовам (1.24). Лінійна комбінація цих приватних рішень також є рішенням рівняння (1.23):
причому
У результаті отримуємо систему рівнянь
містить
Для стійкості досліджуваної різницевої схеми необхідно, щоб за будь-яких значеннях коефіцієнтів
Аналізуючи (1.28) бачимо, що ця нерівність виконується для будь-яких значень параметра
2. Реалізація методу
2.1 Розробка програмного модуля
Поставлено мету: розробити програмний продукт для знаходження наближеного рішення параболічного рівняння:
в області
задовольняє умовам
Розіб'ємо область
де
Замінивши у кожному вузлі похідні кінцево-різницевими відносинами по неявній схемі, отримаємо систему види:
Перетворивши її, отримаємо:
де
У граничних вузлах
У початковий момент
Ця різницева схема стійка при будь-якому
де
Аналогічно
Підставивши значення (1.35) в (1.32) отримаємо:
Звідки
З порівняння (1.35) і (1.37) видно, що
Для
Звідки
або
Звідки, використовуючи (1.35), отримаємо:
Використовуючи даний метод, ми всі обчислення проведемо в наступному порядку для всіх
1) Знаючи значення функції
2) Знайдемо
3) Знайдемо
4) Знайдемо значення шуканої функції на
2.2 Опис логіки програмного модуля
Лістинг програми наведено в додатку 1. Нижче будуть описані функції програмного модуля та його призначення.
Функція main () є базовою. Вона реалізує алгоритм методу сіток, описаного в попередніх розділах роботи.
Функція f (x, y) представляє собою вільну функцію двох змінних диференціального рівняння (1.29). Як аргумент на неї передаються два дійсних числа з плаваючою точкою типу float. На виході функція повертає значення функції
Функції mu_1 (t) і mu_2 (t) представляють собою крайові умови. У них передається по одному аргументу (t) дійсного типу (float).
Функція phi () є відповідальною за початковий умови.
У функції main () визначено такі константи:
Варіюючи
Програма забезпечена трьома механізмами виведення результатів роботи: на екран у вигляді таблиці, в текстовий файл, а також у файл списку математичного пакета Waterloo Maple. Це дозволяє наочно уявити отримане рішення.
Програма написана на мові програмування високого рівня Borland C + + 3.1 у вигляді додатку MS-DOS. Забезпечується повна сумісність програми з усіма широко відомими операційними системами корпорації Майкрософт: MS-DOS 5.x, 6.xx, 7.xx, 8.xx, Windows 9x/Me/2000/NT/XP.
2.3 Приклад роботи програми
В якості прикладу розглянемо чисельне рішення наступного диференціального рівняння параболічного типу:
в області
задовольняє умовам
Задавши прямокутну сітку з кроком осі
2.10 1.91 1.76 1.63 1.53 1.44 1.37 1.31 1.26 1.22 1.18
2.11 1.75 1.23 1.20 1.15 1.10 1.07 1.04 1.04 1.07 1.21
2.12 1.61 0.95 0.96 0.93 0.91 0.90 0.90 0.94 1.03 1.24
2.13 1.51 0.79 0.81 0.81 0.80 0.81 0.83 0.89 1.03 1.27
2.14 1.45 0.69 0.73 0.74 0.74 0.76 0.80 0.88 1.04 1.31
2.15 1.41 0.64 0.69 0.70 0.71 0.74 0.79 0.89 1.05 1.34
У таблиці вісь x розташована горизонтально, а вісь t розташована вертикально і спрямована вниз.
На виконання програми на середньостатистичному персональному комп'ютері витрачається час, що дорівнює декількох мілісекунд, що говорить про високу швидкість алгоритму.
Докладно вихідний файл output.txt, що містить таблицю значень функції
Висновок
У роботі було розглянуто метод сіток рішення параболічних рівнянь в приватних похідних. Розкрито основні поняття методу, апроксимація рівняння та граничних умов, досліджено розв'язність і збіжність одержуваної системи різницевих рівнянь.
На підставі вивченого теоретичного матеріалу була розроблена програмна реалізація методу сіток, проаналізовано її збіжність і швидкодію, проведений тестовий розрахунок, побудований графіки отриманого чисельного рішення.
Список джерел
1. Березін І.С., Жидков Н.П. Методи обчислень. Т.2. - М.: Физматгиз, 1962.
2. Тихонов О.М., Самарський О.А. Рівняння математичної фізики. - М.: Наука, 1972.
3. Пірумов Н.Р. Чисельні методи. - М.: Видавництво МАІ, 1998.
4. Каліткін М.М. Чисельні методи. - М.: Наука, 1976.
Додаток
Текст програми
/ / - / /
# Include <stdio.h>
# Include <conio.h>
# Include <math.h>
void main (void);
float f (float x, float t);
float mu_1 (float t);
float mu_2 (float t);
float phi (float x);
/ / - / /
void main (void)
{
clrscr ();
FILE * myfile;
FILE * plotter;
float a [120] [120];
float b [120] [120];
float u [120] [120];
float T = 0.05;
float l = 1;
float h = 0.1;
float tau = 0.01;
int n, i, j, k;
float s = pow (h, 2) / tau;
n = ceil (l / h);
for (i = 0; i <= 119; i + +)
{
for (j = 1; j <= 119; j + +)
{
u [i] [j] = 0;
a [i] [j] = 0;
b [i] [j] = 0;
}
}
for (i = 0; i <= n; i + +)
{
u [i] [0] = phi (i * h);
}
for (j = 0; j <= floor (T / tau); j + +)
{
u [0] [j] = mu_1 (tau * j);
u [n] [j] = mu_2 (tau * j);
}
for (j = 0; j <= floor (T / tau); j + +)
{
a [1] [j + 1] = 1 / (2 + s);
for (i = 2; i <= n - 1; i + +)
{
a [i] [j + 1] = 1 / (2 + s - a [i - 1] [j + 1]);
}
b [1] [j + 1] = mu_1 ((j + 1) * tau) + s * u [1] [j] + pow (h, 2) * f (h, (j + 1) * tau) ;
for (i = 2; i <= n - 1; i + +)
{
b [i] [j + 1] = a [i - 1] [j + 1] + s * u [i] [j] + pow (h, 2) * f (i * h, (j + 1) * tau);
}
u [n] [j + 1] = mu_2 ((j + 1) * tau);
for (k = 1; k <= n - 1; k + +)
{
u [n - k] [j + 1] = a [n - k] [j + 1] * (b [n - k] [j + 1] + u [n - k + 1] [j + 1] );
}
}
myfile = fopen («output.txt», «w +»);
plotter = fopen («3dplot.txt», «w +»);
fprintf (myfile, «Таблиця значень функції u = u (x, t) в області D = {0 <= X <=% g, 0 <= T <=% g}: \ n», l, T);
printf («Значення функції u (x, t) в області D = {0 <= X <=% g, 0 <= T <=% g}: \ n \ n», l, T);
for (j = 0; j <= floor (T / tau); j + +)
{
for (i = 0; i <= n; i + +)
{
printf ("% .2 f», u [i] [j]);
fprintf (myfile, «u (% g) (% g) =% g; \ n», i * h, j * tau, u [i] [j]);
if (i <n & & j <floor (T / tau))
{
fprintf (plotter, «[[% g,% g,% g], [% g,% g,% g], [% g,% g,% g], [% g,% g,% g]] », i * h, j * tau, u [i] [j], (i + 1) * h, j * tau, u [i + 1] [j], i * h, (j + 1) * tau, u [i] [j + 1], (i + 1) * h, (j + 1) * tau, u [i + 1] [j + 1]);
if (i> = n - 1 & & j> = floor (T / tau) - 1)
{
}
else
{
fprintf (plotter ,»,»);
}
}
}
printf ("\ n");
}
fclose (myfile);
fclose (plotter);
printf ("\ nОсь x розташована горизонтально; вісь t розташована вертикально і спрямована вниз»);
printf ("Крок по осі x дорівнює% g; крок по осі t дорівнює% g. \ n», h, tau);
printf ("\ nДля виходу натисніть ENTER ...»);
while (getch ()! = 13);
}
/ / - / /
float f (float x, float t)
{
return x * t;
}
/ / - / /
float mu_1 (float t)
{
return 2.1 + t;
}
/ / - / /
float mu_2 (float t)
{
return 3.2 * (t + 1 / 2.71828);
}
/ / - / /
float phi (float x)
{
return (1.1 * pow (x, 2) + 2.1) * exp (-x);
}