多約束二階非線性常微分方程拐點的數(shù)值求法

《多約束二階非線性常微分方程拐點的數(shù)值求法》由會員上傳分享，免費在線閱讀，更多相關(guān)內(nèi)容在工程資料-天天文庫。

1、多約束二階非線性常微分方程拐點的數(shù)值求法杜志明馮長根（爆炸災(zāi)害預(yù)防、控制國家重點實驗室，北京理工大學(xué)100081）摘要：本文提出了-?種求解多約束二階非線性常微分方程拐點的數(shù)值解法,并以具有平行反應(yīng)的化學(xué)放熱系統(tǒng)為例，給出了一個具體算例。關(guān)鍵詞：非線性；常微分方程；拐點；數(shù)值計算；熱爆炸二階非線性常微分方程的求解是較為復(fù)雜的問題之一，而其極值點和拐點的確定又比求解方程的普通解要復(fù)雜得多，特別是一些強非線性的微分方程，即使數(shù)值求解也十分困難。在對科學(xué)研究和工程問題的描述過程中，包含有許多約束參數(shù)的非線性微分方程常常是不可避免的，尤其是方程的一些特殊點（如：極

2、值點、拐點等）上的解往往具有十分重要的物理意義。拐點往往對應(yīng)于物理上的轉(zhuǎn)向點、轉(zhuǎn)變點、臨界分歧點等⑴，它是實際系統(tǒng)性質(zhì)發(fā)生根本變化的關(guān)鍵點，因此，確定多約朿微分方程的拐點不僅具有科學(xué)意義，也具有實際意義。1.拐點條件拐點是數(shù)學(xué)上的一類特殊點。對于帶有約束參數(shù)的微分方程而言，微分方程的解與約朿參數(shù)密切相關(guān)。每個約束參數(shù)都有自己的取值范圍，當約束參數(shù)超出其取值范圉時，常常導(dǎo)致方程無解。而當約束參數(shù)取一些特殊值（如在其取值范闈邊界上取值）時，經(jīng)常對應(yīng)于微分方程的特殊解，如極值解、拐點處的解等。顯然，方程的解是約束參數(shù)的函數(shù)。本文中微分方程的拐點，指的就是微分方程

3、的一組特殊解，這組解剛好是解隨約束參數(shù)變化曲線上拐點處的值。下面我們尋求拐點條件。帶有多個約束參數(shù)的二階非線性常微分方程具有以下形式：y"=/（x,y,）/;Z],z2，???z”J=o（1）邊界條件為：%））心0）+03（兀0）=兒（2）式屮，兀是自變量，x0

4、；z1,z2,---,z/J=^i>{x1,77；z1,z2,-sz/J+0/S，〃；Z

5、,Z2,…,zJ—兒=0（4）方程⑴~⑶的極值條件為叫殆gSy）塲=0(5)拐點在物理上相應(yīng)于某一約束參數(shù)z”相應(yīng)于某一狀態(tài)變量的二階導(dǎo)數(shù)等于零，即dSpIdif=0,該條件相當于其中一個約束參數(shù)z°相應(yīng)于另一約束參數(shù)z”取得極值，即dzg/dzp=0o幾何圖形上拐點代表非線性曲線上極大值與極小值合二而一的點，也就是極值消失的點，根據(jù)定態(tài)方程組取得極值的條件以糾，式⑷和式⑸構(gòu)成的方程組取得拐點處解的條件為13.5J0片站%drj由于dF]/dr/=0而式(6)可簡化為

6、上式就是方程(1)~(3)的拐點條件。式⑷、(5)與式⑺聯(lián)立，即可求得原方程(1)~(3)的拐點。1.數(shù)值算法數(shù)值計算吋，首先須進行數(shù)值積分，然后用迭代法求解式(4)、(5)和(7)構(gòu)成的方程組的解。數(shù)值積分可采用Runge-Kutta-Merson法，迭代法選用Newton-Raphson法，其迭代法式為一「卅第，於)?唯境才)其中F=(Mdz“dzP沏dF^乙F=dzpdF3dzp沏為了計算rF,需定義9個輔助函數(shù):(p3=d3y/d?]2dze{屮3=d3y/d?]2dzpQ3=83y/(p=dy/dZg,02=db/d〃dZq,T]=dy/dzp

7、,屮2=d2y/d/jdz/f,Q,=3)73/7，Q2=d2y/d?j2,各輔助函數(shù)及其初值條件可由方程(1)及其初值條件分別對Zq和〃求導(dǎo)得出。3.算例本算例考慮的實際問題是內(nèi)部具有放熱化學(xué)反應(yīng)的系統(tǒng)，用微分方程描述則是帶化學(xué)反應(yīng)項的熱傳導(dǎo)方程。由于化學(xué)反應(yīng)速率與溫度呈強字線性性質(zhì)，微分方程具有非線性性質(zhì)。若考慮系統(tǒng)內(nèi)同時存在多個平行的化學(xué)放熱反應(yīng)，且反應(yīng)遵循Arrhenius定律，當系統(tǒng)達到熱平衡(化學(xué)反應(yīng)放熱速率與系統(tǒng)向環(huán)境散熱速率相等)時可用以下熱平衡方程來描述⑸+~~T~+aiexp[q&/(l+￡&)]=0(9)dp'pdp/0

8、條件可采用Frank-Kamenetskii條件⑹d0/dp=00=0,p=0Qi(10)(H)式中，&表示無量綱溫度，。為無量綱坐標，j是與形狀有關(guān)的參數(shù)，是系統(tǒng)的約束參數(shù)，各參數(shù)的定義及物理意義可參閱文獻2】?！晔窍到y(tǒng)的一個特殊參數(shù)，稱為表觀活化能。當參數(shù)￡超過某一數(shù)值時，系統(tǒng)溫度將隨參數(shù)￡連續(xù)變化，系統(tǒng)溫度不再有突變性質(zhì)，由于系統(tǒng)內(nèi)部熱積累而導(dǎo)致的燃燒爆炸現(xiàn)象將不復(fù)存在。該數(shù)值記為￡“，在熱爆炸理論中稱為轉(zhuǎn)變值，該點即是我們所要求解的拐點MJ色力,則表示系統(tǒng)內(nèi)第，個化學(xué)反應(yīng)的參數(shù)。利用上面介紹的拐點條件和數(shù)值方法，可求方程（9）~（11）的拐點，現(xiàn)將

9、數(shù)值積分吋所用到的微分方程組羅列如下：ro^jef/p=-sf[&