差分方程及微分方程數(shù)值解

1、東南大學(xué)數(shù)學(xué)實(shí)驗(yàn)報(bào)告實(shí)驗(yàn)內(nèi)容：差分方程及微分方程數(shù)值解一 實(shí)驗(yàn)?zāi)康?熟悉迭代法及微分方程數(shù)值方法二 預(yù)備知識(shí) （1）了解差分方程穩(wěn)定性、周期分解、混沌等相關(guān)知識(shí) （2）了解歐拉方法、龍格-庫(kù)特方法。三 實(shí)驗(yàn)內(nèi)容與要求 （一）Volterra方程數(shù)值解方程其中a=1,b=0.1,c=0.5,d=0.04命令與結(jié)果在函數(shù)編輯器中輸入：function dxdt = euler( t,x )dxdt= x(1)*(1-0.1*x(2) x(2)*(-0.5+0.02*x(1);end四階龍格-庫(kù)塔公式：在命令窗口中輸入：tspan=0 15;x0=25;2;t,x=ode45(euler,tspan

2、,x0);plot(t,x(:,1),'r-','LineWidth',0.5);hold on;plot(t,x(:,2),'g-','LineWidth',0.5);hold on;（1） 相圖四階龍格-庫(kù)塔公式：1 / 9axis(0 15 0 125)legend('x(1)','x(2)')grid ontitle('The Numerical Solution Of Calculable way of fourth rank Rounge-kutt')plot(x(:,1

3、),x(:,2)歐拉方法：在命令窗口中輸入：tspan=0 15;x0=25;2;t,x=ode23(euler,tspan,x0);plot(t,x(:,1),'r-','LineWidth',0.5);hold on;plot(t,x(:,2),'g-','LineWidth',0.5);hold on;axis(0 15 0 125)legend('x(1)','x(2)')grid ontitle('The Numerical Solution Of Euler Equation&#

4、39;)>> plot(x(:,1),x(:,2)歐拉公式：（2） 軌線(xiàn)圖四/五階龍格-庫(kù)塔公式：歐拉公式：（2） 差分阻滯增長(zhǎng)模型在t時(shí)刻單位時(shí)間內(nèi)的人口數(shù)量的變化量?jī)H僅與此時(shí)的人口數(shù)量x有關(guān)(等于右邊的值),其中的r表示人口的固有增長(zhǎng)率,N表示能容納的最大人口數(shù).經(jīng)過(guò) 單位時(shí)間，即一個(gè)繁殖周期的種群數(shù)量的增長(zhǎng)量?jī)H僅與前一個(gè)時(shí)期的種群數(shù)量有關(guān)。于是模型為 即這是一個(gè)一階非線(xiàn)性差分方程.對(duì)于給定的初值,我們可以從這個(gè)遞推公式運(yùn)用計(jì)算機(jī)很容易地計(jì)算出一些xk.這是在計(jì)算機(jī)出現(xiàn)以后的一個(gè)新的特點(diǎn).但是我們更關(guān)心的是當(dāng)時(shí)間趨于無(wú)窮時(shí),即k趨于無(wú)窮時(shí), xk的極限如何，即差分方程平衡點(diǎn)的穩(wěn)

5、定性問(wèn)題。求得差分方程(2)的平衡點(diǎn)為0和(b-1)/b.它們分別對(duì)應(yīng)于差分方程(1)的平衡點(diǎn)0和N.我們將這個(gè)差分方程(2)在平衡點(diǎn)附近展開(kāi),有注意到b=1+r>1,平衡點(diǎn)0是不穩(wěn)定的.略去高階項(xiàng)得因此當(dāng)|2-b|<1時(shí),平衡點(diǎn)x* (b-1)/b是穩(wěn)定的. |2-b|>1即 b>3時(shí),平衡點(diǎn)x*是不穩(wěn)定的.（1）穩(wěn)定情況（） 穩(wěn)定性分析b=1.6,在1到3范圍之間穩(wěn)定命令與結(jié)果（圖）function x=block(b,x0)b=1.6;if nargin<2 x0=0.5;endx=x0;n=1:10000;for i=2:10000 x(i)=b*x(i-

6、1)*(1-x(i-1);endplot(n,x,'.');end結(jié)果圖：穩(wěn)定性分析b=2.2,在1到3范圍之間穩(wěn)定命令與結(jié)果（圖）function x=block(b,x0)b=2.2;if nargin<2 x0=0.5;endx=x0;n=1:10000;for i=2:10000 x(i)=b*x(i-1)*(1-x(i-1);endplot(n,x,'.');end結(jié)果圖：（2）不穩(wěn)定情況穩(wěn)定性分析b=3.3在1到3范圍之間不穩(wěn)定命令與結(jié)果（圖）穩(wěn)定性分析b=3.63在1到3范圍之間不穩(wěn)定命令與結(jié)果（圖）（3）倍周期分解 （1） 求穩(wěn)定平衡點(diǎn)可得

7、兩個(gè)平衡點(diǎn)：0.4794 0.8236程序：function x=block(b,x0)b=3.3;if nargin<2 x0=0.5;endx=x0;n=1:100;for i=2:100 x(i)=b*x(i-1)*(1-x(i-1);end（2） 數(shù)值結(jié)果ans = Columns 1 through 160.5000 0.8250 0.4764 0.8232 0.4804 0.8237 0.4792 0.8236 0.4795 0.8236 0.4794 0.8236 0.4794 0.8236 0.4794 0.8236（1）求穩(wěn)定平衡點(diǎn)四個(gè)平衡點(diǎn)：0.5009 0.8750

8、 0.3828 0.8269程序：function x=block(b,x0)b=3.5;if nargin<2 x0=0.5;endx=x0;n=1:100;for i=2:100 x(i)=b*x(i-1)*(1-x(i-1);end（2） 數(shù)值結(jié)果ans = Columns 1 through 160.5000 0.8750 0.3828 0.8269 0.5009 0.8750 0.3828 0.8269 0.5009 0.8750 0.3828 0.8269 0.5009 0.8750 0.3828 0.8269 Columns 17 through 320.4794 0.8236 0.4794 0.8236 0.4794 0.8236 0.4794 0.8236 0.4794 0.8236 0.4794 0.8236 0.4794 0.8236 0.4794 0.8236 Columns 17 through 320.5009 0.8750 0.3828 0.8269 0.