Giải gần đúng phương trình vi phân

Chương 6
GIẢI GẦN ĐÚNG
PHƯƠNG TRÌNH VI PHÂN
I. GIẢI GẦN ĐÚNG PTVP CẤP 1 :
Xét bài toán Cauchy : tìm nghiệm y=y(x) của phương trình vi phân với giá trị ban đầu y0
y` = f(x, y), ?x ? [a,b]
y(a) = y0
Các phương pháp giải gần đúng :
Công thức Euler
Công thức Euler cải tiến
Công thức Runge-Kutta
1. Công thức Euler :
Để tìm nghiệm gần đúng của bài toán Cauchy ta chia đoạn [a,b] thành n đoạn nhỏ bằng nhau với bước h = (b-a)/n
xo= a, x1 = x0 +h, ... , xk = x0 + kh, ... , xn = b
Nghiệm gần đúng của bài toán là dãy {yk} gồm các giá trị gần đúng của hàm tại xk
Ta có yk ? y(xk) , k =0, n
Giả sử bài toán có nghiệm duy nhất y(x) có đạo hàm đến cấp 2 liên tục trên [a,b].
Khai triển Taylor ta có
y(xk+1) = y(xk) + (xk+1-xk) y`(xk) + (xk+1-xk)2 y``(?k)/2 với ?k ? (xk, xk+1)
Ví dụ : Dùng công thức Euler tìm nghiệm gần đúng của bài toán Cauchy
y` = y - x2 +1, 0?x?1
y(0) = 0.5
với n = 5
Tính sai số biết nghiệm chính xác là :
y(x) = (x+1)2 - 0.5ex
giải
ta có h = 0.2
x0 = 0, x1 = 0.2, x2 = 0.4, x3 = 0.6, x4 = 0.8, x5 = 1
Công thức Euler
y0 = 0.5
yk+1 = yk + h f(xk, yk) = yk + 0.2 (yk - xk2 +1)
A = 0
B = 0.5
B = B + 0.2(B - A2 + 1) : A=A+0.2:
(A+1)2-0.5eA:Ans-B
* Nhận xét : công thức Euler đơn gian, nhưng sai số còn lớn nên ít được sử dụng
2. Công thức Euler cải tiến :
yk+1 = yk + (k1+k2)/2 k = 0,1, ..., n-1
k1 = hf(xk, yk),
k2 = hf(xk+h, yk + k1)
với h = xk+1 - xk
Ví dụ : Dùng công thức Euler cải tiến tìm nghiệm gần đúng của bài toán Cauchy
y` = y - x2 +1, 0?x?1
y(0) = 0.5
với n = 5
Tính sai số biết nghiệm chính xác là :
y(x) = (x+1)2 - 0.5ex
giải
ta có h = 0.2
x0 = 0, x1 = 0.2, x2 = 0.4, x3 = 0.6, x4 = 0.8, x5 = 1
A = 0 (xk)
B = 0.5 (yk)
C = 0.2(B - A2 + 1) :
D = 0.2(B + C - (A+0.2)2 + 1):
B=B + (C+D)/2:
A=A+0.2:
(A+1)2-0.5eA:Ans-B
3. Công thức Runge Kutta bậc 4 :
Ví dụ : Xét bài toán Cauchy
y` = 2.7xy + cos (x+2.7y), 1.2?x
y(1.2) = 5.4
Dùng công thức Runge-Kutta tính gần đúng y(1.5) với bước h = 0.3
xo = 1.2, yo = 5.4
y1 = y0 + (K1+ 2K2+ 2K3+ K4) /6
Công thức Runge-Kutta bậc 4
giải
K1= 0.3(2.7xoyo + cos(xo+2.7yo))
K2= 0.3(2.7(xo+0.3/2)(yo+K1/2) +cos(xo+0.3/2 +2.7(yo+K1/2))
K3= 0.3(2.7(xo+0.3/2)(yo+K2/2) +cos(xo+0.3/2 +2.7(yo+K2/2))
K4= 0.3(2.7(xo+0.3)(yo+K3) +cos(xo+0.3 +2.7(yo+K3)
Bấm máy ta được
K1 = 4.949578057 K2 = 8.367054617
K3 = 10.33000627 K4 = 19.41193853
y(1.5) = 15.69260639 ? 15.6926
Ví dụ : Dùng công thức Runge-Kutta tìm nghiệm gần đúng của bài toán Cauchy
y` = y - x2 +1, 0?x?1
y(0) = 0.5
với n = 5
Tính sai số biết nghiệm chính xác là :
y(x) = (x+1)2 - 0.5ex
giải
ta có h = 0.2
x0 = 0, x1 = 0.2, x2 = 0.4, x3 = 0.6, x4 = 0.8, x5 = 1
A = 0 (xk)
B = 0.5 (yk)
C = 0.2(B - A2 + 1) :
D = 0.2(B + C/2 - (A+0.1)2 + 1):
E = 0.2(B + D/2 - (A+0.1)2 + 1):
F = 0.2(B + E - (A+0.2)2 + 1):
B =B + (C+2D+2E+F)/6:
A =A+0.2:
(A+1)2-0.5eA:Ans-B
yk+1 = yk + (K1+ 2K2+ 2K3+ K4) /6
Công thức Runge-Kutta bậc 4
K2 = 0.2 [yk + 0.1(yk - xk2 +1) -(xk+0.1)2 +1 ]
= 0.2(1.1 yk - 1.1xk2 - 0.2xk + 1.09)
K1= 0.2(yk - xk2 +1)
K3 = 0.2[ yk + 0.1(1.1yk - 1.1xk2 - 0.2xk + 1.09)
- (xk+0.1)2 +1 ]
= 0.2(1.11yk - 1.11xk2 - 0.22xk + 1.099)
K4 = 0.2[ yk+0.2(1.11yk-1.11xk2-0.22xk+1.099)
- (xk+0.2)2 +1 ]
= 0.2(1.222yk-1.222xk2-0.444xk+1.1798)
y0 = 0.5
yk+1 = yk+0.2(6.642yk-6.642xk2-1.284xk+6.5578)/6
II. GIẢI GẦN ĐÚNG HỆ PTVP :
Xét hệ phương trình vi phân cấp 1
y`1 = f1(x, y1, y2, ..., ym)
y`2 = f2(x, y1, y2, ..., ym)
. . .
y`m = fm(x, y1, y2, ..., ym)
với a? x ? b và thỏa điều kiện ban đầu
y1(a) = ?1, y2(a) = ?2, .... , ym(a) = ?m
Nghiệm y = (y1, y2, ., ym)
Để tìm nghiệm gần đúng, ta chia đoạn [a,b] thành n đoạn nhỏ bằng nhau với bước h = (b-a)/n và các điểm chia
xo= a, x1 = x0 +h, ... , xk = x0 + kh, ... , xn = b
Nghiệm gần đúng là dãy { yk=(y1 k, y2 k, ., ym k)}
với yi k ? yi(xk)
Công thức Euler cải tiến :
yi k+1 = yi k + (K1 i + K2 i) / 2
K1 i = h fi(xk, y1 k, . , ym k)
K2 i = h fi(xk+h, y1 k+K1 1, . , ym k+K1 m)
?i=1,m; k = 0, n-1
Công thức Runge-Kutta bậc 4 :
yi k+1 = yi k + (K1 i+2K2 i+2K3 i+K4 i) / 6
K1 i = h fi(xk, y1 k, . , ym k)
K2 i = h fi(xk+h/2, y1 k+K11/2, . , ym k+K1 m/2)
K3 i = h fi(xk+h/2, y1 k+K21/2, . , ym k+K2 m/2)
K4 i = h fi(xk+h, y1 k+K31, . , ym k+K3 m)
?i=1,m; k = 0, n-1
Ví dụ : Sử dụng công thức Euler giải gần đúng hệ pt vi phân
y`1 = 3y1 + 2y2 - (2x2 +1)e2x
y`2 = 4y1 + y2 + (x2 +2x -4) e2x
với 0 ?x?0.5
điều kiện ban đầu y1(0)=y2(0)=1
bước h = 0.1
So sánh với nghiệm chính xác
y1(x) = 1/3e5x -1/3e-x+e2x
y2(x) = 1/3e5x +2/3e-x+x2e2x
Công thức Euler
y1 0 = 1
y1 k+1 = y1 k + h (3y1k + 2y2 k - (2xk2 +1)e2xk)
y2 0 = 1
y2 k+1 = y2 k + h (4y1k + y2 k + (xk2 +2xk -4) e2xk)
A=0 (x)
B=1 (y1k)
C=1 (y2k)
D=B + 0.1 (3B + 2C - (2A2 +1)e2A):
C=C + 0.1 (4B + C + (A2 +2A -4) e2A):
B=D:
A=A+0.1
A=0
e5A/3-e-A/3+e2A:
e5A/3+2/3e-A/3+A2e2A:
A=A+0.1
III. GIẢI GẦN ĐÚNG PTVP CẤP CAO:
Xét phương trình vi phân bậc m
y(m)(x) = f(x, y, y`, ... , y(m-1)), a?x?b
với điều kiện ban đầu
y(a) = ?1, y`(a) = ?2, .... , y(m-1)(a) = ?m
Đặt y1 = y, y2 = y`, y3 = y", ... , ym = y(m-1)
Ta chuyển phương trình vi phân bậc m về hệ m phương trình vi phân cấp 1
với điều kiện ban đầu
y1(a) = ?1, y2(a) = ?2, .... , ym(a) = ?m,
Ví dụ : Sử dụng công thức Euler giải gần đúng pt vi phân cấp 2
y " - 2 y` + 2y = sinx e2x , 0?x?0.5
điều kiện ban đầu
y(0) = -0.4, y`(0) = -0.6
với bước h = 0.1
So sánh với nghiệm chính xác biết nghiệm CX y1(x) = 0.2e2x (sinx - 2cosx)
y2(x) = 0.2e2x(4sinx - 3cosx)=y`
đặt y1 = y, y2 = y` chuyển pt về hệ
y`1 = y2
y`2 = sinx e2x- 2 y1 + 2y2
điều kiện y1(0) = -0.4, y2(0) = -0.6
Công thức Euler
y1 0 = -0.4
y1 k+1 = y1 k + 0.1 y2k
y2 0 = -0.6
y2 k+1 = y2 k + 0.1 (sinxke2xk - 2y1k +2y2 k)
A=0
B=-0.4
C=-0.6
D=B+0.1C
C=C+0.1(sinAe2A - 2B + 2C)
B=D
A=A+0.1