第7講-微分方程省公開課金獎(jiǎng)全國(guó)賽課一等獎(jiǎng)微課獲獎(jiǎng)?wù)n件

數(shù)學(xué)建模與數(shù)學(xué)試驗(yàn)微分方程1/19微分方程模型建立微分方程模型要對(duì)研究對(duì)象作詳細(xì)分析，普通有以下三種方法：一是依據(jù)規(guī)律建模；二是用微元法建模；三是用模擬近似法建模。一依據(jù)規(guī)律建模類模型：如后面導(dǎo)彈跟蹤問題；依據(jù)數(shù)學(xué)、力學(xué)、物理、化學(xué)等學(xué)科中已經(jīng)有規(guī)律和定律來進(jìn)行建模；二微元法建模類模型：容器漏水問題（試驗(yàn)P76）此方法也是利用已經(jīng)有規(guī)律和定律來尋求一些微元之間關(guān)系式。三模擬近似法建模：社會(huì)科學(xué)、生物學(xué)、醫(yī)學(xué)、經(jīng)濟(jì)學(xué)等學(xué)科，普通都是用此方法2/19微分方程模型如《數(shù)?！罚簜魅静∧Ｐ?、經(jīng)濟(jì)增加模型、正規(guī)戰(zhàn)與游擊戰(zhàn)模型等等。求解方法：歐拉法（即用差商代替導(dǎo)數(shù)）、改進(jìn)歐拉法（即數(shù)值積分）、使用泰勒公式法動(dòng)態(tài)連續(xù)模型用微分方程模型來建立。下面詳細(xì)介紹3/19（二）建立數(shù)值解法一些路徑1、用差商代替導(dǎo)數(shù)若步長(zhǎng)h較小，則有故有公式：此即歐拉法。4/192、使用數(shù)值積分對(duì)方程y’=f(x,y),兩邊由xi到xi+1積分，并利用梯形公式，有：實(shí)際應(yīng)用時(shí)，與歐拉公式結(jié)合使用：此即改進(jìn)歐拉法。故有公式：5/193、使用泰勒公式以此方法為基礎(chǔ)，有龍格-庫(kù)塔法、線性多步法等方法。4、數(shù)值公式精度當(dāng)一個(gè)數(shù)值公式截?cái)嗾`差可表示為O（hk+1）時(shí)（k為正整數(shù)，h為步長(zhǎng)），稱它是一個(gè)k階公式。k越大，則數(shù)值公式精度越高。歐拉法是一階公式，改進(jìn)歐拉法是二階公式。龍格-庫(kù)塔法有二階公式和四階公式。線性多步法有四階阿達(dá)姆斯外插公式和內(nèi)插公式。返回6/19（三）用Matlab軟件求常微分方程數(shù)值解[t，x]=solver（’f’,ts,x0,options）ode45ode23ode113ode15sode23s由待解方程寫成m-文件名ts=[t0，tf]，t0、tf為自變量初值和終值函數(shù)初值ode23：組合2/3階龍格-庫(kù)塔-芬爾格算法ode45：利用組合4/5階龍格-庫(kù)塔-芬爾格算法自變量值函數(shù)值用于設(shè)定誤差限(缺省時(shí)設(shè)定相對(duì)誤差10-3,絕對(duì)誤差10-6),命令為：options=odeset（’reltol’,rt,’abstol’,at）,rt，at：分別為設(shè)定相對(duì)誤差和絕對(duì)誤差.7/191、在解n個(gè)未知函數(shù)方程組時(shí)，x0和x均為n維向量，m-文件中待解方程組應(yīng)以x分量形式寫成.2、使用Matlab軟件求數(shù)值解時(shí)，高階微分方程必須等價(jià)地變換成一階微分方程組.注意:8/19導(dǎo)彈追蹤問題設(shè)位于坐標(biāo)原點(diǎn)甲艦向位于x軸上點(diǎn)A(1,0)處乙艦發(fā)射導(dǎo)彈，導(dǎo)彈頭一直對(duì)準(zhǔn)乙艦.假如乙艦以最大速度v0(是常數(shù))沿平行于y軸直線行駛，導(dǎo)彈速度是5v0，求導(dǎo)彈運(yùn)行曲線方程.又乙艦行駛多遠(yuǎn)時(shí)，導(dǎo)彈將它擊中？解法一（解析法）9/19由(1),(2)消去t整理得模型:ToMatlab(chase1)軌跡圖見程序chase110/19解法二(數(shù)值解)1.建立m-文件eq1.m

functiondy=eq1(x,y)dy=zeros(2,1);dy(1)=y(2);dy(2)=1/5*sqrt(1+y(1)^2)/(1-x);2.取x0=0，xf=0.9999，建立主程序ff6.m以下：

x0=0，xf=0.9999[x,y]=ode15s('eq1',[x0xf],[00]);plot(x,y(:,1),’b.')holdony=0:0.01:2;plot(1,y,’b*')

結(jié)論:導(dǎo)彈大致在（1，0.2）處擊中乙艦ToMatlab(ff6)令y1=y,y2=y1’，將方程（3）化為一階微分方程組。11/19慢跑者與狗一個(gè)慢跑者在平面上沿橢圓以恒定速率v=1跑步,設(shè)橢圓方程為:x=10+20cost,y=20+5sint.突然有一只狗攻擊他.這只狗從原點(diǎn)出發(fā),以恒定速率w跑向慢跑者,狗運(yùn)動(dòng)方向一直指向慢跑者.分別求出w=20,w=5時(shí)狗運(yùn)動(dòng)軌跡.1.模型建立設(shè)時(shí)刻t慢跑者坐標(biāo)為(X(t)，Y(t))，狗坐標(biāo)為(x(t)，y(t)).則X=10+20cost,Y=20+15sint,狗從(0,0)出發(fā),與導(dǎo)彈追蹤問題類似，建立狗運(yùn)動(dòng)軌跡參數(shù)方程:12/192.模型求解(1)w=20時(shí),建立m-文件eq3.m以下:functiondy=eq3(t,y)dy=zeros(2,1);dy(1)=20*(10+20*cos(t)-y(1))/sqrt((10+20*cos(t)-y(1))^2+(20+15*sin(t)-y(2))^2);dy(2)=20*(20+15*sin(t)-y(2))/sqrt((10+20*cos(t)-y(1))^2+(20+15*sin(t)-y(2))^2);取t0=0，tf=10，建立主程序chase3.m以下：t0=0;tf=10;[t,y]=ode45('eq3',[t0tf],[00]);T=0:0.1:2*pi;X=10+20*cos(T);Y=20+15*sin(T);plot(X,Y,'-')holdonplot(y(:,1),y(:,2),'*')在chase3.m，不停修改tf值,分別取tf=5,2.5,3.5,…,至3.15時(shí),狗剛好追上慢跑者.ToMatlab(chase3)13/19建立m-文件eq4.m以下:functiondy=eq4(t,y)dy=zeros(2,1);dy(1)=5*(10+20*cos(t)-y(1))/sqrt((10+20*cos(t)-y(1))^2+(20+15*sin(t)-y(2))^2);dy(2)=5*(20+15*sin(t)-y(2))/sqrt((10+20*cos(t)-y(1))^2+(20+15*sin(t)-y(2))^2);取t0=0，tf=10，建立主程序chase4.m以下：t0=0;tf=10;[t,y]=ode45('eq4',[t0tf],[00]);T=0:0.1:2*pi;X=10+20*cos(T);Y=20+15*sin(T);plot(X,Y,'-')holdonplot(y(:,1),y(:,2),'*')在chase3.m，不停修改tf值,分別取tf=20,40,80,…,能夠看出,狗永遠(yuǎn)追不上慢跑者.ToMatlab(chase4)(2)w=5時(shí)返回14/19地中海鯊魚問題意大利生物學(xué)家Ancona曾致力于魚類種群相互制約關(guān)系研究，他從第一次世界大戰(zhàn)期間,地中海各港口捕捉幾個(gè)魚類捕捉量百分比資料中，發(fā)覺鯊魚等百分比有顯著增加（見下表），而供其捕食食用魚百分比卻顯著下降.顯然戰(zhàn)爭(zhēng)使打魚量下降，食用魚增加，鯊魚等也隨之增加，但為何鯊魚百分比大幅增加呢？他無法解釋這個(gè)現(xiàn)象，于是求援于著名意大利數(shù)學(xué)家V.Volterra，希望建立一個(gè)食餌—捕食系統(tǒng)數(shù)學(xué)模型，定量地回答這個(gè)問題.15/19

該模型反應(yīng)了在沒有些人工捕捉自然環(huán)境中食餌與捕食者之間制約關(guān)系，沒有考慮食餌和捕食者本身阻滯作用，是Volterra提出最簡(jiǎn)單模型.16/19首先，建立m-文件shier.m以下：functiondx=shier(t,x)dx=zeros(2,1);dx(1)=x(1)*(1-0.1*x(2));dx(2)=x(2)*(-0.5+0.02*x(1));其次，建立主程序shark.m以下：[t,x]=ode45('shier',[015],[252]);plot(t,x(:,1),'-',t,x(:,2),'*')plot(x(:,1),x(:,2))ToMatlab(shark)17/19求解結(jié)果：