Лабораторная работа 1

Численные методы решения

обыкновенных дифференциальных уравнений (4 часа)

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

Обыкновенным дифференциальным уравнением называется равенство

,                                       (1)

в котором  - независимая переменная, изменяющаяся в некотором отрезке , а  - неизвестная функция y
(
x
) и ее первые n производные.

Число называется порядком уравнения.

Задача заключается в нахождении функции y, удовлетворяющей равенству (1). Более того, не оговаривая это отдельно, будем предполагать, что искомое решение обладает той или иной степенью гладкости, необходимой для построения и «законного» применения того или иного метода.

Различают два типа обыкновенных дифференциальных уравнений

 - уравнения без начальных условий

 - уравнения с начальными условиями.

Уравнения без начальных условий - это уравнение вида (1).

Уравнение с начальными условиями - это уравнение вида (1), в котором требуется найти такую функцию , которая при некотором  удовлетворяет следующим условиям:

,

т.е. в точке  функция   и ее первые  производных принимают наперед заданные значения.

Задачи Коши

При изучении способов решения дифференциальных уравнений приближенными методами основной задачей считается задача Коши.

Рассмотрим наиболее популярный метод решения задачи Коши – метод Рунге-Кутта. Этот метод позволяет строить формулы расчета приближенного решения практически любого порядка точности.

Выведем формулы метода Рунге-Кутта второго порядка точности. Для этого решение представим куском ряда Тейлора, отбрасывая члены с порядком выше второго. Тогда приближенное значение искомой функции в точке x
₁
можно записать в виде:

 (2)

Вторую производную y
"(
x
₀
) можно выразить через производную функции f
(
x
,
y
), однако в методе Рунге-Кутта вместо производной используют разность



соответственно подбирая значения параметров

Тогда (2) можно переписать в виде:

y
₁
=
y
₀
+
h
[
β

f
(
x
₀
,
y
₀
) +
α

f
(
x
₀
+
γh
,
y
₀
+
δh
)],                           (3)

где α
,
β
,
γ
и  δ
– некоторые параметры.

Рассматривая правую часть (3) как функцию аргумента h

, разложим ее по степеням h
:

y
₁
=
y
₀
+(
α
+
β
)
h

f
(
x
₀
,
y
₀
) +
αh
²
[
γ

f_x
(
x
₀
,
y
₀
) +
δ

f_y
(
x
₀
,
y
₀
)],


и выберем параметры α
,
β
,
γ
и  δ
так, чтобы это разложение было близко к (2). Отсюда следует, что

α
+
β
=1,
αγ
=0,5,
α

δ
=0,5
f
(
x
₀
,
y
₀
).

С помощью этих уравнений выразим β
,
γ
и  δ
через параметры α
, получим

y
₁
=
y
₀
+
h
[(1 -
α
)
f
(
x
₀
,
y
₀
) +
α

f
(
x
₀
+,
y
₀
+
f
(
x
₀
,
y
₀
)],           (4)

0 <
α
≤ 1.

Теперь, если вместо (x
₀
,
y
₀) в (4) подставить (x
₁
,
y
₁), получим формулу для вычисления y
₂
– приближенного значения искомой функции в точке x
₂
.

В общем случае метод Рунге-Кутта применяется на произвольном разбиении отрезка [x
₀
,
X
] на  n
частей, т.е. с переменным шагом

x₀, x₁, …,x_n;  h_i = x_i+1 – x_i,  x_n = X.                                          (5)

Параметры α
выбирают равными 1 или 0,5. Запишем окончательно расчетные формулы метода Рунге-Кутта второго порядка с переменным шагом для  α
=1:

y_i+1=y_i +h_i f(x_i +

, y_i+

f(x_i, y_i)),                       (6.1)

i
= 0,

1,…,
n
-1.


и α
=0,5:

y_i+1=y_i + [f(x_i , y_i) + f(x_i+ h_i, y_i+ h_i f(x_i , y_i))],                      (6.2)

i
= 0, 1,…,
n
-1.


Наиболее употребляемые формулы метода Рунге-Кутта – формулы четвертого порядка точности:

y_i+1=y_i +

(k₁+ 2k₂+ 2k₃+ k₄),

k₁=f(x_i, y_i),  k₂= f(x_i +

, y_i+

k₁),                                 (7)


k₃= f(x_i +

, y_i+

k₂),  k₄= f(x_i +h, y_i+hk₃).

Для метода Рунге-Кутта применимо правило Рунге для оценки погрешности. Пусть y
(
x
;
h
) – приближенное значение решения в точке x
, полученное по формулам (6.1), (6.2) или (7) с шагом h
, а p
– порядок  точности соответствующей формулы. Тогда погрешность R
(
h
) значения y
(
x
;
h
) можно оценить, используя приближенное значение y
(
x
; 2
h
) решения в точке x
, полученное с шагом 2h
:

                                    (8)

где p
=2 для формул (6.1) и (6.2) и  p
=4 для  (7).

Уточненное решение пишем в виде

.                             (9)

В алгоритмах с автоматическим выбором шага предварительно задают погрешность в виде положительного параметра ε, и на каждом этапе вычисления следующего значения y_i
+1
подбирают такой шаг h
, при котором выполняется неравенство

 ,                              (10)

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

 x  (x₀, X),                                  (11)
y₁(x₀)=y_1,0,  y₂(x₀)=y_2,0,…, y_m(x₀)=y_m,0 .                                   (12)


Приведем для задачи (11), (12) расчетные формулы метода Рунге-Кутта четвертого порядка. Пусть требуется найти систему m
функций y
₁
(
x
),
y
₂
(
x
),…,
y_m
(
x
), удовлетворяющих в интервале (x
₀
,
X
) дифференциальным уравнениям (11), а в точке x
₀
– начальным условиям (12). Предположим, что отрезок [x
₀
,
X
] разбит на N частей:

x_i
=
x
₀
+
i

h_i, 

Тогда каждую l-ю функцию y_l
(
x
) можно приближенно вычислять в точках x_i
+1
по формулам Рунге-Кутта

K_l,1=f_l(x_i, y_1,i, y_2,i,…,y_m,i), i=1, 2, …, m,

K_l,2=f_l(x_i + , y_1,i + K_1,1, y_2,i + K_2,1,…,y_m,i + K_m,1), i=1, 2, …, m,

K_l,3=f_l(x_i + , y_1,i + K_1,2, y_2,i + K_2,2,…,y_m,i + K_m,2), i=1, 2, …, m,
   
         (13)


K_l,4=f_l(x_i + h, y_1,i + hK_1,3, y_2,i + hK_2,3,…,y_m,i + hK_m,3), i=1, 2, …, m,

Y_l,i+1 = y_l,i +( K_l,1 + 2 K_l,2 + 2 K_l,3 + K_l,4), i=1, 2, …, m,

Здесь через y_l
,
i обозначается приближенное значение функции y_l
(
x
) в точке x_i
.

Обратите внимание на порядок вычислений по формулам (13). На каждом шаге сначала вычисляются коэффициенты K_l
,
i
 в следующем порядке:

K_1,1, K_2,1,…, K_m,1,

K_1,2, K_2,2,…, K_m,2,

K
_1,3
, K
_2,3
,…,
K_m
,3
,

K
_1,4
, K
_2,4
,…,
K_m
,4
,

и лишь затем приближенные значения функций y
_1,
i
+1
,
y
_2,
i
+1
,…,
y_m
,
i
+1
.

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

y⁽ⁿ⁾=f(x, y, y', …, y^(n-1)), x  (x₀, X),                                (14)

y(x₀)=y₀,  y'(x₀)=y_1,0, …, y^(n-1)(x₀)=y_n-1,0                                   (15)

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

z
₀
=
y
,
z
₁
=
y
',…
,

z_n-1
=
y^(n-1)
.
                         (16)

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

                             (17)

Начальные условия (15) для функций z_l
переписываются в виде
z₀(x₀)= y₀, z₁(x₀)= y_1,0,…, z_n-1(x₀)= y
_п
-1,0
.                         (18)
Запишем для полученной системы метод Рунге-Кутта:

z_l,i+1= z_l,i +

(K_l,1+ 2K_l,2+ 2K_l,3+ K_l,4),                                    (19)

i
=0, 1, …,
N
, 
l
=0, 1, …,
n
-1.
Для вычисления коэффициентов K_l
,1
,
K_l
,2
,
K_l
,3
и K_l
,4 имеем следующие формулы:

K
_0,1
=
z
_1,
i
,


K_1
,1
=
z_2
,
i,


…………

K_n-1,1= f(x_i, z_0,i, z_1,i,…, z_n-1,i,),

K
_0,2
=
z
_1,
i
+
K
_1,1
,

K_1
,2
=
z_2
,
i
+
K_2
,1
,


…………………

K_n-1,2= f(x_i+ , z_0,i+ K_0,1, z_1,i+ K_1,1,…, z_n-1,i+ K_n-1,1),

K_0,3= z
_1,
i
+
K
_1,
2
,


K_1,3= z_2,i+ K_2,2,

……………………

K_n-1,3= f(x_i+ , z_0,i+ K_0,2, z_1,i+ K_1,2,…, z_n-1,i+ K_n-1,2),

K_0,4= z_1,i+ hK_1,3,

K_1,4= z_2,i+ hK_2,3,

……………………

K_n-1,4= f(x_i+ h, z_0,i+ hK_0,2, z_1,i+ hK_1,2,…, z_n-1,i+ hK_n-1,2).

Задания лабораторной работы 1

1.     Написать файл-функции для решения поставленных далее задач.

2.     Сохранить их в отдельных m-файлах (среда Матлаб)

3. Выполнить и оформить в виде отчета поставленные далее задачи.
Задача №
1. Решить задачу Коши на отрезке [x₀,X] методом Рунге-Кутта четвертого порядка, применяя деление отрезка на N частей. Оценить погрешность.

Варианты заданий в табл.1.

Табл.1.

№ варианта	Уравнение	Начальное условие	[ x 0 , X ]	N
1	y'(x)=sin(xy2)	y(0)=1	[0,2]	10
2	y'(x)=cos(x) + y2	y(0)=2	[0,2]	20
3	y'(x)= cos(xy2)	y(0)=3	[0,2]	30
4	y'(x)=sin	y(0)=1	[0,2]	40
5	y'(x)=tg	y(0)=2	[0,2]	50
6	y'(x)=x + y2	y(1)=3	[1,2]	10
7	y'(x)=	y(1)=1	[1,2]	20
8	y'(x)=cos	y(1)=2	[1,2]	30
9	y'(x)=sin ( x )	y(1)=3	[1,2]	40
10	y'(x)=	y(1)=1	[1,2]	50
11	y'(x)=x ln(1+y2)	y(1)=2	[1,3]	10
12	y'(x)=y cos(x+y2)	y(1)=3	[1,3]	20
13	y'(x)=ex x+y2	y(1)=1	[1,3]	30
14	y'(x)=sin(x(1+y2))	y(1)=2	[1,3]	40
15	y'(x)=lg	y(1)=3	[1,3]	50
16	y'(x)=x+y2 3x	y(-1)=1	[-1,1]	10
17	y'(x)=|x-y|(1+x2+y2)	y(-1)=2	[-1,1]	20
18	y'(x)=	y(-1)=3	[-1,1]	30
19	y'(x)=x+	y(-1)=1	[-1,1]	40
20	y'(x)=	y(-1)=2	[-1,1]	50
21	y'(x)=	y(0)=3	[0,π]	10
22	y'(x)=sin(x) ln(1+y2)	y(0)=1	[0,π]	20
23	y'(x)=sin(y) cos(x+y2)	y(0)=2	[0,π]	30
24	y'(x)=ex sin(y)+x2 ey	y(0)=3	[0,π]	40
25	y'(x)= cos(x) (x+y2)	y(0)=1	[0,π]	50
26	y'(x)=	y( π/2 )=2	[π/2,π]	10
27	y'(x)=x 2y+y 2x	y( π/2 )=1	[π/2,π]	20
28	y'(x)= |x - y| cos(x2 + y2)	y( π/2 )=3	[π/2,π]	30
29	y'(x)=	y( π/2 )=2	[π/2,π]	40
30	y'(x)=(y + x )	y( π/2 )=3	[π/2,π]	50

Задача №
2. Решить задачу Коши для дифференциального уравнения сведением к задачи Коши для системы уравнений первого порядка.

Табл.2.

№ варианта	Дифференциальное уравнение	Начальное условие	[ x 0 , X ]	N
1	y(x)=x y(x)+ sin(x)	y (0)=1, y' (0) =2	[0,2]	10
2	y"'(x)=2x2 y(x) y"(x)	y(0)=2, y' (0) =2, y"(0)=1	[0,2]	20
3	y"(x) – 3cos(x) y(x)=tg(x)	y(0)=3, y' (0) =2	[0,2]	30
4	"'y(x)=x y'(x)	y(0)=1, y' (0) =1, y"(0)=1	[0,2]	40
5	y"'(x)=-cos(x) y"(x) – y(x) sin(x)	y(0)=2, y' (0) =2 , y"(1)=1	[0,2]	50
6	y"(x)– sin(x) y(x)=sin(x)	y(1)=3, y' ( 1 ) =1	[1,2]	10
7	y"(x) – 2x2 y(x)=cos(x)	y(1)=1, y' ( 1 ) =1	[1,2]	20
8	y"'(x)=(x – 1) y(x) + x y"(x)	y(1)=2, y' ( 1 ) =1, y"(1)=1	[1,2]	30
9	y"(x) - sin(x) y(x)=sin3(x)	y(1)=3, y' ( 1 ) =1	[1,2]	40
10	y"'(x)=x y(x) - sin(x) y'(x)	y(1)=1, y' ( 1 ) =1, y"(1)=1	[1,2]	50
11	y"(x)-cos(x) y(x)=x	y(1)=2, y' ( 1 ) =1	[1,3]	10
12	y"'(x) – 2x2 y(x)=x2	y(1)=3, y' (0) =1, y"(0)=1	[1,3]	20
13	y"(x) - lgx y(x)=2x	y(1)=1, y' ( 1 ) =1	[1,3]	30
14	y"'(x) - 2|sin(x)| y'(x)=3x3	y(1)=2, y' ( 1 ) =1, y"(1)=1	[1,3]	40
15	y"(x) – 2lnx y(x)=1+x	y(1)=3, y' ( 1 ) =1	[1,3]	50
16	y"'(x) - |cos(x)| y(x)=x	y(-1)=1, y' ( -1 ) =1, y"(-1)=1	[-1,1]	10
17	y"(x) - 2|x| y(x)=cos2(x)	y(-1)=2, y' ( 1 ) =1	[-1,1]	20
18	y"'(x) - y(x)=e2x	y(-1)=3, y' ( -1 ) =1, y"(-1)=1	[-1,1]	30
19	y"(x) – ln(1+x2) y(x)=sin(2x)	y(-1)=1, y' ( 1 ) =1	[-1,1]	40
20	y"'(x) – sin|x| y(x)=sin(x)	y(-1)=2, y' ( -1 ) =1, y"(-1)=1	[-1,1]	50
21	y"(x) - 2y(x)=sin(x)	y(0)=3, y' (0) =2	[0,π]	10
22	y"'(x)=3y(x)+y"(x) cos(x)	y(0)=1, y' (0) =1, y"(0)=1	[0,π]	20
23	y"(x) - 2x y(x)=x3	y(0)=2, y' (0) =2	[0,π]	30
24	y"'(x) - x y(x)=x4y'(x)	y(0)=3, y' (0) =1, y"(0)=1	[0,π]	40
25	y"(x) - 2x2 y(x)=x2	y(0)=1, y' (0) =2	[0,π]	50
26	y"'(x)=cos(x) y(x)+ex y"(x)	y( 2 )=2, y' ( 2 ) =1, y"(2)=1	[2,π]	10
27	y"(x) - 2x2 y(x)=2x ex	y( 2 )=3, y' (0) =2	[2,π]	20
28	y"'(x) - 5y"(x)=32x	y( 2 )=1, y' ( 2 ) =1, y"(2)=1	[2,π]	30
29	y"(x) - 2sin(x) y(x)=sin(3x)	y( 2 )=2, y' (0) =2	[2,π]	40
30	y"'(x) - lnx y'(x)=1	y( 2 )=3, y' ( 2 ) =1, y"(2)=1	[2,π]	50

Задача №
3.

Найти методом Рунге-Кутта с точностью ε = 10^-8  решение задачи Коши y
'(
x
)=2
x
(1+
y
²
),
y
(0)=0 в точке x
=1.

(Точным решением является функция y
(
x
)=
tg
(
x
²
))
Задача №
4.

Решить методом Эйлера на отрезке [1, 2] задачу Коши

y
'(
x
)=

,
y
(1)=0.

(Точным решением данной задачи является функция y
(
x
)=
tg
(
ln
).


Контрольные вопросы:

1.     Какое уравнение называется обыкновенным дифференциальным уравнением?

2.     Какие методы решения задач для дифференциальных уравнений вы знаете?

3.     В каком случае решение дифференциального уравнения единственно?

4.     Рассказать правило Рунге для оценки погрешности.