Задачи Коши для дифференциальных уравнений высших порядков

Задача Коши для дифференциального уравнения n -го порядка

 

,                   (7.21)

 

           (7.22)

 

легко сводится к задаче Коши для системы дифференциальных уравнений с помощью замены переменных

 

                      (7.23)

 

Учитывая (7.23) из уравнения (7.21) получим систему дифференциальных уравнений

                                     (7.24)

Начальные условия (7.22) для функций z_l переписываются в виде

.                          (7.25)

Запишем для полученной системы метод Рунге-Кутта

.    (7.26)

Для вычисления коэффициентов K_{l ,1}, K_{l ,2}, K_{l ,3} и K_{l ,4} имеем следующие формулы:

 

Пример 7.5. Решить задачу Коши для дифференциального уравнения третьего порядка

 

,

Решение в Excel. Точное решение данной задачи известно . Преобразуем эту задачу с помощью замены  к задаче для системы уравнений:

 

 

Точным решением системы теперь будут функции

 

.

 

Выберем шаг h = 0,1. В программе Excel мы используем функции и алгоритм решения, созданные для решения примера 7.4.

Откроем файл примера 7.4 и внесем необходимые изменения. Результаты расчетов представлены в таблице 7.10.

1) Заменить функции f_1(), f_2() и f_3() на следующие:

 

Function f_1(ByVal x, ByVal y_1, ByVal y_2, ByVal y_3)

f_1 = 2 * x * y_1 / y_2

End Function

Function f_2(ByVal x, ByVal y_1, ByVal y_2, ByVal y_3)

f_2 = 8 * y_3 / (3 * y_2)

End Function

Function f_3(ByVal x, ByVal y_1, ByVal y_2, ByVal y_3)

f_3 = 3 * y_3 / y_1

End Function

 

Функции для вычисления коэффициентов K_{l , I} не изменяем.

2) Изменим начальные условия в ячейках B 3, C 3, D 3 на 0, 1, 0.

В столбце значений переменной x_i вводим числа 0; 0,1; 0,2; …; 1.

В ячейках E 3, F 3, G 3 вводим соответственно формулы =SIN(A3), =COS(A3) и =-SIN(A3).

 

Таблица 7.10

	A	B	C	D	E	F	G
1	h=	0,1			Точные решения
2	xi	y1i	y2i	y3i	y1=sinx	y2=cosx	y3=-sinx
3	0	0	1	0	0	1	0
4	0,1	0,099833	0,995004	-0,09983	0,09983	0,995	-0,099833
5	0,2	0,198669	0,980067	-0,19867	0,19867	0,98007	-0,198669
6	0,3	0,29552	0,955337	-0,29552	0,29552	0,95534	-0,29552
7	0,4	0,389418	0,921061	-0,38942	0,38942	0,92106	-0,389418
8	0,5	0,479425	0,877583	-0,47943	0,47943	0,87758	-0,479426
9	0,6	0,564642	0,825336	-0,56464	0,56464	0,82534	-0,564642
10	0,7	0,644217	0,764843	-0,64422	0,64422	0,76484	-0,644218
11	0,8	0,717356	0,696707	-0,71736	0,71736	0,69671	-0,717356
12	0,9	0,783326	0,621611	-0,78333	0,78333	0,62161	-0,783327
13	1	0,84147	0,540303	-0,84147	0,84147	0,5403	-0,841471

 

Сравнивая столбцы B и E, видим, что метод Рунге-Кутта дает хорошие приближения к значениям точного решения y1=sinx.

Краевые задачи для обыкновенных дифференциальных уравнений

Примерами краевых задач для обыкновенных дифференциальных уравнений являются следующие задачи:

1) Найти функцию u (x), которая удовлетворяет на отрезке [ a, b ] уравнению

,                                             (7.27)

а на концах отрезка условиям

 

u (a) = u (b) = 0.                                              (7.28)

 

Задача (7.27), (7.28) имеет следующее содержание. Между точками
x = a и x = b натянута упругая струна, находящаяся под действием внешней изгибающей нагрузки.   f (x) — величина нагрузки, а u (x) — прогиб струны в безразмерных единицах.

2) Двухточечная краевая задача для линейного дифференциального уравнения второго порядка

 

                              (7.29)

 

                                           (7.30)

 

Метод прогонки

Рассмотрим частный случай задачи (7.29), (7.30)

 

                              (7.31)

 

                                        (7.32)

 

Введем сетку: . Обозначим через u_i приближенные значения u (x) в узлах сетки. Рассмотрим уравнение (7.31) во внутренних узлах  и заменим производную второго порядка разностной формулой (5.11):

 

.                                   (7.33)

 

Тогда из (7.31), (7.32) получим для определения u_i систему линейных уравнений

 

                    (7.34)

 

Система (7.34) при p (x) ≥ 0 имеет решение. Доказательство этого утверждения можно найти в [1].

Система (7.34) представляет собой систему линейных уравнений с трехдиагональной матрицей

 

                    (7.35)

 

и для её решения применим метод прогонки, который фактически является методом исключения неизвестных Гаусса.

1. Прямой ход прогонки. Запишем первое уравнение (7.35) в виде

 

                      (7.36)

 

Подставив (7.36) во второе уравнение (7.35) и упростив выражения, увидим, что можно для u_i получить формулы

 

(7.37)

 

Из последнего уравнения (7.35), учитывая (7.37) при , получим

 

                     (7.38)

 

2. Обратный ход прогонки. После вычисления прогоночных коэффициентов можно найти значения решения задачи. Формула (7.38) дает значение :

 

                                                      (7.39)

 

Остальные значения вычисляем в обратном порядке

 

                              (7.40)

 

Пример 7. 6. Решить краевую задачу методом прогонки

 

.

 

Решение в программе Excel. Выберем N = 10, т.е. шаг h = π/10 ≈ 0,31415926.

Составим таблицу значений аргумента x, функций p (x) и f (x), вычислим прогоночные коэффициенты и приближенные значения искомого решения по формулам (7.36) — (7.40). В таблице 7.11 приведен фрагмент листа Excel с результатами вычислений.

1) Заполним столбцы A 1: A 15, B 1: B 15 и строку A 4: H 4 как в таблице 7.11. При заполнении столбца значений переменной x можно ввести 0 в ячейку B 5, 3,1415926 в ячейку B 6, затем выделить B 5: B 6 и маркером заполнения протянуть до ячейки B 15.

2) Присвоим имена a, b и h ячейкам B 1, B 2, B 3 соответственно с помощью команды «Вставка — Имя — Присвоить».

Приведем содержимое для остальных ячеек, которые содержат формулы.

3) В столбце для значений p (x), т.е. в ячейках C 5: C 16 записываем единицы, так как p (x) = 1.

4) В диапазоне D 5: D 16 вводим формулу =-2*SIN(x),  так как
f (x) = – 2sin x. Для этого достаточно ввести эту формулу в ячейку D 5 и протянуть D 5 маркером заполнения до D 15.

5) В ячейку E 5 введем число 0, а в ячейку E 6 — формулу
=1/(2+C6*h^2-E5) и протянем ячейку E 6 маркером заполнения до E 13. В ячейку E 14 запишем число 0.

6) В ячейку F 5 введем число 0, в ячейку F 6 вводим формулу
=(a-D6*h^2)/(2+C6*h^2-E5), и в ячейку F 7 введём =(F6-D7*h^2)/(2+C7*h^2-E6) и протянем ячейку F7 маркером заполнения до F 14.

7) Вычислим решение. В ячейку G 14 введем =F14, в ячейку G 13 введем формулу =E13*G14+F13 и протянем G 13 вверх до ячейки G 6.

8) Для проверки правильности вычислений вычислим в столбце H 5: H 15 значения точного решения задачи . В ячейку H 5 введем формулу =sin(x) и протянем ячейку H 5 до H 15.

Как видим из таблицы, полученные приближенные значения имеют абсолютную погрешность не более 0,01.

Таблица 7.11

	A	B	C	D	E	F	G	H
1	a=	0
2	b=	0
3	h=	0,31415926
4	i	x	p(x)	f(x)	alfa	beta	ui	sin(x)
5	0	0	1	0	0	0		0
6	1	0,31415926	1	-0,61803	0,476486	0,029064	0,310289	0,309017
7	2	0,62831852	1	-1,17557	0,616443	0,089439	0,590204	0,587785
8	3	0,94247778	1	-1,61803	0,674649	0,168077	0,812347	0,809017
9	4	1,25663704	1	-1,90211	0,702224	0,249857	0,954971	0,951057
10	5	1,5707963	1	-2	0,71609	0,320271	1,004116	1
11	6	1,88495556	1	-1,90211	0,723272	0,367423	0,954971	0,951057
12	7	2,19911482	1	-1,61803	0,727048	0,383239	0,812347	0,809017
13	8	2,51327408	1	-1,17557	0,72905	0,363988	0,590204	0,587785
14	9	2,82743334	1	-0,61803	0	0,310289	0,310289	0,309017
15	10	3,1415926	1	-1,1E-07				5,36E-08

Рассмотрим ещё один пример простой краевой задачи, решение которой получим с помощью той же таблицы Excel.

Пример 7. 7. Решить краевую задачу методом прогонки

 

 

Решение в программе Excel. Выберем шаг h = 0,1, т.е. разделим отрезок [1, 2] на N = 10 частей. Откроем файл с решением примера 7.4 и внесем в таблицу изменения. Эти изменения не касаются столбцов E, F и G, в которых записаны формулы (7.27) — (7.31) метода прогонки.

Так как a = 2, b = 12 и h = 0,1, в ячейки B 1: B 3 вводим числа 2, 12 и 0,1.

Отметим те изменения, которые не видны явно в таблице 7.12.

В столбце для значений p (x), т.е. в ячейках C 5: C 16 записываем нули, так как p (x) = 0.

В ячейку D 5 вводим формулу =6*x+2,  так как f (x) = 6 x + 2, и протянем ячейку D 5 до D 15.

В ячейку H 5 введем формулу =x^3+x^2 и протянем ячейку H 5 до H 15.

Результаты приведены в таблице 7.12.

Таблица 7.12

	A	B	C	D	E	F	G	H
1	a=	2
2	b=	12
3	h=	0,1
4	i	x	p(x)	f(x)	alfa	beta	ui	x^3+x^2
5	0	1	0	8	0	0		2
6	1	1,1	0	8,6	0,5	0,957	2,54100000	2,54100000
7	2	1,2	0	9,2	0,666667	0,576667	3,16800000	3,16800000
8	3	1,3	0	9,8	0,75	0,359	3,88700000	3,88700000
9	4	1,4	0	10,4	0,8	0,204	4,70400000	4,70400000
10	5	1,5	0	11	0,833333	0,078333	5,62500000	5,62500000
11	6	1,6	0	11,6	0,857143	-0,03229	6,65600000	6,65600000
12	7	1,7	0	12,2	0,875	-0,135	7,80300000	7,80300000
13	8	1,8	0	12,8	0,888889	-0,23378	9,07200000	9,07200000
14	9	1,9	0	13,4	0	10,469	10,46900000	10,46900000
15	10	2	0	14				12

 

Точное совпадение решения системы уравнений (7.34) со значениями решения краевой задачи (7.31), (7.32) объяснятся следующим обстоятельством. Разностная формула для производной второго порядка (7.33) дает точное значение (!) для кубических функций.