Mathcad - ЛР6

ìy¢ = f (t, y, z)

ï

= g(t, y, z)

Задана система двух ОДУ 1-го порядка:

ïz¢

.

í

ïy(t0 ) = y0

ï

îz(t0 ) = z0

Решением системы на сетке ti = t0 + h ×i ,

h

- заданный шаг, является набор чисел (yi , zi ) ,

которые могут быть найдены, исходя из следующих итерационных соотношений:

ì

i+1 = yi + hf (ti , yi , zi )

y

ï

= yi + hg(ti , yi , zi )

ïzi+1

ï

h

íy

i+1

= y

+

( f (t

, y

, z

) + f (t

i+1

, y

i+1

, zi+1 )) .

ï

i

2

i

i

i

ï

h

ïzi+1

= zi

+

(g(ti , yi , zi ) + g(ti+1 , yi+1, zi+1 ))

2

î

Исходные данные:

f(t) := 14×t e- t

A := 0

B := 4.3

b1 := 2

b2 := -0.4

p := 3

q := -4 t0 := 0

i := 0 .. N

ti := t0 + i×h

t0 := 0

f1(t , y, z) := z

f2(t , y, z) := 14×t×e- t

- p×z - q×y

Решение задачи:

Аналитическое решение:

y(t) := - 73 ×e- t ×t - 187 ×e- t - 145 ×e- 4×t× -31519 + 145 ×et × 111140

Решение методом Эйлера-Коши:

EK(f , g, y0, z0 , t0 , B , N) :=

h ¬

B - t0

N

y0 ¬ y0

z0 ¬ z0

for

i Î 0 .. N

y_i+1 ¬ yi + h×f(t0 + i×h, yi , zi)

z_i+1 ¬ zi + h×g(t0 + i×h , yi , zi)

yi+1 ¬ yi +

h

×éf(t0 + i×h, yi , zi)

+ fét0 + (i + 1)×h, y_i+1

, z_i+1ùù

2

ë

ë

ûû

zi+1 ¬ zi +

h

×ég(t0 + i×h, yi , zi) + gét0 + (i + 1)×h, y_i+1

, z_i+1ùù

2

ë

ë

ûû

y

Решение при h=(B-A)/10 :

	N := 10 i := 0 .. N h :=	B - A	ti := t0 + i×h
		N

ek10 := EK(f1 , f2, b1 , b2 , t0 , B , N)

200

150

y(ti)

ek10i 100

50

0

2

4

6

ti

Вычисление N, при котором рассчитанная по правилу Рунге погрешность будет меньше 10- 3

get_N :=

N ¬ 10

e ¬ 1

while e > 10- 3

ek1 ¬ EK(f1, f2 , b1 , b2 , t0 , B , N)

ek2 ¬ EK(f1, f2 , b1 , b2 , t0 , B , 2×N)

for

j Î 0 .. N

Ej

¬

ek22× j - ek1j

3

e ¬ max(E)

N ¬ 2×N

N

d := 64

N := get_N

i := 0 .. N

N = 2560

h :=

B - A

N2 :=

N

N2 = 40 j := 0 .. N2

N

d

EK_res := EK(f1 , f2 , b1 , b2 , t0 , B , N)			ti := t0 + i×h		h2 :=	B - A	t2j := t0 +		j×h2
						N2
200
150
y(t2j)
100
EK_resj×d
50
0	0.5	1	1.5	2	2.5	3	3.5	4	4.5
					t2j

Расчет фактической абсолютной погрешности:

max(y, s, N) := E ¬ 0 max ¬ 0

for i Î 0 .. N

E ¬ si - y(ti)

max ¬ E if E > max max

max(y, EK_res , N) = 3.91066975566901 ´ 10− 4

Вывод: В ходе проведенного численного эксперимента было решено дифференциальное уравнение второго порядка с постоянными коэффициентами методом сведения его к системе линейных дифференциальных уравнений с последующим её решением по методу Эйлера-Коши. Была достигнута практическая абсолютная

погрешность 3.91066975566901 ´ 10− 4 при заданной 10− 3.