Xem mẫu

  1. ĐẠ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í
  2. Đề 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.
  3. 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
  4. 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
  5. 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ó:
  6. 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:
  7. 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
  8. 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
  9. 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