【正文】 酒精從系統(tǒng) I向系統(tǒng) II的轉(zhuǎn)移的速 率系數(shù)，及向體外的排出的速率系數(shù)，與該系統(tǒng)的酒精濃度成正比，這兩個(gè)速率系數(shù) 1k 、 2k 是由人體的身體機(jī)能所決定的常數(shù)。具體程序如下： 設(shè)方程（ ）中的 2 2 1( ) , , 2 , ( ) .N t y k k N t b??? ? ? ?于是有下面程序： s dso l v e ( D y a * b * K 2 a * b * y a * K 2 * y a * ^ 2y?? ? ? ? ?,y(0)=0`) S= (K2*exp(t*a*b+t*a*K2+log(b/K2)/(b+K2)*b+log(b/K2)/(b+K2)*K2)b)/ (1+exp(t*a*b+t*a*K2+log(b/K2)/(b+K2)*b+log(b/K2)/(b+K2)*K2)) 司機(jī)飲酒駕車防避模型的數(shù)值解法 在 2021年全國(guó)大學(xué)生數(shù)學(xué)建模競(jìng)賽題中有一個(gè)關(guān)于司機(jī)飲酒駕車模型。本節(jié)經(jīng)過(guò)細(xì)致的分析，建立數(shù)學(xué)模型，對(duì)這一現(xiàn)象做出了科學(xué)的解釋。具體地說(shuō)，對(duì)一階線性微分方程組 ( ) ( 8 . 5 )dy A y xdx ? ? ? 其中 y ， mR?? , A為 m 階方陣。當(dāng)在區(qū)間上使用的步長(zhǎng)太大時(shí)，會(huì)出現(xiàn)不穩(wěn)定，有時(shí)表現(xiàn)為計(jì)算解的振蕩性。 對(duì)于多步法，常用的預(yù)估 — 校正公式有 MilneHamning 方法、 MilneSimpon方法、 4 階隱預(yù)估 — 校正公式等。 階龍格 庫(kù)塔（ Runge— Kutta）方法公式 要進(jìn)一步提高精度，必須取更多的點(diǎn)，如取 4 點(diǎn)構(gòu)造如下形式的公式： 常微分方程的高精度求解方法 安徽大學(xué)江淮學(xué)院 07 計(jì)算機(jī) (1)班 15 51 1 1 2 2 3 3 4 412 1 1 13 2 2 1 3 24 3 4 1 2 6 3()( , )( , ) ( 6 .5 )( , )( , )nnnnnnnnnny y h k k k kk f x yk f x h y h kk f x h y h k h kk f x h y h k h k h k? ? ? ???? ? ?? ? ? ??? ? ? ? ? ????? ? ? ??? ? ? ? ??? ? ? ? ? ?? 其中待定系數(shù) i? ， i? ， i? 共 13 個(gè)，經(jīng)過(guò)與推導(dǎo) 2 階龍格 — 庫(kù)塔公式類似、但更復(fù)雜的計(jì)算，得到使局部誤差 511( ) ( )nny x y o h????的 11 個(gè)方程。 如上分析啟示我們，在區(qū)間 ? ?1,nnxx? 內(nèi)多取幾個(gè)點(diǎn)，將它們的斜率加權(quán)平均作為 K ，就有可能構(gòu)造出精度更高的計(jì)算公式。如果 N 固定，則理論上可以推導(dǎo)出步長(zhǎng) h ，使之滿足任意想要的最終全局誤差。對(duì) 1()yt 求解方程（ ） 1010( ) ( ) ( , ( ) )tty t y t f t y t d t?? ? （ ） 然后可用數(shù)值積分方法逼近（ ）中的定積分，如果采用步長(zhǎng) 10h t t?? 的梯形公式，則結(jié)果為 1 0 0 0 1 1( ) ( ) ( ( , ( ) ) ( , ( ) ) )2hy t y t f t y t f t y t? ? ? （ ） 注意公式（ ）的右端包含了待定值 1()yt 。 梯形公式也是隱式格式，一般需用迭代法求解，迭代公式為 ( 0 )1( 1 ) ( )1 1 1( , )( , ) ( , ) ( )2( 0 , 1 , 2 , )n n n nkkn n n n n ny y hf x yhy y f x y f x yk??? ? ?? ??????? ? ?? ?????? 由于函數(shù) ( , )f xy 關(guān)于 y滿足 Lipschitz條件，容易看出 ( 1 ) ( ) ( ) ( 1 )1 1 1 12k k k kn n n nhLy y y y??? ? ? ?? ? ? 其中 L為 Lipschitz常數(shù)。向后 Euler公式的右端含有 1ny? ，因此是隱式公式，一般要用迭代法求解，迭代公式通常為 常微分方程的高精度求解方法 安徽大學(xué)江淮學(xué)院 07 計(jì)算機(jī) (1)班 8 ( 0 )1( 1 ) ( )1 1 1( , ) ( 3 . 2 )( , ) ( 0 , 1 , 2 , )n n n nkkn n n ny y h f x yy y h f x y k? ?? ? ?? ????? ? ??? Euler方法的誤差估計(jì) 對(duì)于向前 Euler 公式（ ）我們看到，當(dāng) 1,2,n? 時(shí)公式右端的 ny 都是近似的，所以用它計(jì)算的 1ny? 會(huì)有累積誤差，分析累積誤差比較復(fù)雜，這里先討論比較簡(jiǎn)單的，所謂局部截?cái)嗾`差。 getchar()。 for(i=0。 for (i=0。 } /********************* changing part end**************************/ void cal_error() //計(jì)算誤差值以對(duì)各種方法進(jìn)行比較 { int i。 （ iii） Taylor 多項(xiàng)式近似 將函數(shù) ()yx 在 nx 處展開(kāi)，取一次 Taylor 多項(xiàng)式近似，則得 1( ) ( ) ( ) ( , ( , ( ) )n n n n n ny x y x h y x h f x y x y x? ?? ? ? 再將 ()nyx 的近似值 ny 代入上式右端，所得結(jié)果作為 1()nyx? 的近似值 1ny? ，得到離散化的計(jì)算公式 1 ( , )n n n ny y hf x y? ?? 以上三種方法都是將微分方程離散化的常用方法，每一類方法又可導(dǎo)出不同形式的計(jì)算公式。 常微分方程的高精度求解方法 安徽大學(xué)江淮學(xué)院 07 計(jì)算機(jī) (1)班 4 常溫分方程的離散化 下面主要討論 一階常微分方程的初值問(wèn)題，其一般形式是 0( , ) ( )()dy f x y a x bdxy a y? ? ? ???? ?? 在下面的討論 我們總假定函數(shù) f (x, y) 連續(xù)，且關(guān)于 y 滿足李普希茲(Lipschitz)條 件，即存在常數(shù) L ，使得 ( , ) ( , )f x y f x y L y y? ? ? 這樣，由常微分方程理論知，初值問(wèn)題 (1)的解必定存在唯一。y39。據(jù)此我們可找到溫度與時(shí)間之間的函數(shù)關(guān)系式，這事實(shí)上就是一個(gè)微分方程的建立問(wèn)題。 Huon method。文章首先 案例引入微分方程概念，然后 給出了微分方程的基本概念。最后，再討論高階常微分方程和一階常微分方程組：一般的高階常微分方程都可以通過(guò)相應(yīng)的變量代換轉(zhuǎn)化為一階常微分方程組，一階常微分方程的初值問(wèn)題求數(shù)值解與一階常微分方程的初值問(wèn)題求數(shù)值解的方法基本相同。 案例：一次謀殺案，在某天下午四點(diǎn)發(fā)現(xiàn)尸體，尸體的體溫為 30℃，假設(shè)當(dāng)時(shí)屋內(nèi)空間的溫度保護(hù) 20℃不變，現(xiàn)判斷謀殺是何時(shí)發(fā)生的？ 解決此問(wèn)題首先必須要從尸體溫度的變化尋求關(guān)系式，這就需要知道物理學(xué)中的加熱與冷卻規(guī)律。yxy39。 1 e e 39。 需要說(shuō)明的是，用不同的差商近似導(dǎo)數(shù)，將得到不同的計(jì)算公式。 double X [MAX], Y [MAX], Z[MAX], F[MAX], G[MAX], CY[MAX], E[MAX]。 } return。 for(i=0。 for(i=0。 如果在微分方程離散化時(shí)，用向后差商代替導(dǎo)數(shù)，即 11 ( ) ( )() nnn y x y xyx h?? ?? ?，則得計(jì)算公式 1 1 10( , ) ( 0 , 1 , ) ( 3 . 1 )()n n n ny y h f x y ny y a? ? ?? ? ??? ?? 用這組公式求問(wèn)題（ ）的數(shù)值解稱為向后 Euler 公式。 改進(jìn)的 Euler 方法 梯形公式 利用數(shù)值積分方法將微分方程離散化時(shí)，若用梯形公式計(jì)算式 ()中之右端積分，即 ? ?111( , ( ) ) ( , ( ) ) ( , ( ) )2nnx n n n nx hf x y x d x f x y x f x y x? ????? 常微分方程的高精度求解方法 安徽大學(xué)江淮學(xué)院 07 計(jì)算機(jī) (1)班 9 這就是求解初值問(wèn)題（ ）的梯形公式。 為便于編制程序上機(jī)，式（ ）常改寫(xiě)成 1( , )( , ) ( )1 ()2p n n nq n n pn p qy y hf x yy y hf x h yy y y??? ????? ? ???? ???? 常微分方程的高精度求解方法 安徽大學(xué)江淮學(xué)院 07 計(jì)算機(jī) (1)班 10 改進(jìn) Euler 法是二階方法。求導(dǎo)公式可以遞歸得計(jì)算： ()y t f? ? t y t yy f f y f f f?? ?? ? ? ? 22 ( )tt ty y y yy f f y f y f y??? ? ?? ?? ? ? ? 22 ( )tt ty y y y tf f f f f f f f f? ? ? ? ? ( 4 ) 23 3 ( ) 3tt t tt y ty y tyy f f y f y f y? ? ??? ? ? ? 33 ( )y yy yyyf y f y y f y??? ? ?? ?? ? ? 2 3 2( 3 3 ) ( 2 )tt t tt y ty y y y y y tt ty yf f f f f f f f f f f f f? ? ? ? ? ? ? 23 ( ) ( ) ( )t y ty y y y t yf f f f f f f f f f? ? ? ? ? () 并且一般有， ( ) ( 1 )( ) ( , ( ) )NNy t P f t y t?? （ ） 其中 P 為導(dǎo)數(shù)算子 ()Pfty???? 區(qū)間 0[ , ]Ntt 上的初值 ( ) ( , )y t f t y? ? 的近似數(shù)值解可由各子區(qū)間 1[ , ]kktt? 上的公式（ ）來(lái)推導(dǎo)?？墒墙o出一種斜率 K ,()式就對(duì)應(yīng)地導(dǎo)出一種算法。由于（ ）式有 4 個(gè)未知數(shù)而只有 3 個(gè)方程，所以解不唯一。 多步 法的一個(gè)優(yōu)點(diǎn)是，可以確定它的局部截?cái)嗾`差，并可以包含一個(gè)校正項(xiàng)， 用于在 每 一步計(jì)算中改變解的精確度。 正確的步長(zhǎng) 用預(yù)估 — 校正方法在 大區(qū)間上求解初值問(wèn)題 00( , ), ( )y f t y y t y? ??時(shí)，有時(shí)也出現(xiàn)問(wèn)題。 高階微分方程的數(shù)值解法 高階微分方程的初值問(wèn)題可以通過(guò)變量代換化為一階微分方程組初值問(wèn)題。據(jù)此在時(shí)間一銷售坐標(biāo)系給出的曲線稱為產(chǎn)品的生命曲線，其形狀呈鐘型。 其中 1()Nt， ? 表示外部信息使“創(chuàng)新型”顧客購(gòu)買新產(chǎn)品的比率；112 K?表示口傳信息使“創(chuàng)新型”顧客購(gòu)買新產(chǎn)品的比率； 12K? 表示口傳信息使“模仿型’’顧客購(gòu)買新產(chǎn)品的比率。 參考數(shù)據(jù) 65％至 70％，其中血液只占體重的 7％左右：而藥物 (包括酒精 )在血液中的含量與在體液中的含量大體是一樣的。 符號(hào)說(shuō)明： 0f ：酒精進(jìn)入胃腔的速率， 設(shè)為常數(shù) t ：測(cè)試時(shí)間 (小時(shí) ) 0t ：飲酒時(shí)間 (小時(shí) ) 1()xt： t 時(shí)刻人體胃腔中的酒精含量 (毫克／百毫。 3．怎樣估計(jì)血液中 的酒精含量在什么時(shí)間最高。綜上，斯蒂芬斯一莫賽模型是一常微分方程組的初值問(wèn)題模型： ? ? ? ?? ? ? ?? ?11 1 122 2 1 212() ( ) ( ) ,() ( ) ( ) ( ) , ( )( 0) 0 , 0 0dN t K N t N tdtdN t K N t N t N tdtNN???? ? ? ???? ? ? ???? ????