實驗報告七常微分方程初值問題的數值解法

1、課程名稱數值計算方法實驗項目名稱常微分方程初值問題的數值解法實驗成績 指導老師(簽名 ) 日期2015/12/16一.實驗目的和要求1 .用Matlab軟件掌握求微分方程數值解的歐拉方法和龍格-庫塔方法；2 .通過實例學習用微分方程模型解決簡化的實際問題。二.實驗內容和原理編程題2-1要求寫出Matlab源程序(m文件)，并有適當的注釋語句；分析應用題2-2, 2-3 ,2-4 , 2-5要求將問題的分析過程、Matlab源程序和運行結果和結果的解釋、算法的分析寫在實驗報告上。2-1編程編寫用向前歐拉公式和改進歐拉公式求微分方程數值解的Matlab程序，問題如下：在區(qū)間a,b內(N 1)個等距

2、點處，逼近下列初值問題的解，并對程序的每一句添上注釋語 句。y f(x, y) a x byEuler 法 y=euler(a,b,n,y0,f,f1,b1)改進 Euler 法 y=eulerpro(a,b,n,y0,f,f1,b1)2-2分析應用題假設等分區(qū)間數 n 100,用歐拉法和改進歐拉法在區(qū)間t 0,10內求解初值問題y(t) y(t) 20y(0) 10并作出解的曲線圖形，同時將方程的解析解也畫在同一張圖上，并作比較，分析這兩種方法的精度。2-3分析應用題用以下三種不同的方法求下述微分方程的數值解，取h 10y y 2x y(0) 1畫出解的圖形，與精確值比較并進行分析。1)歐拉

3、法；2)改進歐拉法；3)龍格一庫塔方法；2-4分析應用題考慮一個涉及到社會上與眾不同的人的繁衍問題模型。假設在時刻t (單位為年)，社會上有人口 x(t)人，又假設所有與眾不同的人與別的與眾不同的人結婚后所生后代也是與眾 不同的人。而固定比例為r的所有其他的后代也是與眾不同的人。如果對所有人來說出生率假定為常數b,又如果普通的人和與眾不同的人的婚配是任意的，則此問題可以用微分方程表示為：dp。 rb(1 p(t)dt其中變量p(t) Xi(t”x(t)表示在時刻t社會上與眾不同的人的比例，x(t)表示在時刻t人口中與眾不同的人的數量。1)假定p(0) 0.01,b 0.02和r 0.1,當步長

4、為h 1年時，求從t 0到t 50解p(t)的近 似值，并作出近似解的曲線圖形。2)精確求出微分方程的解p(t),并將你當t 50時在分題(b)中得到的結果與此時的精確值進行比較?！綧ATLABf關函數求微分方程的解析解及其數值的代入dsolve(egn1 ,egn2 , L x)subs (expr, x,y, x1,y1,)其中egni表示第i個方程，x表示微分方程中的自變量，默認時自變量為t。subs命令中的expr、x、y為符合型表達式，x、y分別用數值x1、x2代入。 syms x y z subs(x+y+z,x,y,z,1,2,3)ans =6 syms x subs(xA2,x

5、,2)ans =4s=dsolve( Dy 1 y 2 , y(0) 1 , x)ans =tan(x 1 4 pi) syms x subs(s,x,2)ans =右端函數f(x, y)的自動生成f= inline(expr , varl ,var2 ,)其中expr表示函數的表達式，van ,var2 表示函數表達式中的變量，運行該函數，生成一個新的函數表達式為f (val var2, )。 f=inline(x+3*y,x,y)f =Inline function:f(x,y) = x+3*y f(2,3)ans =114, 5階龍格-庫塔方法求解微分方程數值解t,x=ode45(f,t

6、s,x0,options)其中f是由待解方程寫成的m文件名；x0為函數的初值；t,x分別為輸出的自變量和函數彳1(列向量)，t的步長是程序根據誤差限自動選定的。若ts呻0,t1,t2,tf,則輸出在自變量指定值，等步長時用ts=t0:k:tf ,輸出在等分點；options用于設定誤差限(可以缺省，缺省時設定為相對誤差10 3 ,絕對誤差10 6),程序為：options=odeset( reltol ,rt, abstol ,at),這里 rt,at 分別為設定的相對誤差 和絕對誤差。常用選項見下表。選項名功能可選值省缺值AbsTol設定絕對誤差正數1e 6RelTol設定相對誤差正數1e

7、3InitialStep設定初始步長正數自動MaxStep僅7E J K-t界正數tspan/10MaxOrder設定ode15s的最高階數1,2,3,4,5 15Stats顯示計算成本統(tǒng)計on,offoffBDF設定ode15s是否用反向差分on,offoff例：解微分方程y y 2-0 t 4yy(0) 1在命令窗口執(zhí)行odefun = inline ( y 2* t/y , t , y );t, y ode45(odefun,0,4,1) ； t,yans =0plot( t , y , o- ,)%解函數圖形表示ode45(odefun,0,4,1)% 不用輸出變量，則直接輸出圖形t,

8、 y ode45(odefun,0 : 4,1) ； t,yans =0三 . 操作方法與實驗步驟(包括實驗數據記錄和處理)2-1 編程編寫用向前歐拉公式和改進歐拉公式求微分方程數值解的 Matlab 程序，問題如下：在區(qū)間 a,b 內 (N 1)個等距點處，逼近下列初值問題的解，并對程序的每一句添上注釋語句。y f (x, y) a x b y(a)y0Euler 法 y=euler(a,b,n,y0,f,f1,b1)改進 Euler 法 y=eulerpro(a,b,n,y0,f,f1,b1)Euler 法y=euler(a,b,n,y0,f,f1,b1)y=zeros(1,n+1);y(

9、1)=y0;h=(b-a)/n;x=a:h:b;for i=1:n;y(i+1)=y(i)+h*f(x(i),y(i);end plot(x,y)hold on% 求微分方程的精確解x1=linspace(a,b,100); 精確解為 s=dsolve(f1,b1,x)syms xy1=zeros(1,100);fori=1:100y1(i)=subs(s,x,x1(i);endplot(x1,y1,r)title( 紅色代表精確解)改進 Euler 法y=eulerpro(a,b,n,y0,f,f1,b1)% 求微分方程的數值解y=zeros(1,n+1);y(1)=y0;h=(b-a)/n

10、;x=a:h:b;fori=1:n;T1=f(x(i),y(i);T2=f(x(i+1),y(i)+h*T1);y(i+1)=y(i)+(h/2)*(T1+T2);endplot(x,y) hold on% 求微分方程的精確解x1=linspace(a,b,100); 精確解為 s=dsolve(f1,b1,x)syms xy1=zeros(1,100);for i=1:100y1(i)=subs(s,x,x1(i);endplot(x1,y1,r)title( 紅色代表精確解)2-2 分析應用題假設等分區(qū)間數n 100 ，用歐拉法和改進歐拉法在區(qū)間 t 0,10 內求解初值問題y (t) y

11、(t) 20 y(0) 10并作出解的曲線圖形， 同時將方程的解析解也畫在同一張圖上， 并作比較， 分析這兩種方法 的精度。(1) 向前歐拉法 euler(0,10,100,10,inline(y-20,x,y),Dy=y-20,y(0)=10)ans =精確解為s =20 - 10*exp(x)ans =+005 *Columns 1 through 8Columns 9 through 16Columns 17 through 24Columns 25 through 32Columns 33 through 40Columns 41 through 48Columns 49 throug

12、h 56Columns 57 through 64Columns 65 through 72Columns 73 through 80Columns 81 through 88Columns 89 through 96Columns 97 through 1012)改進歐拉法 eulerpro(0,10,100,10,inline(y-20,x,y),Dy=y-20,y(0)=10) ans =精確解為s =20 - 10*exp(x)ans =+005 *Columns 1 through 8Columns 9 through 16Columns 17 through 24Columns 2

13、5 through 32Columns 33 through 40Columns 41 through 48Columns 49 through 56Columns 57 through 64Columns 65 through 72Columns 73 through 80Columns 81 through 88Columns 89 through 96Columns 97 through 101改進歐拉法的精度比向前歐拉法更高。2-3分析應用題用以下三種不同的方法求下述微分方程的數值解，取 h 10y y 2x 0x1 y(0) 1畫出解的圖形，與精確值比較并進行分析。1)歐拉法；2)改

14、進歐拉法； 2-4分析應用題考慮一個涉及到社會上與眾不同的人的繁衍問題模型。假設在時刻t (單位為年)，社會上有人口 x(t)人，又假設所有與眾不同的人與別的與眾不同的人結婚后所生后代也是與眾 不同的人。而固定比例為r的所有其他的后代也是與眾不同的人。如果對所有人來說出生率假定為常數b,又如果普通的人和與眾不同的人的婚配是任意的，則此問題可以用微分方程表示為：dp(t) rb(1 P(t) dt其中變量p(t) x(t)/x(t)表示在時刻t社會上與眾不同的人的比例，x(t)表示在時刻t人口中與眾不同的人的數量。1)假定p(0) 0.01,b 0.02和r 0.1,當步長為h 1年時，求從t 0到t 50解p(t)的近 似值，并作出近似解的曲線圖形。2)精確求出微分方程的解p(t),并將你當t 50時在分題(b)中得到的結果與此時的精確值進行比較。1)euler(0,50,50,inline(,t,p),Dp=,p(0)= 1)ans =精確解為1 - 99/(100*exp(x/500)ans =Columns 1 thro