Xem mẫu
- ĐẠI HỌC QUỐC GIA THÀNH PHỐ HỒ CHÍ MINH
TRƯỜNG ĐẠI HỌC BÁCH KHOA
***
KHOA ĐIỆN – ĐIỆN TỬ
Giáo viên hướng dẫn: Lê Thị Quỳnh Hà
Các thành viên trong nhóm:
40901457: Nguyễn Phước Lộc
40902767: Võ Nhựt Tiến
40902741: Cù Văn Tiến
40903057: Huỳnh Trung Trực
40902907: Nguyễn Kim Triển
40902935: Phan Đức Trí
- Đề tài: Viết chương trình giải phương trình bằng phương pháp Runge-Kutta bậc 4.
Vẽ đồ thị hàm nhận được.
- CÔ SÔÛ LYÙ THUYEÁT
BAØI TOAÙN CAUCHY
M t ph ng trình vi phân c p 1 có th vi t d i d ng gi i c
yf (x,y) mà ta có th tìm c hàm y t o hàm c a nó. T n t i vô s
nghi m tho mãn ph ng trình trên. M i nghi m ph thu c vào m t h ng s
tu ý. Khi cho tr c giá tr ban u c a y là yo t i giá tr u xo ta nh n c
m t nghi m riêng c a ph ng trình. Bài toán Cauchy (hay bài toán có i u
ki n u) tóm l i nh sau: cho x sao cho b x a, tìm y(x) tho mãn i u
ki n:
y ( x) f ( x , y )
(1)
y(a )
Ng i ta ch ng minh r ng bài toán này có m t nghi m duy nh t n u f
tho mãn i u ki n Lipschitz:
f( x , y 1 ) f( x , y 2 ) L y 1 y 2
v i L là m t h ng s d ng.
Ng i ta c ng ch ng minh r ng n u f y ( o hàm c a f theo y ) là liên
t c và b ch n thì f tho mãn i u ki n Lipschitz.
M t cách t ng quát h n, ng i ta nh ngh a h ph ng trình b c 1:
y 1 f1 ( x , y 1 , y 2 ,..., y n )
y2 f2 ( x , y 1 , y 2 ,..., y n )
yn fn ( x , y 1 , y 2 ,..., y n )
Ta ph i tìm nghi m y1, y2,..., yn sao cho:
Y ( x) f ( x , X )
Y(a )
v i:
y1 f1 y1
y2 f2 y2
F Y
Y .. .. ..
.. .. ..
fn yn
yn
- N u ph ng trình vi phân có b c cao h n (n), nghi m s ph thu c vào
n h ng s tu ý. nh n c m t nghi m riêng, ta ph i cho n i u ki n u.
Bài toán s có giá tr u n u v i giá tr xo ã cho ta cho y(xo), y (xo), y (xo),....
M t ph ng trình vi phân b c n có th a v thành m t h ph ng
trình vi phân c p 1. Ví d n u ta có ph ng trình vi phân c p 2:
y f( x , y , y )
y(a ) , y (a )
t u = y và v = y ta nh n c h ph ng trình vi phân c p 1:
Khi
uv
v g( x , u , v )
v i i u ki n u: u(a) = và v(a) =
Các ph ng pháp gi i ph ng trình vi phân c trình bày trong
ch ng này là các ph ng pháp r i r c: o n [a, b] c chia thành n o n
nh b ng nhau c g i là các b c tích phân h = ( b a) / n.
§2. PH NG PHÁP EULER
- PHÖÔNG PHAÙP RUNGE KUTTA
ch a chính xác i v i các bài toán th c t .
Xét bài toán Cauchy (1). Gi s ta ã tìm c giá tr g n úng yi c a
y(xi) và mu n tính yi+1 c a y(xi+1). Tr c h t ta vi t công th c Taylor:
h2 h m (m) h m 1 ( m 1)
y( x i 1 ) y( x i ) hy ( x i ) y (x i ) y (x i ) y (c) (11)
2 m! m!
v i c (xi, xi+1) và:
y ( x i ) f x i , y( x i )
dk 1
( k)
y (xi ) f x i , y(x i )
dx k 1
Ta vi t l i (11) d i d ng:
h2 h m (m) h m 1 ( m 1)
y i 1 y i hy ( x i ) y (x i ) y (x i ) y ( c) (12)
2 m! m!
Ta ã kéo dài khai tri n Taylor k t qu chính xác h n. tính y i, y i v.v. ta
có th dùng ph ng pháp Runge Kutta b ng cách t:
y i 1 y i r1 k (1i ) r2 k (2i ) r3 k (3i ) r4 k (4i ) (13)
trong ó:
k (1i ) hf( x i , y i )
k (2i ) k (1i ) )
hf( x i ah , y i
(14)
k (3i ) k (1i ) k (2i ) )
hf( x i bh , y i
.......
và ta c n xác nh các h s a, b,..; , , ,...; r1, r2,.. sao cho v ph i c a (13)
khác v i v ph i c a (12) m t vô cùng bé c p cao nh t có th có i v i h.
Khi dùng công th c Runge Kutta b c hai ta có:
k (1i ) hf( x i , y i )
(15)
k (2i ) hf( x i ah , y i k (1i ) )
y i 1 y i r1 k 1i ) r2 k (2i )
(
và (16)
Ta có:
y (x) = f[x,y(x)]
y ( x) fx x , y( x) fy x , y( x)
................
Do ó v ph i c a (12) là:
h2
hf( x i , y i ) fx (x i , y i ) fy ( x i , y i ) y ( x) (17)
2
M t khác theo (15) và theo công th c Taylor ta có:
- k 1i )
(
hf ( x i , y i ) hy i
k (2i ) k (1i ) fy ( x i , y i )
h[f( x i , y i ) ahfx ( x i , y i ) ]
Do ó v ph i c a (16) là:
h(r1 r2 )f( x i , y i ) h 2 [ar2 fx ( x i , y i ) r2 y i fy ( x i , y i )] (18)
Bây gi cho (17) và (18) khác nhau m t vô cùng bé c p O(h3) ta tìm c các
h s ch a bi t khi cân b ng các s h ng ch a h và ch a h2:
r1 + r2 = 1
a.r1 = 1/ 2
.r2 = 1
Nh v y: = a, r1 = (2a 1)/ 2a, r2 = 1/ 2a v i a c ch n b t kì.
N u a = 1 / 2 thì r1 = 0 và r2 = 1. Lúc này ta nh n c công th c Euler. N u
a=1 thì r1 = 1 / 2 và r2 = 1/2. Lúc này ta nh n c công th c Euler c i ti n.
M t cách t ng t chúng ta nh n c công th c Runge Kutta b c 4.
Công th c này hay c dùng trong tính toán th c t :
k1 = h.f(xi, yi)
k2 = h.f(xi+h/ 2, yi + k1/ 2)
k3 = h.f(xi+h/ 2, yi + k2/ 2)
k4 = h.f(xi+h, yi + k3)
yi+1 = yi + (k1 + 2k2 + 2k3 + k4) / 6
Ta xây d ng hàm rungekutta() th c hi n công th c Runge Kutta b c 4:
- 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
- Xây dựng hàm rk4 trong matlab để giải phương trình vi phân theo phương pháp trên.
function[x,y]=rk4(f,x0,x1,y0,h)
if nargin [X,Y]=rk4(inline('y-x^2+1','x','y'),0,1,0.5)
Định nghĩa hàm f trong của sổ Script rồi lưu thành file f.m
function [dy] = f(a,b)
dy = b-a^2 +1;
end
>> [X,Y]=rk4(@f,0,1,0.5)
Vẽ đồ thị từ các giá trị xk,yk (giá trị trả về của hàm rk4):
plot(X,Y,'rd-', 'LineWidth',3)
xlabel('Truc X')
ylabel('Truc Y')
title('DO THI CUA HAM DA CHO')
grid on
- Kết quả cửa sổ Command Windows:
>> [X,Y]=rk4(inline('y-x^2+1','x','y'),0,1,0.5)
Nhap so doan chia n = 5
n=
5
X=
0 0.2000 0.4000 0.6000 0.8000 1.0000
Y=
0.5000 0.8293 1.2141 1.6489 2.1272 2.6408
>> ve
>>
nguon tai.lieu . vn