淺談常微分方程的數(shù)值解法及其應(yīng)用[畢業(yè)論文]

1、<p>　　本科畢業(yè)論文（設(shè)計(jì)）</p><p><b>　?。?0 屆）</b></p><p>　　淺談常微分方程的數(shù)值解法及其應(yīng)用</p><p>　　所在學(xué)院 </p><p>　　專業(yè)班級 信息與計(jì)算科學(xué) </p

2、><p>　　學(xué)生姓名 學(xué)號 </p><p>　　指導(dǎo)教師 職稱 </p><p>　　完成日期 年 月 </p><p>　　摘要：對于一些簡單而典型的微分方程模型，譬如線性方程、某些特殊的一階非線性方程等是

3、可以設(shè)法求出其解析解的，但在數(shù)學(xué)建模中碰到的常微分方程初值問題模型，通常很難，甚至根本無法求出其解析解，而只能求其近似解。因此，研究其數(shù)值方法，以便快速求得數(shù)值有其重大意義。針對于此，本文對常微分方程初值問題模型現(xiàn)有的數(shù)值解法問題進(jìn)行了綜述研究，探討了一些數(shù)值解法的應(yīng)用。</p><p>　　關(guān)鍵詞：常微分方程 ，數(shù)值解法，應(yīng)用</p><p>　　Numerical Solution o

4、f Differential Equation and Its Applications</p><p>　　Abstract: for some simple and typical differential equation model,such as linear equation, some special first-order nonlinear equation can be managed to

5、find out its analytical solution,but in mathematical modeling of ordinary differential equation met in initial value problem model,it is often hard to,even can't find out the analytical solution,but only for its appr

6、oximate solution.Therefore, study the numerical method for quick get numerical has its great significance. Based on this,the paper </p><p>　　Keywords: ordinary differential equation, the numerical solution,

7、applications</p><p><b>　　目 錄</b></p><p><b>　　1 緒 論1</b></p><p>　　1.1問題的背景1</p><p>　　1.2問題的意義1</p><p>　　2 常微分方程概念介紹3</p>

8、;<p>　　2.1常微分方程概況3</p><p>　　2.2常微分方程初值問題描述3</p><p>　　2.3初值問題(1)解的存在惟一性定理3</p><p>　　2.4常微分微分方程產(chǎn)生的歷史背景以及發(fā)展4</p><p>　　3 常微分方程的數(shù)值解法6</p><p>　　3.1常

9、微分方程求解的數(shù)學(xué)思想6</p><p>　　3.2常微分方程的數(shù)值解法6</p><p>　　3.2.1 Euler 法6</p><p>　　3.2.2 泰勒級數(shù)法8</p><p>　　3.2.3 龍格—庫塔方法9</p><p>　　3.2.4 預(yù)報(bào)—校正方法11</p><p

10、>　　4 常微分方程的數(shù)值解法的應(yīng)用13</p><p>　　4.1初軌計(jì)算13</p><p>　　4.2司機(jī)飲酒駕車防避模型16</p><p>　　5 結(jié)論錯(cuò)誤!未定義書簽。</p><p><b>　　參考文獻(xiàn)24</b></p><p><b>　　致 謝

11、26</b></p><p><b>　　1 緒 論</b></p><p>　　1.1 問題的背景</p><p>　　微分方程差不多是和微積分同時(shí)產(chǎn)生的，它的形成和發(fā)展是和力學(xué)、天文學(xué)、物理學(xué)，以及其他科學(xué)技術(shù)的發(fā)展密切相關(guān)的。蘇格蘭數(shù)學(xué)家耐普爾創(chuàng)立對數(shù)的時(shí)候，就討論過微分方程的近似解。牛頓在建立微積分的同時(shí)，對簡單的微分方程

12、用級數(shù)來求解。后來瑞士數(shù)學(xué)家雅各布·利、歐拉、法國數(shù)學(xué)家克雷洛、達(dá)朗貝爾、拉格朗日等人又不斷地研究和豐富了微分方程的理論。[1]</p><p>　　常微分方程在很多學(xué)科領(lǐng)域內(nèi)有著重要的作用，自動(dòng)控制、各種電子學(xué)裝置的設(shè)計(jì)、彈道的計(jì)算、飛機(jī)和導(dǎo)彈飛行的穩(wěn)定性的研究、化學(xué)反應(yīng)過程穩(wěn)定性的研究等等，這些問題都可以化為求微分方程的解，或者化為研究解的性質(zhì)的問題。[2]</p><p>

13、<b>　　1.2問題的意義</b></p><p>　　常微分方程的概念、解法和相關(guān)理論很多。求通解在歷史上曾作為微分方程的主要目標(biāo)，不過能夠求出通解的情況不多，在實(shí)際應(yīng)用中多是求滿足某種指定條件的特解。我們知道，自然界中很多事物的運(yùn)動(dòng)規(guī)律可用微分方程來刻畫。常微分方程是研究自然科學(xué)和社會(huì)科學(xué)中的事物、物體和現(xiàn)象運(yùn)動(dòng)、演化和變化規(guī)律的最為基本的數(shù)學(xué)理論和方法。物理、化學(xué)、生物、工程、航空航

14、天、醫(yī)學(xué)、經(jīng)濟(jì)和金融領(lǐng)域中的許多原理和規(guī)律都可以描述成適當(dāng)?shù)某Ｎ⒎址匠?，如牛頓的運(yùn)動(dòng)定律、萬有引力定律、機(jī)械能守恒定律，能量守恒定律、人口發(fā)展規(guī)律、生態(tài)種群競爭、疾病傳染、遺傳基因變異、股票的漲伏趨勢、利率的浮動(dòng)、市場均衡價(jià)格的變化等，對這些規(guī)律的描述、認(rèn)識和分析就歸結(jié)為對相應(yīng)的常微分方程描述的數(shù)學(xué)模型的研究。因此，常微分方程的理論和方法不僅廣泛應(yīng)用于自然科學(xué)，而且越來越多的應(yīng)用于社會(huì)科學(xué)的各個(gè)領(lǐng)域。它的學(xué)術(shù)價(jià)值是無價(jià)的，應(yīng)用價(jià)值是立竿

15、見影的。求一階常微分方程的解是數(shù)學(xué)工作者的一項(xiàng)基本的且重要的工作。由于國內(nèi)外眾多數(shù)學(xué)家的努力，使此學(xué)科基本上形成了一套完美的學(xué)科體系；由于該問題比較復(fù)雜且涉及的面廣，使得有些問題的解析解很難求出，而對于一些典型的微分方程(如</p><p>　　“常微分方程”是理學(xué)院數(shù)學(xué)系所有專業(yè)學(xué)生的重要專業(yè)基礎(chǔ)課之一，也是工科、經(jīng)濟(jì)等專業(yè)必學(xué)內(nèi)容之一。由此可見，“常微分方程”在我們的理論學(xué)習(xí)和實(shí)際應(yīng)用中具有很重要的地位和實(shí)用

16、性，成為我們發(fā)現(xiàn)定理和解決問題的重要方法之一。[4]</p><p>　　2常微分方程概念介紹</p><p>　　2.1 常微分方程概況</p><p>　　我們知道微分方程就是聯(lián)系著自變量、未知函數(shù)及其導(dǎo)數(shù)的關(guān)系式。如果在微分方程中，自變量的個(gè)數(shù)只有一個(gè)，我們稱這種微分方程為常微分方程；自變量的個(gè)數(shù)為兩個(gè)或兩個(gè)以上的微分方程稱為偏微分方程。</p>

17、<p><b>　　方程</b></p><p>　　， （1）</p><p>　　為常微分方程。其中出現(xiàn)的最高階導(dǎo)數(shù)的階數(shù)，叫做常微分方程的階。例如  ， ，是一階常微分方程。 是二階常微分方程。設(shè)定義于區(qū)間上，有直到階的導(dǎo)數(shù)，將它代入（1），使（1）變成關(guān)于的恒等式，即</p><p>

18、<b>　　， </b></p><p>　　就稱＝為（1）的一個(gè)定義于上的解，并稱為該解的定義區(qū)間。[5]</p><p>　　2.2 常微分方程初值問題描述</p><p>　　在自然科學(xué)和經(jīng)濟(jì)的許多領(lǐng)域中。常常會(huì)遇到一階常微分方程的初值問題</p><p><b>　　(1)</b><

19、/p><p>　　這里是充分光滑，即關(guān)于或，滿足李普希茨條件的二元函數(shù)，是給定的初始值，稱為初始條件[5]。</p><p>　　2.3初值問題(1)解的存在惟一性定理</p><p>　　一個(gè)常微分方程是不是有特解呢?如果有，又有幾個(gè)呢?這是微分方程論中一個(gè)基本的問題，數(shù)學(xué)家把它歸納成基本定理，叫做存在和唯一性定理。因?yàn)槿绻麤]有解，而我們要去求解，那是沒有意義的；如果

20、有解而又不是唯一的，那又不好確定。因此，存在和唯一性定理對于微分方程的求解是十分重要的。這個(gè)重要的存在和唯一性就是下面列出的著名的存在惟一性定理。</p><p>　　定理 如果在矩形區(qū)域上連續(xù)且關(guān)于滿足利普希茨條件，即存在正常數(shù)，使得</p><p>　　對所有的以及所有、都成立，則（1）存在惟一的連續(xù)可微解[6]。</p><p>　　2.4常微分微分方程產(chǎn)生

21、的歷史背景以及發(fā)展</p><p>　　微分方程的形成與發(fā)展是和力學(xué)、天文學(xué)、物理學(xué)，以及其他科學(xué)技術(shù)的發(fā)展密切相關(guān)的。常微分方程是屬于數(shù)學(xué)分析的一支，是數(shù)學(xué)中與應(yīng)用密切的相關(guān)的基礎(chǔ)學(xué)科，其自身也在不斷發(fā)展中，學(xué)好常微分方程基本的理論與方法對進(jìn)一步學(xué)習(xí)研究數(shù)學(xué)理論和實(shí)際應(yīng)用均非常重要。而數(shù)學(xué)的其他分支的新發(fā)展，如復(fù)變函數(shù)、李群、組合拓?fù)鋵W(xué)等，都對微分方程的發(fā)展產(chǎn)生了深刻的影響，當(dāng)前計(jì)算機(jī)的發(fā)展更是為常微分方程的應(yīng)

22、用及理論研究提供了非常有力的工具。 [7]</p><p>　　常微分方程發(fā)展的初期是對具體的常微分方程希望能用初等函數(shù)或超越函數(shù)表示其解，屬于“求通解”時(shí)代。萊布尼茨成專門研究利用變量變換解決一階微分方程的求解問題，而歐拉則試圖用積分因子統(tǒng)一處理，伯努利、里卡蒂微分方程就是在研究初等積分時(shí)提出后人以他們的名字命名的方程。[8]</p><p>　　早期的常微分方程的求解熱潮被劉維爾在18

23、41年證明里卡蒂方程不存在一般的初等解而中斷。加上柯西初值問題的提出，常微分方程從“求通解”，轉(zhuǎn)向“求定解”時(shí)代。同時(shí)，由于天文計(jì)算的需要促進(jìn)了常微分方程攝動(dòng)理論以及小參數(shù)冪級數(shù)等近似方法的研究。[8]</p><p>　　19世紀(jì)末。天體力學(xué)中的太陽系穩(wěn)定性問題需研究常微分方程解的大范圍性態(tài)，從而使常微分方程的研究從“求定解問題”轉(zhuǎn)化為“求所有解時(shí)代”。[8]</p><p>　　20世

24、紀(jì)六七十年代以后，常微分方程由于計(jì)算機(jī)技術(shù)的發(fā)展迎來了新的時(shí)期，從“求所有解”轉(zhuǎn)入求“求特殊解”時(shí)代，發(fā)現(xiàn)了具有新性質(zhì)的特殊的解和方程，如混沌（解）、奇異吸引子等。塞蒙斯成如此評價(jià)微分方程在數(shù)學(xué)中的地位：“300年來分析是數(shù)學(xué)里首要的分支，而微分方程又是分析的心臟，這是初等微積分的天然后的繼課，又是為了解物理科學(xué)的一門最重要的數(shù)學(xué)，而且在它所產(chǎn)生的較深的問題中，它又是高等分析里大部分思想和理論的根源?！?[9]</p>&

25、lt;p>　　3 常微分方程的數(shù)值解法</p><p>　　3.1 常微分方程求解的數(shù)學(xué)思想</p><p>　　從常微分發(fā)展歷程可以看出，化歸是常微分方程的重要數(shù)學(xué)思想方法，常數(shù)變易法、代換法、級數(shù)解法、逐次逼近法、算子法、相平面分析法等，都是用聯(lián)系、變化的觀點(diǎn)，有意識地將問題化繁為簡，化歸解決的。非齊次方程問題化為齊次方程問題，一階線性方程組化為一階線性方程問題， 高階方程問題

26、化為低階方程問題，在常微分方程發(fā)展的各個(gè)階段包含著這種化歸范例。</p><p>　　常系數(shù)非齊次線性微分方程，經(jīng)采用歐拉的待定指數(shù)函數(shù)法，將求解問題化歸為代數(shù)方程根的問題，從而省去了積分運(yùn)算，這是十分引人入勝的。皮卡逼近法，將微分方程的解問題化歸為積分方程的解問題，進(jìn)而化歸為一致收斂的函數(shù)列問題，完全符合化難為易，化未知為已知，化繁為簡的化歸原則。拉普拉斯變換將常系數(shù)線性非齊次微分方程的邊值問題，化歸為關(guān)于未知

27、函數(shù)的拉氏變換像函數(shù)的代數(shù)方程問題。</p><p>　　3.2 常微分方程的數(shù)值解法</p><p>　　3.2.1 Euler 法</p><p>　　Euler法是最簡單的數(shù)值方法，為求解良態(tài)初值問題，的區(qū)間。實(shí)際上，下面的過程不是要找到滿足該初值問題的可微函數(shù)，而是要生成點(diǎn)集，并且將這些點(diǎn)作為近似解，即。如何構(gòu)造“近似滿足微方程”的“點(diǎn)集”呢？首先為這些點(diǎn)

28、選擇橫坐標(biāo)，為方便起見，將區(qū)間劃分為個(gè)等距子區(qū)間，并選擇網(wǎng)絡(luò)點(diǎn)</p><p>　　， ，其中 （1）</p><p>　　值稱為步長。然后近似解</p><p>　　在上， （2）</p><p>　　設(shè)，和連續(xù)，利用泰勒定理將在處展開，對每個(gè)值，存

29、 在一個(gè)和之間的值，使得</p><p>　　， （3）</p><p>　　將和代人等式（3），得到的表示：</p><p>　　， （4）</p><p>　　如果步長足夠小，則可以忽略 2 次項(xiàng)（包含的項(xiàng)），得到</p><p>　　，

30、 （5）</p><p><b>　　這就是歐拉近似。</b></p><p>　　重復(fù)該過程，就能得到近似解曲線的一個(gè)點(diǎn)序列。歐拉方法的一般步驟是</p><p>　　， 其中 [8]（6）</p><p>　　例1 用歐拉方法求解初值問題，取步長，計(jì)算過程保留四位數(shù)字。</p>

31、<p>　　解 ，，首先建立歐拉迭代格式</p><p>　?。?）當(dāng)，時(shí)，已知，，有 ；</p><p>　?。?）當(dāng)，時(shí)，已知，，有 </p><p><b>　?。?lt;/b></p><p>　?。?）當(dāng)，時(shí)，已知，，有</p><p><b>　　；</b

32、></p><p>　　Euler算法的程序設(shè)計(jì)</p><p>　?。?） 送初值，打印</p><p><b>　?。?） </b></p><p><b>　　（3） </b></p><p><b>　?。?） </b></p

33、><p><b>　?。?） 打印</b></p><p><b>　?。?） </b></p><p>　?。?） 是轉(zhuǎn)（8）</p><p>　　(8） stop</p><p>　　3.2.2 泰勒級數(shù)法</p><p>　　泰勒級數(shù)

34、法有著廣泛的應(yīng)用，并且是比較求解初值問題的各種不同數(shù)值方法的標(biāo)準(zhǔn)，它可設(shè)計(jì)為任意指定的精度。下面首先將泰勒定理用新的公式表示，使之適合于求解微分方程。</p><p>　　定理（泰勒定理）設(shè) ，且在不動(dòng)點(diǎn)處有次泰勒級數(shù)展開：</p><p>　　， （1）</p><p><b>　　其中，</b></p><

35、;p>　　， （2）</p><p>　　表示函數(shù)關(guān)于的次全導(dǎo)數(shù)。求導(dǎo)公式可以遞歸地計(jì)算：</p><p><b>　?。?）</b></p><p><b>　　并且一般有</b></p><p>　　，

36、 （4）</p><p><b>　　其中為導(dǎo)數(shù)算子</b></p><p><b>　　，</b></p><p>　　區(qū)間上的初值問題的近似數(shù)值解可由各子區(qū)間上的公式（1）來推導(dǎo)。次泰勒方法的一般步驟為</p><p>　　， （5） </p&

37、gt;<p><b>　　其中在各步有。</b></p><p>　　次泰勒方法的最終全局誤差是階的，因此可選擇所需大小的，使得誤差足夠小。如果是固定，則理論上可以推導(dǎo)出步長，使之滿足任意的最終全局誤差。然而在實(shí)際運(yùn)算中，通常用和計(jì)算兩個(gè)近似結(jié)果集，然后比較其結(jié)果[9]。</p><p>　　3.2.3 龍格—庫塔方法</p><p&

38、gt;　　泰勒方法的優(yōu)點(diǎn)是最終全局誤差的階為，并且可以通過選擇較大的來得到較小的誤差。然而泰勒方法的缺點(diǎn)是，需要先確定，并且要計(jì)算高階導(dǎo)數(shù)，它們可能十分復(fù)雜。每個(gè)龍格一庫塔（Runge-Kutta ）方法都由一個(gè)合適的泰勒方法推導(dǎo)而來，使得其最終全局誤差為。一種折中方法是每步進(jìn)行若干次函數(shù)求值，從而省去高階導(dǎo)數(shù)計(jì)算。這種方法可構(gòu)造任意 階精度的近似公式。最常用的是的龍格一庫塔方法，它適用于一般的應(yīng)用，因?yàn)樗浅＞_、穩(wěn)定，且易于編程。許

39、多專家聲稱，沒有必要使用更高階的方法，因?yàn)樘岣叩木扰c增加的計(jì)算量相抵消。如果需要更高的精度，則應(yīng)該使用更小的步長或某種自適應(yīng)方法。</p><p>　　4階龍格一庫塔方法（RK4）可模擬的泰勒方法的精度。它基于如下對，的計(jì)算：</p><p>　　， （1）</p><p><b>　　其中，，和形如</b></p>

40、<p><b>　?。?）</b></p><p>　　通過與 階的泰勒級數(shù)方法的系數(shù)匹配，使得局部誤差為，龍格和庫塔得出了如下方程組：</p><p><b>　?。?）</b></p><p>　　該方程組有11個(gè)方程和13個(gè)未知量，必須補(bǔ)充兩個(gè)條件才可以求解。最有用的選擇是 </p>&l

41、t;p>　　， （4）</p><p><b>　　其余變量的解為</b></p><p><b>　?。?）</b></p><p>　　將式（4）和（5）中的值代入式（2）和式（1），得到標(biāo)準(zhǔn)的階龍格—庫塔方法，其描述如下。自初始點(diǎn)開始，利用</p>

42、<p>　　， （6）</p><p>　　生成近似值序列，其中[10]</p><p>　　3.2.4 預(yù)報(bào)—校正方法</p><p>　　歐拉方法、休恩方法、泰勒方法以及龍格一庫塔方法都稱為單步長方法，因?yàn)樗鼈冎焕们耙粋€(gè)點(diǎn)的信息來計(jì)算下一個(gè)點(diǎn)，即計(jì)算時(shí)只使用了初始點(diǎn)。一般地，只有用來。當(dāng)計(jì)算出若干個(gè)點(diǎn)之后，就可以利用幾個(gè)已

43、計(jì)算出的點(diǎn)來計(jì)算下一個(gè)點(diǎn)。以亞當(dāng)斯一巴什福斯4步法的推導(dǎo)為例，計(jì)算需要，，和。該方法不是自啟動(dòng)的，要生成點(diǎn)，必須先給出其4個(gè)初始點(diǎn)，，， （可用前面各節(jié)中的方法完成）。多步法的一個(gè)優(yōu)點(diǎn)是，可以確定它的局部截?cái)嗾`差（local truncation error ，簡稱 L.T.E.），并可以包含一個(gè)校正項(xiàng)，用于在每一步計(jì)算中改善解的精確度。該方法還可以確定步長是否小到能得到的精確值，同時(shí)又大到能夠免除不必要的和費(fèi)時(shí)的計(jì)算。使用預(yù)報(bào)子和校正

44、子的組合在每一步只需要進(jìn)行兩次函數(shù)求值[11]。</p><p>　　2.4.4.1亞當(dāng)斯一巴什福斯一莫爾頓方法</p><p>　　亞當(dāng)斯一巴什福斯一莫爾頓方法（Adams—Bashforth—Moulton）是由基本微積分定理推導(dǎo)出的多步法：</p><p>　　， (1)</p><p>　　預(yù)報(bào)子使用基于點(diǎn)和的

45、的拉格朗日多項(xiàng)式逼近值，并在區(qū)間上對式（1）積分，這個(gè)過程產(chǎn)生亞當(dāng)斯一巴什福斯預(yù)報(bào)子：</p><p>　　， （2）</p><p>　　校正子的推導(dǎo)類似。這時(shí)可以實(shí)用剛剛計(jì)算出的值。基于點(diǎn)，和新的點(diǎn)構(gòu)造的一個(gè)新的拉格朗日多項(xiàng)式逼近，然后在區(qū)間上對該多項(xiàng)式積分，即可得到亞當(dāng)斯一莫爾頓校正子[12]：</p><p>　　。

46、 （3）</p><p>　　2.4.4.2 米爾恩—辛普森方法</p><p>　　米爾恩—辛普森方法是預(yù)報(bào)子基于區(qū)間上的對的積分：</p><p>　　， （4）</p><p>　　預(yù)報(bào)子使用基于和的拉格朗日多項(xiàng)式逼近，在區(qū)間上對它積分，得到米爾恩預(yù)報(bào)子：</p><p>　　

47、， (5)</p><p>　　校正子的推導(dǎo)類似。此時(shí)值已知，基于點(diǎn)，和新點(diǎn)構(gòu)造的新的拉格朗日多項(xiàng)式，然后在區(qū)間上對該多項(xiàng)式積分，結(jié)果為大家所熟悉的辛普森公式[13]：</p><p>　　。 （6）</p><p>　　4 常微分方程的數(shù)值解法的應(yīng)用</p><p><b&g

48、t;　　4.1 初軌計(jì)算</b></p><p>　　根據(jù)理論力學(xué)，若將地球看作一個(gè)密度均勻分布的球體，則它對衛(wèi)星的吸引可等效于一個(gè)質(zhì)點(diǎn)，這樣地球和衛(wèi)星就構(gòu)成一個(gè)二體系統(tǒng)。在地心慣性坐標(biāo)系下考慮衛(wèi)星相對地球的運(yùn)動(dòng)，記衛(wèi)星位置矢量為，衛(wèi)星速度矢量，萬有引力常數(shù)為，地球質(zhì)量為，衛(wèi)星質(zhì)量為，根據(jù)萬有引力定律，衛(wèi)星受地球引力 </p>&

49、lt;p>　　? ，</p><p><b>　　其中。</b></p><p>　　根據(jù)牛頓第二定律，衛(wèi)星的運(yùn)動(dòng)方程可寫成</p><p><b>　　，</b></p><p>　　式中，是地球引力常數(shù)，為了計(jì)算方便，可選取軌道計(jì)算單

50、位，即取，這樣方程可寫成</p><p>　　， （1）</p><p>　　衛(wèi)星的運(yùn)動(dòng)狀態(tài)可由這6個(gè)量來表示，稱為軌道根數(shù)。其中稱為半長軸，稱為軌道偏心率，這2個(gè)量決定了軌道的大小和形狀，稱為軌道傾角，稱為升交點(diǎn)赤經(jīng)，這2個(gè)量確定了衛(wèi)星軌道平面的空間方向，稱為近地點(diǎn)幅角，稱為平近點(diǎn)角。</p><p>　　若衛(wèi)星某一時(shí)刻的，，

51、給定，相應(yīng)的軌道根數(shù)就完全確定，而任意時(shí)刻的衛(wèi)星位置矢量，速度矢量又完全由軌</p><p>　　道根數(shù) 來確定，因此、應(yīng)該可以用時(shí)刻的，來表達(dá)，即 </p><p><b>　　，</b></p><p>　　當(dāng)不太大時(shí)，可以在處作泰勒展開為的冪級數(shù)</p><p><b>　　，</b><

52、/p><p>　　根據(jù)二體問題運(yùn)動(dòng)方程（），的二次導(dǎo)數(shù)可化為不含微商的形式，的三次導(dǎo)數(shù)可化為及其一次導(dǎo)數(shù)的線性組合，推廣到任意次導(dǎo)數(shù)亦這樣，以分量為例</p><p><b>　　，</b></p><p>　　可以看到總能用形式來替換，故各高階導(dǎo)數(shù)總可以表示成不含和的形式： </p><p><b>　　，<

53、;/b></p><p>　　分量也可以表示這樣的形式，所以的時(shí)間冪級數(shù)又可以表示對，的展開式形式： </p><p><b>　　，</b></p><p><b>　　其中</b></p><p><b>　　，</b></p>

54、<p><b>　　另外、有封閉表達(dá)式</b></p><p>　　， （2）</p><p><b>　　其中，而滿足</b></p><p><b>　　， （3）</b></p><p>　　可用牛頓迭代法求解。</p><p

55、>　　對于也可以用相同的方法得到</p><p><b>　　，</b></p><p>　　當(dāng)有兩點(diǎn)、時(shí)刻的位置矢量、時(shí)，令=，則</p><p><b>　　，</b></p><p><b>　　于是可得</b></p><p><b&

56、gt;　　，</b></p><p>　　由于在計(jì)算時(shí)還有未知量，所以在計(jì)算上式時(shí)還有一個(gè)迭代過程，迭代第一步可取</p><p><b>　　，</b></p><p>　　在計(jì)算出、后即可得到計(jì)算的未知量，在計(jì)算出，進(jìn)行迭代，直到、與前一次計(jì)算的、的差達(dá)到精度要求為止。得到了，，，相應(yīng)的軌道根數(shù)就完全確定。</p>

57、<p>　　當(dāng)測量數(shù)據(jù)為測站地平坐標(biāo)系下的時(shí)，可直接轉(zhuǎn)換為地心坐標(biāo)系下的位置矢量，于是可以得到一系列的?？梢韵扔蒙厦娴姆椒ㄊ褂脙牲c(diǎn)數(shù)據(jù)計(jì)算出一組初軌，然后使用多點(diǎn)數(shù)據(jù)采用最小二乘法對其修正，得到較精確的初軌，步驟如下：</p><p>　　(1)將測站地平坐標(biāo)系下的數(shù)據(jù)轉(zhuǎn)換為地心坐標(biāo)系下的位置矢量；</p><p>　　(2)使用兩點(diǎn)數(shù)據(jù)計(jì)算出；</p><

58、p>　　(3)由計(jì)算出軌道根數(shù)；</p><p>　　(4)對式(4.3)計(jì)算出時(shí)刻的；</p><p>　　(5)利用式(4.2)求出時(shí)刻的；</p><p>　　(6)計(jì)算下列值（其中 為數(shù)據(jù)總數(shù)）；</p><p><b>　　， ，</b></p><p><b>

59、　　， ， </b></p><p><b>　　， ；</b></p><p><b>　?。?）求，</b></p><p><b>　　， ；</b></p><p>　?。?）返回(3)進(jìn)行循環(huán)迭代，直到與前一次計(jì)算出的的差達(dá)到精度要求為

60、止； </p><p>　?。?）由計(jì)算出軌道根數(shù)。 </p><p>　　4.2 司機(jī)飲酒駕車防避模型</p><p>　　在2004年全國大學(xué)生數(shù)學(xué)建模競賽題中有一個(gè)關(guān)于司機(jī)飲酒駕車模型。</p><p><b>　　1 問題的提出</b></p><p>　　《車輛駕駛?cè)藛T血液、呼氣酒精含量

61、閩值與檢驗(yàn)》國家新標(biāo)準(zhǔn)規(guī)定，車輛駕駛?cè)藛T血液中的酒精含量大于或等于20毫克／百毫升，小于80毫克／百毫升為飲酒駕車，血液中的酒精含量大于或等于80毫克／百毫升為醉酒駕車。大李在中午12點(diǎn)喝了一瓶啤酒，下午6點(diǎn)檢查時(shí)符合新的駕車標(biāo)準(zhǔn)，緊接著他在吃晚飯時(shí)又喝了一瓶啤酒，為了保險(xiǎn)起見他呆到凌晨2點(diǎn)才駕車回家，又一次遭遇檢查時(shí)卻被定為飲酒駕車，這讓他既懊惱又困惑，為什么喝同樣多的酒，兩次檢查結(jié)果會(huì)不一樣呢?</p><p&g

62、t;　　請你參考下面給出的數(shù)據(jù)(或自己收集資料)建立飲酒后血液中酒精含量的數(shù)學(xué)模型，并討論以下問題：</p><p>　　<1>．對大李碰到的情況做出解釋：</p><p>　　<2>．在喝了3瓶啤酒或者半斤低度白酒后多長時(shí)間內(nèi)駕車就會(huì)違反上述標(biāo)準(zhǔn)，在以下情況下回答：</p><p>　　1)、酒是在很短時(shí)間內(nèi)喝的；</p>&

63、lt;p>　　2)、酒是在較長一段時(shí)間(比如2小時(shí))內(nèi)喝的。</p><p>　　<3>．怎樣估計(jì)血液中的酒精含量在什么時(shí)間最高。</p><p>　　<4>．根據(jù)你的模型論證：如果天天喝酒，是否還能開車?</p><p>　　<5>．根據(jù)你做的模型并結(jié)合新的國家標(biāo)準(zhǔn)寫一篇短文，給想喝一點(diǎn)酒的司機(jī)如何駕車提出忠告。</

64、p><p><b>　　參考數(shù)據(jù)</b></p><p>　　1)、人的體液占人的體重的65％至70％，其中血液只占體重的7％左右：而藥物(包</p><p>　　括酒精)在血液中的含量與在體液中的含量大體是一樣的。</p><p>　　2)、體重約70kg的某人在短時(shí)間內(nèi)喝下2瓶啤滔后，隔一定時(shí)間測量他的血液中滔</

65、p><p>　　精含量(毫克／百毫升)，得到數(shù)據(jù)如下； </p><p><b>　　2 模型假設(shè)：</b></p><p>　　<l>.駕駛司機(jī)沒有其他疾病，消化系統(tǒng)良好，屬于健康人群，其體重為70kg左右。</p><p>　　<2>.酒精在人體內(nèi)的循環(huán)系統(tǒng)分為胃腔系統(tǒng)(系統(tǒng)I)和體液系統(tǒng)(系統(tǒng)I

66、I)，兩個(gè)系統(tǒng)的容積(即血液體積或酒精分布容積)在過程中保持不變。</p><p>　　<3>.酒精從系統(tǒng)I向系統(tǒng)II的轉(zhuǎn)移的速率系數(shù)，及向體外的排出的速率系數(shù)，與該系統(tǒng)的酒精濃度成正比，這兩個(gè)速率系數(shù)也是由人體的身體機(jī)能所決定的常數(shù)。</p><p>　　<4>.循環(huán)過程只考慮由體外進(jìn)入胃腔，再由胃腔進(jìn)入體液，最后由體液排除體外。不考慮人體其他機(jī)體對酒精的吸收，體

67、液的變化可以忽略而保持一定。</p><p><b>　　<5>.符號說明：</b></p><p>　　：酒精進(jìn)入胃腔的速率，設(shè)為常數(shù)</p><p><b>　?。簻y試時(shí)間(小時(shí))</b></p><p><b>　　：飲酒時(shí)間(小時(shí))</b></p>

68、;<p>　　：時(shí)刻人體胃腔中的酒精含量(毫克／百毫升)</p><p>　?。何盖恢谐跏季凭?毫克)</p><p>　　：剛喝完酒時(shí)胃腔中的酒精量(毫克)</p><p>　?。壕凭晌盖晦D(zhuǎn)移至體液的速率系數(shù)</p><p>　　：酒精由體液排出體外的速率系數(shù)</p><p>　?。壕凭晌盖晦D(zhuǎn)移至

69、體液的轉(zhuǎn)移速率(毫克／小時(shí))</p><p>　?。簳r(shí)刻人體體液中的酒精含量(毫克／百毫升)</p><p>　?。簳r(shí)刻人體體液中酒精濃度(毫克／百毫升)</p><p>　?。喝梭w體液的體積(百毫升)</p><p>　?。后w液系統(tǒng)中初始酒精濃度(毫克／百毫升)記為</p><p>　　：酒精排出體外的速率(毫克／小

70、時(shí))</p><p><b>　　3 模型建立：</b></p><p>　　由酒精在人體內(nèi)吸收、轉(zhuǎn)移規(guī)律的特點(diǎn)，應(yīng)用藥物動(dòng)力學(xué)原理建立人體內(nèi)胃腔與體液循環(huán)系統(tǒng)模型，可用微分方程描述其動(dòng)態(tài)過程。一般情長時(shí)間飲酒，原身體內(nèi)有殘余酒精</p><p><b>　　胃腔系統(tǒng)過程：</b></p><p>

71、<b>　　酒精含量，初始。</b></p><p><b>　　當(dāng)時(shí)，有</b></p><p>　　， （1）</p><p><b>　　求解得到</b></p><p><b>　　，</b></p><p>&

72、lt;b>　　當(dāng)時(shí)，有 </b></p><p>　　， （2）</p><p><b>　　解得</b></p><p><b>　　，</b></p><p>　　綜合（1）、（2）得到</p><p>　　，

73、 （3）</p><p><b>　　于是</b></p><p><b>　　。</b></p><p>　　2） 體液系統(tǒng)過程：</p><p>　　轉(zhuǎn)移因子，速率 轉(zhuǎn)移因子，速率</p><p>　　酒精含量，體液，酒精濃度，初始= ，則有

74、 </p><p><b>　　， ，</b></p><p>　　， （4）</p><p><b>　　當(dāng)時(shí)，有</b></p><p><b>　　，</b></p><p><b>　　令，，則

75、</b></p><p><b>　　， （5）</b></p><p><b>　　當(dāng)時(shí)，有</b></p><p><b>　　，</b></p><p><b>　　令，則</b></p><p><b&g

76、t;　　，</b></p><p><b>　　于是</b></p><p><b>　　，</b></p><p>　　對于短時(shí)間飲酒，體內(nèi)殘余酒精可以解釋如下：</p><p><b>　　，從而 </b></p><p><b&g

77、t;　　，</b></p><p>　　即 </p><p><b>　　。</b></p><p><b>　　4 模型求解 </b></p><p>　　上述模型的表達(dá)式(1)，(2)，(4)均可歸結(jié)為常微分方程初值問題，對于其解可用上面介紹的數(shù)值解法的方法給出。

78、在這里給出了一個(gè)計(jì)算數(shù)值解的程序。在模型中考慮長時(shí)間飲酒的情況，用MATLAB計(jì)算出，當(dāng)大李飲酒的時(shí)間達(dá)到1．865個(gè)小時(shí)，檢測時(shí)其酒精含量是20．0276毫克／百毫升，正好超標(biāo)。</p><p>　　大李短時(shí)間繼續(xù)飲酒8小時(shí)后體內(nèi)酒精含量(對于上式(6)的求解程序)</p><p>　　k1=1.8100;</p><p>　　k2=0.2100;</p&g

79、t;<p>　　DO=51200/2+O.0132;</p><p>　　c0=18.3404;</p><p><b>　　t=8;</b></p><p>　　a=(k1*D0)／(v*(kl—k2));</p><p>　　c=c0*exp(一k2*t)+a*(exp(一k2*t)一exp(一kl*t

80、))</p><p><b>　　運(yùn)行結(jié)果：</b></p><p><b>　　c=15.4695</b></p><p>　　大李長時(shí)間飲酒后體內(nèi)酒精含量(對于上式(5)的求解程序)</p><p><b>　　kl=1.8lOO</b></p><p&

81、gt;<b>　　k2=0.2lOO</b></p><p><b>　　v=447.867</b></p><p><b>　　tO=1.865</b></p><p><b>　　t=8</b></p><p><b>　　dO=0.0132

82、</b></p><p>　　f0=51200/4</p><p>　　cO=18.3404</p><p>　　a=fO/(v*k2)</p><p>　　b=(fO—k1)/(v*(kl—k2))</p><p>　　dl=f0/k1+(dO—f0/k1)*exp(一k1*t)</p>&

83、lt;p>　　cl=cO+a+b*exp(一k1*t)+(cO—a—b)*exp(一k2*t0)</p><p>　　al=d1/(y*(k1一k2))</p><p>　　c=(c1+a1)*exp(-k2*(t—tO))一a1*exp(一k1*(t—tO))</p><p><b>　　運(yùn)行結(jié)果：</b></p><

84、;p><b>　　c=20.0276</b></p><p><b>　　5 模型評價(jià)</b></p><p>　　<1>.本模型成功剖析了一部分想喝酒駕車的司機(jī)人員的心理。他們總想僥幸，然而事實(shí)不允許他們這么做，我們所做的工作讓他們的這種心理無跡可遁，對促進(jìn)交通安全也不無貢獻(xiàn)。</p><p>　　&l

85、t;2>.缺點(diǎn)：沒有考慮其他可能的因素給飲酒駕車問題帶來的影響，比如人的體重、司機(jī)的健康狀況、交警檢驗(yàn)程序不夠科學(xué)等。求得的方案也許并不是最優(yōu)的，但是相比之下比較滿意的。</p><p><b>　　6 模型推廣</b></p><p>　　嚴(yán)禁酒后駕車?現(xiàn)有動(dòng)力系統(tǒng)模型基本解決駕駛員飲酒量與停駕時(shí)間量化分析的交通難題，對駕駛員掌握駕駛時(shí)機(jī)有重要意義；模型的實(shí)際

86、應(yīng)用是當(dāng)今社會(huì)非常急需，酒后駕車者被視為公路第一殺手；應(yīng)用課題：如駕駛員飲酒量與停駕時(shí)間量化分析， 駕駛員理論培訓(xùn)．肇事時(shí)血液中酒精濃度的反推算，車保賠償?shù)鹊难芯?。我們將研究初步結(jié)果送到相關(guān)單位專家手中，聽取他們的意見。他們是本項(xiàng)目涉及到的實(shí)際應(yīng)用領(lǐng)域的執(zhí)行者和評判者。確切地說，他們的意見對我們進(jìn)一步如何完善模型是非常有積極意義的。根據(jù)他們對該研究初步結(jié)果提出的寶貴意見：</p><p>　　<1>

87、.對于酒后駕駛的安全性，保險(xiǎn)對酒后肇事的賠付等有著指導(dǎo)作用。</p><p>　　<2>.對于法醫(yī)學(xué)中所用的血中乙醇濃度反推生前飲酒量有意義。</p><p>　　<3>.實(shí)驗(yàn)嚴(yán)謹(jǐn)，結(jié)論有明顯的對比性．對于酒精在人體內(nèi)的代謝濃度，有較完整數(shù)據(jù)。</p><p>　　<4>.在“嚴(yán)禁酒后駕車”、“酒后駕車肇事不予賠償”的規(guī)定和現(xiàn)實(shí)之

88、間尋求一種合情合理又合法的新途徑，提出了“安全飲酒”的新概念。</p><p>　　<5>.“酒后安全駕車時(shí)刻表”，對于有效地預(yù)防和避免交通事故的發(fā)生有者一定的積極意義。</p><p>　　<6>.研究提供了更科學(xué)、數(shù)字化地判斷駕駛員是否應(yīng)該駕車的依據(jù)，有利于解決駕駛員飲酒量與停駕時(shí)間量化分析的交通執(zhí)法難題。對于上述寶貴意見，筆者綜合分析后，找到進(jìn)一步深化、提高、

89、拓廣研究的途徑。</p><p>　　由于模型研究中存在一個(gè)重要的假設(shè)，即酒精在各人體中的吸收、消除速率基本相同，該假設(shè)為小概率事件。事實(shí)上人對酒精的吸收與代謝的各項(xiàng)個(gè)體差異顯著，尤其是乙醇脫氫酶的個(gè)體差異非常顯著；提出以下研究設(shè)想。</p><p>　　動(dòng)力系統(tǒng)模型中考慮乙醇脫氫酶因素；探討駕駛員反應(yīng)曲線與安全駕車的關(guān)系。根據(jù)專家意見，我們對在交通執(zhí)法中的一大難題，肇事時(shí)體內(nèi)酒精濃度的反

90、推算問題方面再作了一些深入的研究。此外，由于呼吸式酒精測試儀的局限性和誤差，我們大膽提出應(yīng)用于機(jī)動(dòng)車輛的手持感應(yīng)酒精測試儀的研制設(shè)想，并提出將該模型方法用于其它手工操作的機(jī)械業(yè)，以及在醫(yī)學(xué)、農(nóng)學(xué)、養(yǎng)殖業(yè)等其他領(lǐng)域的相關(guān)研究設(shè)想，進(jìn)一步完善模型的相關(guān)內(nèi)容。</p><p><b>　　5 結(jié)論</b></p><p>　　本文從傳統(tǒng)的數(shù)值方法的簡單介紹出發(fā)，介紹了常微

91、分方程數(shù)值解法一些常用的方法，并例舉了常微分?jǐn)?shù)值解法在數(shù)學(xué)建模中的應(yīng)用。從這些應(yīng)用中我們不難看出，常微分?jǐn)?shù)值解法在實(shí)際應(yīng)用中的重要作用和意義。同時(shí)，我們應(yīng)在探索常微分?jǐn)?shù)值解法的理論研究的前提下，不斷發(fā)展相應(yīng)的軟件來為其服務(wù)。</p><p><b>　　參考文獻(xiàn)</b></p><p>　　[1] 張良勇,董曉芳. 常微分方程的起源與發(fā)展[J]. 高等函授學(xué)報(bào)(自然科

92、學(xué)版),2006(3): 34-38.</p><p>　　[2] 黃賽.常微分方程發(fā)展的主要?dú)v史沿革[J].教育與職業(yè).2006(2):36-41.</p><p>　　[3] 黃煥福.常微分方程課程建設(shè)初探[J].教育與職業(yè),2007(14):146-148.</p><p>　　[4] 林建平.常微分方程早期發(fā)展概觀[J]. 南京工程學(xué)院學(xué)報(bào),2001(2):

93、1-4.</p><p>　　[5] 錢明忠,陳友朋.常微分方程的通解[J].高等數(shù)學(xué)研究:2007(4):106-108.</p><p>　　[6] 于全祠,徐成浩.淺析常微分方程的數(shù)學(xué)思想方法[J]. 高校理科研究,2005(6)：86-87.</p><p>　　[7] 周仲旺.幾類特殊的常微分方程[J].濰坊學(xué)院學(xué)報(bào).2003(6): 7.</p&g

94、t;<p>　　[8] 邵曉鋒,徐衛(wèi)衛(wèi),李龍星.試析常微分方程模型的歸結(jié)方法[J]. 黃岡職業(yè)技術(shù)學(xué)院報(bào), 2007(1):37-39.</p><p>　　[9] 張同斌,管軍軍,焦萬堂. 基本微分方程的應(yīng)用研究[J].河南科學(xué),2007(6):888-890.</p><p>　　[10] 常廣平.常微分方程的思想方法與應(yīng)用[J].北京聯(lián)合大學(xué)學(xué)報(bào)(自然科學(xué)版), 200

95、5(60): 45-47.</p><p>　　[11] 歐陽瑞,孫要偉.常微分方程在數(shù)學(xué)建模中的應(yīng)用[J]. 宿州教育學(xué)院學(xué)報(bào),2008(2): 146 -147.</p><p>　　[12] 張學(xué)元.一類常微分方程的參數(shù)解及其應(yīng)用[J]. 湖南工程學(xué)院學(xué)報(bào),2004(4): 80-85.</p><p>　　[13] 哈迪拉·馬合買提. 論常微分方程

96、應(yīng)用題解法[J].新疆教育學(xué)院學(xué)報(bào)2003(2):41-44.</p><p>　　[14] 張善卿.一些直接代數(shù)方法的幾何解釋及應(yīng)用[J].物理學(xué)報(bào),2008(3):1335-1338. </p><p>　　[15] 王高雄,周之銘等編.?；辗址匠蘙M].高等教育出版社.1992(2).</p><p>　　[16] 王五生,歐陽云. 在常微分方程教學(xué)中滲透數(shù)

97、學(xué)建模[J]. 河池學(xué)院學(xué)報(bào),2008(2): 129- 131.</p><p>　　[17] 李大潛.傳染病動(dòng)力學(xué)的一個(gè)偏微分方程模型[J].高校應(yīng)用數(shù)學(xué)學(xué)報(bào),1986(1):17-26.</p><p>　　[18] W. F. 盧卡斯主編,朱煜民、周宇虹譯.微分方程模型[M].國防科技大學(xué)出版社,1988.</p><p>　　[19] Wolfgang W

98、alter. Ordinary Differential Equations[M].世界圖書出版公司.2003(6): 27-35.</p><p>　　[20] Howard,Irl Bivens and Stephen Davis. Calculus[M]. Eighth Edition. Beijing:高等教育出版社.2008(1):420-424.</p><p>　　[21]