數(shù)學(xué)建模-常微分方程模型及差分模型

微分方程模型洛陽理工學(xué)院數(shù)理部1精選課件微分方程模型人口增長(zhǎng)的預(yù)測(cè)傳染病模型種群模型2精選課件動(dòng)態(tài)模型

描述對(duì)象特征隨時(shí)間(空間)的演變過程.

分析對(duì)象特征的變化規(guī)律.

預(yù)報(bào)對(duì)象特征的未來性態(tài).

研究控制對(duì)象特征的手段.

根據(jù)函數(shù)及其變化率之間的關(guān)系確定函數(shù).微分方程建模

根據(jù)建模目的和問題分析作出簡(jiǎn)化假設(shè).

按照內(nèi)在規(guī)律或用類比法建立微分方程.3精選課件對(duì)微分方程的研究方法解在很廣泛的條件下存在，但能用有限解析式表達(dá)者很少.另辟它徑：1、求數(shù)值解（近似解）；2、定性方法分析.4精選課件背景

年1625183019301960197419871999人口(億)5102030405060世界人口增長(zhǎng)概況中國(guó)人口增長(zhǎng)概況

年19011929195319651982199019952000人口(億)4.265.486.027.2510.3211.3012.0012.95研究人口變化規(guī)律控制人口過快增長(zhǎng)人口增長(zhǎng)的預(yù)測(cè)5精選課件指數(shù)增長(zhǎng)模型——馬爾薩斯提出(1798)常用的計(jì)算公式x(t)~時(shí)刻t的人口基本假設(shè)

:人口(相對(duì))增長(zhǎng)率r

是常數(shù)今年人口x0,年增長(zhǎng)率rk年后人口隨著時(shí)間增加，人口按指數(shù)規(guī)律無限增長(zhǎng)6精選課件年17901800181018201830184018501860人口3.95.37.29.612.917.123.231.4年18701880189019001910192019301940人口38.650.262.976.092.0106.5123.2131.7年195019601970198019902000人口150.7179.3204.0226.5251.4281.4美國(guó)人口統(tǒng)計(jì)數(shù)據(jù)2012-3-177精選課件指數(shù)增長(zhǎng)模型的應(yīng)用及局限性

與19世紀(jì)以前歐洲一些地區(qū)人口統(tǒng)計(jì)數(shù)據(jù)吻合

適用于19世紀(jì)后遷往加拿大的歐洲移民后代

可用于短期人口增長(zhǎng)預(yù)測(cè)

不符合19世紀(jì)后多數(shù)地區(qū)人口增長(zhǎng)規(guī)律

不能預(yù)測(cè)較長(zhǎng)期的人口增長(zhǎng)過程19世紀(jì)后人口數(shù)據(jù)人口增長(zhǎng)率r不是常數(shù)(逐漸下降)2012-3-178精選課件阻滯增長(zhǎng)模型(Logistic模型)人口增長(zhǎng)到一定數(shù)量后，增長(zhǎng)率下降的原因：資源、環(huán)境等因素對(duì)人口增長(zhǎng)的阻滯作用且阻滯作用隨人口數(shù)量增加而變大假設(shè)r~固有增長(zhǎng)率(x很小時(shí))xm~人口容量（資源、環(huán)境能容納的最大數(shù)量）r是x的減函數(shù)2012-3-179精選課件dx/dtx0xmxm/2xmtx0x(t)~S形曲線,x增加先快后慢x0xm/2阻滯增長(zhǎng)模型(Logistic模型)2012-3-1710精選課件參數(shù)估計(jì)用指數(shù)增長(zhǎng)模型或阻滯增長(zhǎng)模型作人口預(yù)報(bào)，必須先估計(jì)模型參數(shù)r或r,xm

利用統(tǒng)計(jì)數(shù)據(jù)用最小二乘法作擬合例：美國(guó)人口數(shù)據(jù)（單位~百萬）186018701880……196019701980199031.438.650.2……179.3204.0226.5251.4阻滯增長(zhǎng)模型(Logistic模型)r=0.2557,xm=392.12012-3-1711精選課件2012-3-1712精選課件模型檢驗(yàn)用模型計(jì)算2000年美國(guó)人口，與實(shí)際數(shù)據(jù)比較實(shí)際為281.4(百萬)模型應(yīng)用——預(yù)報(bào)美國(guó)2010年的人口加入2000年人口數(shù)據(jù)后重新估計(jì)模型參數(shù)Logistic模型在經(jīng)濟(jì)領(lǐng)域中的應(yīng)用(如耐用消費(fèi)品的售量)阻滯增長(zhǎng)模型(Logistic模型)r=0.2490,xm=434.0x(2010)=306.02012-3-1713精選課件傳染病模型問題

描述傳染病的傳播過程

分析受感染人數(shù)的變化規(guī)律

預(yù)報(bào)傳染病高潮到來的時(shí)刻

預(yù)防傳染病蔓延的手段

按照傳播過程的一般規(guī)律，用機(jī)理分析方法建立模型2012-3-1714精選課件

已感染人數(shù)(病人)i(t)

每個(gè)病人每天有效接觸(足以使人致病)人數(shù)為模型1假設(shè)若有效接觸的是病人，則不能使病人數(shù)增加必須區(qū)分已感染者(病人)和未感染者(健康人)建模?2012-3-1715精選課件模型2區(qū)分已感染者(病人)和未感染者(健康人)假設(shè)1）總?cè)藬?shù)N不變，病人和健康人的比例分別為2）每個(gè)病人每天有效接觸人數(shù)為,且使接觸的健康人致病建模~日接觸率SI模型2012-3-1716精選課件模型21/2tmii010ttm~傳染病高潮到來時(shí)刻(日接觸率)tmLogistic模型病人可以治愈！?t=tm,di/dt最大2012-3-1717精選課件模型3傳染病無免疫性——病人治愈成為健康人，健康人可再次被感染增加假設(shè)SIS模型3）病人每天治愈的比例為~日治愈率建模~日接觸率1/~平均感染期~一個(gè)感染期內(nèi)每個(gè)病人的有效接觸人數(shù)，稱為接觸數(shù)。2012-3-1718精選課件模型3i0i0接觸數(shù)=1~閾值感染期內(nèi)有效接觸感染的健康者人數(shù)不超過病人數(shù)1-1/i0模型2(SI模型)如何看作模型3(SIS模型)的特例idi/dt01>10ti>11-1/i0t1di/dt<02012-3-1719精選課件模型4傳染病有免疫性——病人治愈后即移出感染系統(tǒng)，稱移出者SIR模型假設(shè)1）總?cè)藬?shù)N不變，病人、健康人和移出者的比例分別為2）病人的日接觸率

,日治愈率,

接觸數(shù)=/建模需建立的兩個(gè)方程2012-3-1720精選課件模型4SIR模型無法求出的解析解在相平面上研究解的性質(zhì)2012-3-1721精選課件模型4消去dtSIR模型相軌線的定義域相軌線11si0D在D內(nèi)作相軌線的圖形，進(jìn)行分析2012-3-1722精選課件si101D模型4SIR模型相軌線及其分析傳染病蔓延傳染病不蔓延s(t)單調(diào)減相軌線的方向P1s0imP1:s0>1/i(t)先升后降至0P2:s0<1/

i(t)單調(diào)降至01/~閾值P3P4P2S02012-3-1723精選課件模型4SIR模型預(yù)防傳染病蔓延的手段(日接觸率)衛(wèi)生水平(日治愈率)醫(yī)療水平傳染病不蔓延的條件——s0<1/的估計(jì)

降低s0提高r0

提高閾值1/

降低(=/),群體免疫2012-3-1724精選課件模型4SIR模型被傳染人數(shù)的估計(jì)記被傳染人數(shù)比例x<<s0i0P1i00,s01小,s01提高閾值1/降低被傳染人數(shù)比例xs0-1/=2012-3-1725精選課件

種群的弱肉強(qiáng)食(食餌-捕食者模型)

種群甲靠豐富的天然資源生存，種群乙靠捕食甲為生，形成食餌-捕食者系統(tǒng)，如食用魚和鯊魚，美洲兔和山貓，害蟲和益蟲。

模型的歷史背景——一次世界大戰(zhàn)期間地中海漁業(yè)的捕撈量下降(食用魚和鯊魚同時(shí)捕撈)，但是其中鯊魚的比例卻增加，為什么？精選課件食餌（甲）數(shù)量x(t),

捕食者（乙）數(shù)量

y(t)甲獨(dú)立生存的增長(zhǎng)率r乙使甲的增長(zhǎng)率減小，減小量與

y成正比乙獨(dú)立生存的死亡率d甲使乙的死亡率減小，減小量與x成正比方程(1),(2)無解析解食餌-捕食者模型(Volterra)a~捕食者掠取食餌能力b~食餌供養(yǎng)捕食者能力精選課件Volterra模型的平衡點(diǎn)及其穩(wěn)定性平衡點(diǎn)穩(wěn)定性分析P點(diǎn)穩(wěn)定性不能用近似線性方程分析p=0,q>0P:臨界狀態(tài)q<0P′不穩(wěn)定精選課件tx(t)y(t)020.00004.00000.100021.24063.96510.200022.56493.94050.300023.97633.9269………5.10009.616216.72355.20009.017316.2064………9.500018.47504.04479.600019.61363.99689.700020.83113.9587用數(shù)學(xué)軟件MATLAB求微分方程數(shù)值解x~y平面上的相軌線精選課件計(jì)算結(jié)果（數(shù)值，圖形）x(t),y(t)是周期函數(shù)，相圖(x,y)是封閉曲線觀察，猜測(cè)x(t),y(t)的周期約為9.6xmax65.5,xmin6,ymax20.5,ymin3.9用數(shù)值積分可算出

x(t),y(t)一周期的平均值：x(t)的平均值約為25,

y(t)的平均值約為10。食餌-捕食者模型(Volterra)精選課件

消去dt用相軌線分析點(diǎn)穩(wěn)定性c由初始條件確定取指數(shù)精選課件x0fmf(x)x0g(y)gmy0y0在相平面上討論相軌線的圖形用相軌線分析點(diǎn)穩(wěn)定性相軌線時(shí)無相軌線以下設(shè)精選課件y2y1xQ3Q4qy1y2x1x2pyy0xx0P0x1x2Q1Q2Q1(x1,y0),Q2(x2,y0)Q3(x,y1),Q4(x,y2)相軌線退化為P點(diǎn)

存在x1<x0<x2,使f(x1)=f(x2)=p存在y1<y0<y2,使g(y1)=g(y2)=q相軌線是封閉曲線族xQ3Q4f(x)xx0fm0g(y)gmy0y0相軌線P~中心精選課件相軌線是封閉曲線x(t),y(t)是周期函數(shù)(周期記T)求x(t),y(t)在一周期的平均值軌線中心用相軌線分析點(diǎn)穩(wěn)定性精選課件?T2T3T4T1PT1

T2

T3

T4x(t)的“相位”領(lǐng)先y(t)模型解釋初值相軌線的方向精選課件模型解釋r~食餌增長(zhǎng)率d~捕食者死亡率b~食餌供養(yǎng)捕食者能力捕食者數(shù)量食餌數(shù)量Pr/ad/ba~捕食者掠取食餌能力捕食者數(shù)量與r成正比,與a成反比食餌數(shù)量與d成正比,與b成反比精選課件模型解釋一次大戰(zhàn)期間地中海漁業(yè)的捕撈量下降，但是其中鯊魚的比例卻在增加，為什么？rr-1,dd+1捕撈戰(zhàn)時(shí)捕撈rr-2,dd+2,2<1???xy食餌(魚)減少，捕食者(鯊魚)增加自然環(huán)境

還表明：對(duì)害蟲(食餌)—益蟲(捕食者)系統(tǒng)，使用滅兩種蟲的殺蟲劑,會(huì)使害蟲增加，益蟲減少。精選課件食餌-捕食者模型(Volterra)的缺點(diǎn)與改進(jìn)Volterra模型改寫多數(shù)食餌—捕食者系統(tǒng)觀察不到周期震蕩,而是趨向某個(gè)平衡狀態(tài),即存在穩(wěn)定平衡點(diǎn)加Logistic項(xiàng)有穩(wěn)定平衡點(diǎn)精選課件

相軌線是封閉曲線，結(jié)構(gòu)不穩(wěn)定——一旦離開某一條閉軌線，就進(jìn)入另一條閉軌線，不恢復(fù)原狀。

自然界存在的周期性平衡生態(tài)系統(tǒng)是結(jié)構(gòu)穩(wěn)定的，即偏離周期軌道后，內(nèi)部制約使系統(tǒng)恢復(fù)原狀。食餌-捕食者模型(Volterra)的缺點(diǎn)與改進(jìn)r1=1,N1=20,1=0.1,w=0.2,r2=0.5,2=0.18相軌線趨向極限環(huán)結(jié)構(gòu)穩(wěn)定精選課件相關(guān)Matlab知識(shí)命令：常微分方程的符號(hào)解函數(shù)：dsolve格式：r=dsolve(‘eq1,eq2,…’,’cond1,cond2,…’,’v’)說明：對(duì)給定的常微分方程（組）eq1,eq2,…中指定的符號(hào)自變量v，與給定的邊界條件和初始條件cond1,cond2,…，求符號(hào)解（即解析解）r；若沒有指定變量v，則默認(rèn)變量為t；在微分方程（組）的表達(dá)式eq中，大寫字母D表示對(duì)自變量（設(shè)為x）的微分算子：D=s/sx,D2=d2/x2,….微分算子D后面的字母則表示為因變量，即待求解的未知函數(shù).2012-3-1740精選課件相關(guān)Matlab知識(shí)初始和邊界條件由字符串表示：y(a)=b,Dy(c)=d,D2y(e)=f等，分別表示若邊界條件少于方程（組）的階數(shù)，則返回的結(jié)果r中會(huì)出現(xiàn)任意常數(shù)C1,C2,…；dsolve命令最多可以接受12個(gè)輸入變量（包括方程組與定解條件個(gè)數(shù)，當(dāng)然可以做到輸入的方程個(gè)數(shù)多于12個(gè)，只要將多個(gè)方程置于一個(gè)字符串內(nèi)即可）.若沒有給定輸出參量，則在命令窗口顯示解列表.2012-3-1741精選課件相關(guān)Matlab知識(shí)若該命令找不到解析解，則返回一警告信息，同時(shí)返回一空的sym對(duì)象.這時(shí)，用戶可以用命令ode23或ode45求解方程組的數(shù)值解.例1：解常微分方程：（注：【程序】『輸出結(jié)果』）（1）求的通解；【>>s=dsolve('Dy=a*y+b')】『s=-b/a+exp(a*t)*C1』2012-3-1742精選課件相關(guān)Matlab知識(shí)（2）初值問題【>>dsolve('Dy=y-2*t/y','y(0)=1')】『ans=(2*t+1)^(1/2)』（3）高階方程【>>dsolve('D2y=cos(2*x)-y','y(0)=1','Dy(0)=0','x')】『ans=4/3*cos(x)-1/3*cos(2*x)』2012-3-1743精選課件相關(guān)Matlab知識(shí)（4）邊值問題【>>dsolve('x*D2y-3*Dy=x^2','y(1)=0','y(5)=0','x')】『ans=31/468*x^4-1/3*x^3+125/468』（5）方程組問題【>>S=dsolve('Df=f+g','Dg=-f+g','f(0)=1','g(0)=2');>>%S是一個(gè)結(jié)構(gòu)>>S.f,S.g】2012-3-1744精選課件相關(guān)Matlab知識(shí)『ans=exp(t)*(2*sin(t)+cos(t))ans=exp(t)*(2*cos(t)-sin(t))』（6）無解析解問題【>>dsolve('Dy=x+y^2','y(0)=0','x')】『ans=(3^(1/2)*AiryAi(1,-x)+AiryB