數(shù)學建模實例

微分方程模型二、微分方程模型三、微分方程案例分析一、微分方程建模簡介四、微分方程的MATLAB求解五、微分方程綜合案例分析5/8/20241數(shù)學建模實例微分方程是研究變化規(guī)律的有力工具，在科技、工程、經(jīng)濟管理、生態(tài)、環(huán)境、人口和交通各個領(lǐng)域中有廣泛的應用。不少實際問題當我們采用微觀眼光觀察時都遵循著下面的模式：凈變化率＝輸入率－輸出率（守恒原理）一、微分方程模型簡介5/8/20242數(shù)學建模實例引例一在凌晨1時警察發(fā)現(xiàn)一具尸體，測得尸體溫度是29oC，當時環(huán)境的溫度是21oC。1h后尸體溫度下降到27oC，若人體的正常溫度是37oC，估計死者的死亡時間。解：設(shè)T(t)為死者在被殺害后t時刻尸體的溫度；k為比例系數(shù)。由牛頓冷卻定律，得則通解為5/8/20243數(shù)學建模實例由已知，由因此死者大約是在前一天的夜晚10:35被害的。可得微分方程的特解：，代入解得圖1尸體的溫度下降曲線5/8/20244數(shù)學建模實例建立微分方程的常用方法1、按變化規(guī)律直接列方程，如：利用人們熟悉的力學、數(shù)學、物理、化學等學科中的規(guī)律，如牛頓第二定律，放射性物質(zhì)的放射規(guī)律等。對某些實際問題直接列出微分方程．2、利用微元分析方法建模根據(jù)已知的規(guī)律或定理，通過尋求微元之間的關(guān)系式得出微分方程。3、模擬近似法，如：在生物、經(jīng)濟等學科中，許多現(xiàn)象所滿足的規(guī)律并不很清楚，而且現(xiàn)象也相當復雜，因而需根據(jù)實際資料或大量的實驗數(shù)據(jù)，提出各種假設(shè)，在一定的假設(shè)下，給出實際現(xiàn)象所滿足的規(guī)律，然后利用適當?shù)臄?shù)學方法得出微分方程。5/8/20245數(shù)學建模實例微分方程的建模步驟1、翻譯或轉(zhuǎn)化：在實際問題中許多表示導數(shù)的常用詞，如“速率”、‘增長”(在生物學以及人口問題研究中)，“衰變”(在放射性問題中)，以及“邊際的”(在經(jīng)濟學中)等．2、建立瞬時表達式：根據(jù)自變量有微小改變△t時，因變量的增量△W，建立起在時段△t上的增量表達式，令△t→0，即得到的表達式．二、微分方程模型5/8/20246數(shù)學建模實例3、配備物理單位：在建模中應注意每一項采用同樣的物理單位．4、確定條件：這些條件是關(guān)于系統(tǒng)在某一特定時刻或邊界上的信息，它們獨立于微分方程而成立，用以確定有關(guān)的常數(shù)。為了完整充分地給出問題的數(shù)學陳述，應將這些給定的條件和微分方程一起列出。5/8/20247數(shù)學建模實例案例1：一位女士每天攝入2500cal食物，1200cal用于基本新陳代謝(即自動消耗)，并以每千克體重消耗16cal用于日常鍛煉，其他的熱量轉(zhuǎn)變?yōu)樯眢w的脂肪（設(shè)10000cal可轉(zhuǎn)換成1kg脂肪）。星期天晚上，該女士的體重是57.1526kg，星期四那天她飽餐了一頓，共攝入了3500cal的食物，要求建立一個通過時間預測體重的數(shù)學模型，并用它估計：（1）星期六該女士的體重？（2）為了不增重，每天她最多的攝入量是多少？（3）若不進食，N周后她的體重是多少？二、微分方程案例分析5/8/20248數(shù)學建模實例解1、翻譯或轉(zhuǎn)化：2、配備物理單位：3、建立表達式：4、確定條件：5/8/20249數(shù)學建模實例1、“每天”：體重的變化＝輸入一輸出其中輸入指扣除了基本新陳代謝之后的凈重量吸收；輸出是進行健身訓練時的消耗．2、上述陳述更好的表示結(jié)構(gòu)式：取天為計時單位，記W(t)為t天時體重(kg)，則：每天的凈吸收量＝2500–1200＝1300（cal)每天的凈輸出量＝16(cal)×W＝16W(cal)轉(zhuǎn)換成脂肪量＝1300–16W（cal）3、體重的變化／天＝(千克／天)5/8/202410數(shù)學建模實例1、翻譯或轉(zhuǎn)化：2、配備物理單位：3、建立表達式：4、確定條件：5/8/202411數(shù)學建模實例有些量是用能量(cal)的形式給出的，而另外一些量是用重量的形式(cal)給出，考慮單位的匹配，利用單位匹配5/8/202412數(shù)學建模實例1、翻譯或轉(zhuǎn)化：2、配備物理單位：3、建立表達式：4、確定條件：5/8/202413數(shù)學建模實例建立表達式積分后可求得其通解為：（1）當時，每天體重的變化：初始條件為：，代入解出則5/8/202414數(shù)學建模實例積分后可求得其通解為：（2）當時，每天體重的變化：初始條件為：，代入解出則5/8/202415數(shù)學建模實例積分后可求得其通解為：（2）當時，食物的攝入量恢復正常初始條件為：，代入解出則5/8/202416數(shù)學建模實例最后得到不同階段的微分方程是：5/8/202417數(shù)學建模實例（1）代入對應方程，求得現(xiàn)回答上述問題（2）要滿足體重不增，即所以因此每天總卡路里攝取量是1200+914＝2114cal(cal)（3）由于每天不攝取能量，所以解得因此，n周后的體重為5/8/202418數(shù)學建模實例案例2在一個巴基斯坦洞穴里，發(fā)現(xiàn)了具有古代尼安德特人特征的人骨碎片，科學家們把它們帶到實驗室，作碳14年代測定。分析表明C14與C12的比例僅僅是活組織內(nèi)的6.24％，此人生活在多少年前？（碳14年代測定：活體中的碳有一小部分是放射性同位素C14。這種放射性碳是由于宇宙射線在高層大氣中的撞擊引起的，經(jīng)過一系列交換過程進入活組織中，直到在生物體中達到平衡濃度。這意味著在活體中，C14的數(shù)量與穩(wěn)定的C12的數(shù)量成定比。生物體死亡后，交換過程就停止了，放射性碳便以每年八千分之一的速度減少）5/8/202419數(shù)學建模實例（1）問題分析與模型的建立1、放射性衰變的這種性質(zhì)還可描述為“放射性物質(zhì)在任意時刻的衰變速度都與該物質(zhì)現(xiàn)存的數(shù)量成比例”。而C14的比例數(shù)為每年八千分之一。2、碳14年代測定可計算出生物體的死亡時間；所以，我們問題實際上就是：“這人死去多久了？”若設(shè)t為死后年數(shù)，y(t)為比例數(shù)，則y(t)=C14/C12(mgC14/mgC12)，則上文中最后一句話就給出了我們的微分方程，單位為mgC14/mgC12/yr(與關(guān)鍵詞“速率”相當)5/8/202420數(shù)學建模實例（2）解微分方程的通解為：由初始條件，故有由問題，當，代入原方程5/8/202421數(shù)學建模實例案例3、追線問題我緝私艦雷達發(fā)現(xiàn)，在其正西方距c海里處有一艘走私船正以勻速度a沿直線向北行駛，緝私艦立即以最大的速度b追趕，若用雷達進行跟蹤，保持船的瞬時速度方向始終指向走私船，試求緝私艦追逐路線和追上的時間。5/8/202422數(shù)學建模實例圖2走私船與緝私艦的位置關(guān)系(c,0)xD(x,y)走私船R(0,at)緝私艇O5/8/202423數(shù)學建模實例幾何關(guān)系5/8/202424數(shù)學建模實例如何消去時間t?1、求導：2、速度與路程的關(guān)系：3、分解得：（這里有負號是因為s隨x的減小而增大）4、將第2、3步代入第1步，可得模型5/8/202425數(shù)學建模實例追線模型：模型的解：5/8/202426數(shù)學建模實例解的進一步討論（1）若a<b，從而k<1，由積分式得當x=0時，即走私船被緝私艦捕捉前所花的時間為所跑過的距離為（2）若a=b，即k=1，由積分式得顯然x不能取零值，即緝私艦不可能追上走私船。（3）若a>b，即k>1，顯然緝私艦也不可能追上走私船。5/8/202427數(shù)學建模實例如圖所示一個容量為2000m3的小湖的示意圖，通過小河A水以0.1m3/s的速度流入，以相同的流量湖水通過B流出。在上午11:05時，因交通事故一個盛有毒性化學物質(zhì)的容器傾翻，圖中X點處注入湖中。在采取緊急措施后，于11:35事故得到控制，但數(shù)量不詳案例4湖泊污染問題的化學物質(zhì)Z已瀉入湖中，初步估計Z的量在5~20m3之間。建立一個模型，通過它來估計湖水污染程度隨時間的變化并估計：（1）湖水何時到達污染高峰；（2）何時污染程度可降至安全水平(<0.05%)圖3小湖示意圖5/8/202428數(shù)學建模實例湖泊污染問題分析設(shè)湖水在t時的污染程度為C(t)，即每立方米受污染的水中含有Cm3的化學物質(zhì)和(1-C)m3的清潔水。用分鐘作為時間t的單位。在0<t<30的時間內(nèi)，污染物流入湖中小湖示意圖的速率是Z/30(m3/min)，而排出湖外的污染物的速率是60×0.1C(m3/min)，因為每立方流走的水中含有Cm3的污染物，而湖水始終保持2000m3的容積不變，所以可列方程：5/8/202429數(shù)學建模實例由初始條件：，可得微分方程的特解為顯然，t=30時，污染達到高峰，所以因污染源被截斷，故微分方程變?yōu)樗奶亟鉃椋汉泻廴疚锏乃矔r變化率=污染物流入量－污染物排出量5/8/202430數(shù)學建模實例當達到安全水平，即C(t)=0.0005時，可求出此時的t=T，即解得Z取不同值時的濃度C(30)和時間T51015200.002390.004780.007170.0095655273891810145/8/202431數(shù)學建模實例何為房室系統(tǒng)？在用微分方程研究實際問題時，人們常常采用一種叫“房室系統(tǒng)”的觀點來考察問題。根據(jù)研究對象的特征或研究的不同精度要求，我們把研究對象看成一個整體（單房室系統(tǒng)）或?qū)⑵淦史殖扇舾蓚€相互存在著某種聯(lián)系的部分（多房室系統(tǒng)）。房室具有以下特征：它由考察對象均勻分布而成，（注：考察對象一般并非均勻分布，這里采用了一種簡化方法一集中參數(shù)法）；房室中考察對象的數(shù)量或濃度（密度）的變化率與外部環(huán)境有關(guān)，這種關(guān)系被稱為“交換”且交換滿足著總量守衡。在本節(jié)中，我們將用房室系統(tǒng)的方法來研究藥物在體內(nèi)的分布。在下一節(jié)中，我們將用多房室系統(tǒng)的方法來研究另一問題。兩者都很簡單，意圖在于介紹建模方法。交換環(huán)境內(nèi)部單房室系統(tǒng)均勻分布案例5藥物在體內(nèi)的變化（房室模型）5/8/202432數(shù)學建模實例藥物的分解與排泄（輸出）速率通常被認為是與藥物當前的濃度成正比的，即：藥物分布的單房室模型單房室模型是最簡單的模型，它假設(shè)：體內(nèi)藥物在任一時刻都是均勻分布的，設(shè)t時刻體內(nèi)藥物的總量為x(t)；系統(tǒng)處于一種動態(tài)平衡中，即成立著關(guān)系式：藥物的輸入規(guī)律與給藥的方式有關(guān)。下面，我們來研究一下在幾種常見的給藥方式下體內(nèi)藥體的變化規(guī)律。機體環(huán)境藥物總量圖3-8

假設(shè)藥物均勻分布5/8/202433數(shù)學建模實例情況1快速靜脈注射機體環(huán)境只輸出不輸入房室其解為：藥物的濃度：

與放射性物質(zhì)類似，醫(yī)學上將血漿藥物濃度衰減一半所需的時間稱為藥物的血漿半衰期：負增長率的Malthus模型在快速靜脈注射時，總量為D的藥物在瞬間被注入體內(nèi)。設(shè)機體的體積為V，則我們可以近似地將系統(tǒng)看成初始總量為D，濃度為D/V，只輸出不輸入的房室，即系統(tǒng)可看成近似地滿足微分方程：（3.12）

5/8/202434數(shù)學建模實例情況2恒速靜脈點滴機體環(huán)境恒定速率輸入房室藥物似恒速點滴方式進入體內(nèi)，即:則體內(nèi)藥物總量滿足：（x(0)=0）（3.13）這是一個一階常系數(shù)線性方程，其解為：或易見：稱為穩(wěn)態(tài)血藥濃度對于多次點滴，設(shè)點滴時間為T1，兩次點滴之間的間隔時間設(shè)為T2，則在第一次點滴結(jié)束時病人體內(nèi)的藥物濃度可由上式得出。其后T2時間內(nèi)為情況1。故：（第一次）0≤t≤T1

T1≤t≤T1

+T2

類似可討論以后各次點滴時的情況，區(qū)別只在初值上的不同。第二次點滴起，患者體內(nèi)的初始藥物濃度不為零。5/8/202435數(shù)學建模實例情況3口服藥或肌注y(t)x(t)k1ykx環(huán)境機體外部藥物口服藥或肌肉注射時，藥物的吸收方式與點滴時不同，藥物雖然瞬間進入了體內(nèi)，但它一般都集中與身體的某一部位，靠其表面與肌體接觸而逐步被吸收。設(shè)藥物被吸收的速率與存量藥物的數(shù)量成正比，記比例系數(shù)為K1，即若記t時刻殘留藥物量為y(t)，則y滿足：D為口服或肌注藥物總量因而：所以：解得：從而藥物濃度：在通常情況下，總有k1>k（藥物未吸收完前，輸入速率通?？偞笥诜纸馀c排泄速率），但也有例外的可能（與藥物性質(zhì)及機體對該藥物的吸收、分解能力有關(guān)）。當k1>k時，體內(nèi)藥物量均很小，這種情況在醫(yī)學上被稱為觸發(fā)翻轉(zhuǎn)（flip-flop）。當k1=k時，對固定的t，令k→k1取極限（應用羅比達法則），可得出在這種情況下的血藥濃度為：5/8/202436數(shù)學建模實例如下圖給出了上述三種情況下體內(nèi)血藥濃度的變化曲線。容易看出，快速靜脈注射能使血藥濃度立即達到峰值，常用于急救等緊急情況；口服、肌注與點滴也有一定的差異，主要表現(xiàn)在血藥濃度的峰值出現(xiàn)在不同的時刻，血藥的有效濃度保持時間也不盡相同，（注：為達到治療目的，血藥濃度應達到某一有效濃度，并使之維持一特定的時間長度）。房室系統(tǒng)我們已求得三種常見給藥方式下的血藥濃度C(t)，當然也容易求得血藥濃度的峰值及出現(xiàn)峰值的時間，因而，也不難根據(jù)不同疾病的治療要求找出最佳治療方案。5/8/202437數(shù)學建模實例國家質(zhì)量監(jiān)督檢驗檢疫局2004年5月31日發(fā)布了新的《車輛駕駛?cè)藛T血液、呼氣酒精含量閾值與檢驗》國家標準，新標準規(guī)定，車輛駕駛?cè)藛T血液中的酒精含量大于或等于20毫克／百毫升，小于80毫克／百毫升為飲酒駕車（原標準是小于100毫克／百毫升），血液中的酒精含量大于或等于80毫克／百毫升為醉酒駕車（原標準是大于或等于100毫克／百毫升）。五、微分方程綜合案例分析5/8/202438數(shù)學建模實例大李在中午12點喝了一瓶啤酒，下午6點檢查時符合新的駕車標準，緊接著他在吃晚飯時又喝了一瓶啤酒，為了保險起見他呆到凌晨2點才駕車回家，又一次遭遇檢查時卻被定為飲酒駕車，這讓他既懊惱又困惑，為什么喝同樣多的酒，兩次檢查結(jié)果會不一樣呢？請你參考下面給出的數(shù)據(jù)（或自己收集資料）建立飲酒后血液中酒精含量的數(shù)學模型，對大李碰到的情況做出解釋.5/8/202439數(shù)學建模實例參考數(shù)據(jù)1.人的體液占人的體重的65%至70%，其中血液只占體重的7%左右；而藥物（包括酒精）在血液中的含量與在體液中的含量大體是一樣的。2.體重約70kg的某人在短時間內(nèi)喝下2瓶啤酒后，隔一定時間測量他的血液中酒精含量（毫克／百毫升），得到數(shù)據(jù)如下：時間(小時)0.250.50.7511.522.533.544.55酒精含量306875828277686858515041時間(小時)678910111213141516酒精含量38352825181512107745/8/202440數(shù)學建模實例問題分析一個人的血中酒精含量取決于他喝了多少酒、他體內(nèi)原有的酒精含量以及喝酒方式等。由科普知識知道，酒精是經(jīng)胃腸（主要是肝臟）的吸收與分解進體液的。因此本文把酒精的從胃腸（含肝臟）向體液轉(zhuǎn)移情況用如下簡圖直觀地表示：k11為酒精從胃腸滲透到（除體液外）其它地方的速率系數(shù)；k12為酒精從胃腸進入體液的速率系數(shù)；k21為酒精在體液中消耗（向外排除或分解或吸收）的速率系數(shù)；f(t)為酒精進入胃腸的速率。5/8/202441數(shù)學建模實例由題意，參照房室模型，可建立如下微分方程組：（1）大李在中午12點喝一瓶啤酒時，即在t=0時，胃腸中的酒精量x1(0)為一瓶酒中的酒精a與飲酒瓶數(shù)N的乘積Na，而此時體液中的酒精量y1(0)為0。因此初始條件為體液（或血液）中的酒精的濃度為5/8/202442數(shù)學建模實例（2）大李第二次喝酒時胃腸和體液中已經(jīng)有酒精，所以在第二次喝酒即t=0時胃腸中的酒精量x2(0)為N瓶酒中的酒精質(zhì)量Na與第一次喝酒后殘留在胃腸中的酒精質(zhì)量x1(T1)之和，而此時體液中的酒精量y1(0)為第一次喝酒后殘留在胃腸中的酒精質(zhì)量y1(T1)，因此大李第二次喝酒的模型如下：5/8/202443數(shù)學建模實例解以上微分方程組，得令，，解可轉(zhuǎn)化為5/8/202444數(shù)學建模實例N＝2，運用最小二乘擬合法，求解得作圖如下：5/8/202445數(shù)學建模實例將以上數(shù)據(jù)代入問題一的模型中，可求得大李在中午12點飲一瓶啤酒，即N=1時，到下午6點第一次檢查時體液中的酒精含量（即血液中的酒精含量）所以大李通過了第一次檢查。5/8/202446數(shù)學建模實例大李第二次喝酒模型的方程解為：考慮到大李在下午6點接受檢查，之后由于停車等待等原因耽誤了大約半個小時，假設(shè)大李從第一次檢驗到第二次喝酒之間間隔0.5小時，代入數(shù)據(jù)計算可得第二次檢驗時，大李血液中酒精含量為：20.2448(毫克/百毫升)。這就解釋了大李在第一次喝酒通過檢查，第二次喝同樣的酒且經(jīng)過更長的時間檢查卻被定為飲酒駕車的情況，因為第二次喝酒時有第一次喝酒的殘留量。5/8/202447數(shù)學建模實例5/8/202448數(shù)學建模實例問題及其背景降落傘（parachute）是由柔性紡織物制成的傘狀氣動力減速器，可分為人用傘和物用傘，在軍事、搶險救災等方面有著廣泛的用途。人用傘使用時通常離開飛行器后馬上張開，但有時因為特殊需要，要求跳傘者從高空跳下，降落適當距離后才在半空中手動開傘，一般要求跳傘者在空中停留的時間盡可能短，又能以安全的速度落地。那么，何時張傘好呢？降落傘何時張開好5/8/202449數(shù)學建模實例分析我們主要關(guān)心什么呢？是跳傘者的落

地速度和在空中的停留時間。因此我們首先要考慮跳傘者的降落速度，它是時間的函數(shù)。跳傘者（包括降落傘，下同）在降落過程中主要受到重力和空氣阻力的作用以及氣流運動的影響，一般所受到的空氣阻力與降落速度成正比。因為我們主要關(guān)心一般情況下降落速度的垂直分量變化情況，可以忽略水平分量，不考慮氣流運動的影響，只考慮其作垂直降落運動。5/8/202450數(shù)學建模實例根據(jù)牛頓運動學第二定律以及加

速度是速度關(guān)于時間的導數(shù)，我們

就能列出降落速度滿足的微分方程。雖然張傘前后跳傘者所受到的空氣阻力的情況差別較大，但可認為僅是所受空氣阻力與降落速度的比例系數(shù)不同而已，這樣，如果忽略張傘時間，可以分張傘前和張傘后來建立類似的模型。下面主要建立張傘后的模型。5/8/202451數(shù)學建模實例模型的假設(shè)跳傘者（包括降落傘）在降落過程中只受到重力和空氣阻力的作用，只作垂直降落運動。所受到的空氣阻力的大小與降落速率成正比，比例系數(shù)是與時間無關(guān)的常數(shù)，設(shè)為k。張傘時刻為t=0，此時降落速率為v0。5/8/202452數(shù)學建模實例模型的建立設(shè)跳傘者（包括降落傘）的質(zhì)量為重力加速度為，降落速度為由Newton力學第二定律，可得：，即

這就是跳傘者的降落速度滿足的數(shù)學模型，這是一個常微分方程的初值問題。5/8/202453數(shù)學建模實例模型的求解如何求解這個模型呢？注意到即，我們有其中為任意常數(shù)。將定解條件代入，即得因此，模型的解為：5/8/202454數(shù)學建模實例模型解的分析和應用因為，隨著時間的增大，降落速度將很快趨于常值落地有足夠的時間，那么落地時的降落。。如果從張傘到速度約等于5/8/202455數(shù)學建模實例的大小與傘張開時傘面的形狀和有效面積有關(guān)，考慮安全等原因，張傘經(jīng)過10秒后約下降速度約6.0000076米/秒，通常設(shè)計降落傘使得，若米/秒，117.55米，此時降落已非常接近6米/秒的速度。所以，一個經(jīng)過一定訓練的跳傘者，若從離地面8000米的高空跳傘，即使離地面只有幾百米時才張開傘，也能安全著地，而這和剛跳落時就張傘相比，空中滯留時間將大大減少。5/8/202456數(shù)學建模實例進一步的考慮

上述模型可根據(jù)需要作進一步的改進。比如，因為空氣的稀薄程度與海拔高度有關(guān)，因此的大小也與海拔有一定關(guān)系，可認為此時可轉(zhuǎn)而考慮降落速度與的關(guān)系。高度，設(shè)時開始張傘，此時速度為因為，5/8/202457數(shù)學建模實例所以原模型可以改寫為：同學們可以根據(jù)具體情況分析求解此模型。5/8/202458數(shù)學建模實例求微分方程（組）解析解的命令:dsolve(‘方程1’,‘方程2’,…,‘方程n’,‘初始條件’,‘自變量’)ToMATLAB（ff1）結(jié)果：u=tan(t-c)五、微分方程的MATLAB求解ezplot5/8/202459數(shù)學建模實例解輸入命令:y=dsolve('D2y+4*Dy+29*y=0','y(0)=0,Dy(0)=15','x')結(jié)果為:y=3e-2xsin（5x）ToMATLAB（ff2）5/8/202460數(shù)學建模實例解輸入命令：

[x,y,z]=dsolve('Dx=2*x-3*y+3*z','Dy=4*x-5*y+3*z','Dz=4*x-4*y+2*z','t')；x=simple(x)%將x化簡y=simple(y)z=simple(z)結(jié)果為：x=(c1-c2+c3+c2e-3t-c3e-3t)e2t

y=-c1e-4t+c2e-4t+c2e-3t-c3e-3t+c1-c2+c3)e2t

z=(-c1e-4t+c2e-4t+c1-c2+c3)e2t

ToMATLAB（ff3）返回5/8/202461數(shù)學建模實例微分方程的數(shù)值解（一）常微分方程數(shù)值解的定義在生產(chǎn)和科研中所處理的微分方程往往很復雜,且大多得不出一般解．而實際中的對初值問題，一般是要求得到解在若干個點上滿足規(guī)定精確度的近似值，或者得到一個滿足精確度要求的便于計算的表達式．因此，研究常微分方程的數(shù)值解法是十分必要的．返回5/8/202462數(shù)學建模實例（二）用MAT