第2章非線性方程求解

第2章非線性方程求解

2.1化工實際問題的提出2.2實根的對分法2.3直接迭代法2.4松弛迭代法2.5韋格斯坦法2.6牛頓迭代法2.7割線法2.8非線性方程組的牛頓方法2.9化工生產中非線性方程組求解應用實例

目錄2.1化工實際問題的提出

求解非線性方程是化工設計及模擬計算中必須解決的一個問題。與線性方程相比，非線性方程問題無論是從理論上還是從計算公式上，都要比線性方程復雜的多。對于一般的非線性f(x)=0，計算方程的根既無一定章程可循也無直接方法可言。例如，求解高次方程組7x6-x3+x-1.5=0的根，求解含有指數和正弦函數的超越方程ex-sin(x)=0的零點。解非線性方程或非線性方程組也是計算方法中的一個主題。一般地，我們用符號f(x)來表示方程左端的函數，方程的一般形式表示為f(x)=0，方程的解稱為方程的根或函數的零點。通常，非線性方程的根不止一個，而任何一種方法只能算出一個根。因此，在求解非線性方程時，要給定初始值或求解范圍。而對于具體的化工問題，初值和求解范圍常?？筛鶕唧w的化工知識來決定。常見的雷諾數和摩擦系數關系方程在雷諾數低于4000時有以下關系式：

2.12.82.72.52.32.2總目錄2.92.62.42.1化工實際問題的提出這是一個典型的非線性方程。我們在管路設計中經常碰到。當我們已知雷諾數Re，如何根據公式（2-1）求出摩擦系數λ，這是我們在管路設計中必須首先解決的問題。對于方程（2-1）而言，無法用解析的方法求出摩擦系數，只能用數值求解的方法。如用在下面即將介紹的松弛迭代法，假設：

則利用松弛迭代公式可得：經11次迭代可得摩擦系數為0.07593。同樣，在n個組分的等溫閃蒸計算中，通過物料和相平衡計算，我們可得到如下非線性方程：

2.12.82.72.52.32.2總目錄2.92.62.42.1化工實際問題的提出在方程（2-3）中只有α是未知數，ki為相平衡常數，zi為進料組分的摩爾濃度，均為已知數。和上面的情況一樣，方程（2-3）也無法直接解析求解，必須利用數值的方法，借助于計算機方可精確的計算。對于這個問題的求解，可利用我們下面介紹的牛頓迭代法進行計算，也可利用其他迭代公式進行計算，如采用牛頓迭代公式，則可以得到如下的具體迭代公式：

（2-4）

飽和蒸氣壓是我們經常要用到的數據，雖然我們可以通過實驗測量來獲取飽和蒸氣壓的數據，但我們通常利用前人已經測量得到的數據或回歸的公式來獲取，這可以減輕我們大量的基礎實驗工作。公式（2-5）是一種常用的飽和蒸氣壓計算公式：其中p為飽和蒸氣壓，單位為mmHg，T為溫度，單位為K，A、B、C、D為已知系數。要想得到某一溫度下的飽和蒸氣壓，直接利用公式（2-5）是無法得到的。因為公式（2-5）兩邊都有未知變量，并且無法用解析的方法求解，必須用數值計算的方法求解。通過上面的一些例子，我們可以發(fā)現，如果沒有適當的手段和辦法來求解非線性方程，那么化學化工中的許多研究、設計等工作將無法展開，這勢必影響化學化工的發(fā)展，下面我們將介紹一些實用的非線性方程求解方法，并提供計算機程序。

2.12.82.72.52.32.2總目錄2.92.62.42.2實根的對分法

2.2.1使用對分法的條件

2.2.2對分法求根算法

2.2.3對分法VB程序清單

2.12.82.72.52.32.2總目錄2.92.62.42.2.1使用對分法的條件對分法或稱二分法是求方程近似解的一種簡單直觀的方法。設函數f(x)在[a,b]上連續(xù)，且f(a)f(b)<0，則f(x)在[a,b]上至少有一零點，這是微積分中的介值定理，也是使用對分法的前提條件。計算中通過對分區(qū)間，逐步縮小區(qū)間范圍的步驟搜索零點的位置。如果我們所要求解的方程從物理意義上來講確實存在實根，但又不滿足f(a)f(b)<0，這時，我們必須通過改變a和b的值來滿足二分法的應用條件。

2.12.82.72.52.32.2總目錄2.92.62.42.2.2對分法求根算法

計算f(x)=0的一般計算步驟如下：

1、輸入求根區(qū)間[a,b]和誤差控制量ε，定義函數f(x)。

2、判斷:如果f(a)f(b)<0則轉下，否則，重新輸入a和b的值。

3、計算中點x=(a+b)/2以及f(x)的值分情況處理（1）|f(x)|<ε：停止計算x*=x，轉向步驟4（2）f(a)f(x)<0：修正區(qū)間[a,x]→[a,b]，重復3（3）f(x)f(b)<0：修正區(qū)間[x,b]→[a,b]，重復34、輸出近似根x*。

右圖給出對分法的示意圖。

x3=(x0+x2)/2

x2=

(x0+x1)/2

x0

x3

x1

x1

2.12.82.72.52.32.2總目錄2.92.62.42.2.3對分法VB程序清單

PrivateSubCommand1_Click()Dimx1,x2,x,y1,y2,y,eer80x1=InputBox("x1")x2=InputBox("x2")y1=f(x1)y2=f(x2)Ify1*y2<0Then

GoTo100ElsePrint"pleaserepeatinputx1andx2"

GoTo80EndIf100x=(x1+x2)/2y=f(x)IfAbs(y)<=0.001ThenPrint"thefunctionrootis";xPrint"y=";yElseIfy1*y<0Thenx2=xy2=yGoTo100Elsex1=xy1=yGoTo100EndIfEndIfEndSubPublicFunctionf(x)Dimyy=x^3+x^2-1f=yEndFunction2.12.82.72.52.32.2總目錄2.92.62.42.2.3對分法求解實例用對分法求在區(qū)間[1,2]之間的根。解：(1)f(1)=-2.8，f(2)=0.3，由介值定理可得有根區(qū)間[a,b]=[1,2]。

(2)計算x2=(1+2)/2=1.5，f(1.5)=-0.45，有根區(qū)間[a,b]=[1.5,2]。

(3)計算x3=(1.5+2)/2=1.75，f(1.75)=0.078125，有根區(qū)間[a,b]=[1.5,1.75]。

一直做到|f(xn)|<ε(計算前給定的精度)或|a-b|<ε時停止。詳細計算結果見表2-1。對分法的算法簡單，然而，若f(x)在[a,b]上有幾個零點時，如不作特殊處理只能算出其中一個零點；另一方面，即使f(x)在[a,b]上有零點，也未必有f(a)f(b)<0。這就限制了對分法的使用范圍。對分法只能計算方程f(x)=0的實根。

對于多個零點的方程，我們可以通過將給定的區(qū)間[a,b]進行細分，然后在細分后的區(qū)間內用二分法分別求解，從而得到多個零點。例如求方程在0-30內的所有根。需要對二分法進行以下處理：即先給定一個a，本例中為0，然后不斷增加，直到找到一個b，使f(a)f(b)<0，調用二分法，計算在[a，b]范圍內的根，然后將b作為a，重復上面的工作，直到計算范圍超出30為止。

2.12.82.72.52.32.2總目錄2.92.62.4VB調用2.3直接迭代法

對給定的方程f(x)=0，將它轉換成等價形式：。給定初值x0，由此來構造迭代序列，k=1,2,…，如果迭代收斂，即有,則就是方程f(x)=0的根。在計算中當小于給定的精度控制量時，取為方程的根。例如，代數方程x3-2x-10=0的三種等價形式及其迭代格式如下：對于方程構造的多種迭代格式，怎樣判斷構造的迭代格式是否收斂？收斂是否與迭代的初值有關？根據數學知識，我們可以直接利用以下收斂條件：

1、

當有

2、

在[a,b]上可導，并且存在正數L<1，使任意的，有則在[a,b]上有唯一的點滿足，稱為的不動點。而且迭代格式對任意初值均收斂于的不動點，并有下面誤差估計式：

(2-6)

要構造滿足收斂條件的等價形式一般比較困難。事實上，如果為f(x)的零點，若能構造等價形式，而，由的邊疆性，一定存在的鄰域，其上有，這時若初值迭代也就收斂了。由此構造收斂迭代格式有兩個要素，其一，等價形式應滿足；其二，初值必須取自的充分小鄰域，其大小決定于函數f(x)，及做出的等價形式。2.12.82.72.52.32.2總目錄2.92.62.42.3直接迭代法例：求代數方程x3-2x-5=0，在x0=2附近的零點。

解：1）x3=2x+5構造的迭代序列收斂。取x0=2，則

準確的解是x=2.09455148150。2）將迭代格式寫為

迭代格式不能保證收斂,但并不一定不收斂。VB程序界面：

2.12.82.72.52.32.2總目錄2.92.62.42.4松弛迭代法

有些非線性方程或方程組當用上一節(jié)中的直接迭代法求解時,迭代過程是發(fā)散的。這時可引入松弛因子,利用松弛迭代法。通過選擇合適的松弛因子,就可以使迭代過程收斂。松弛法的迭代公式如下:

（2-7）由上式可知，當松弛因子ω等于1時，松弛迭代變?yōu)橹苯拥?。當松弛因子ω大?時松弛法使迭代步長加大,可加速迭代，但有可能使原來收斂的迭代變成發(fā)散。當0<ω<1時,松弛法使迭代步長減小,這適合于迭代發(fā)散或振蕩收斂的情況,可使振蕩收斂過程加速。當ω<0時,將使迭代反方向進行,可使一些迭代發(fā)散過程收斂。松弛法是否有效的關鍵因子是松弛因子ω的值能否正確選定。如果ω值選用適當,能使迭代過程加速,或使原來不收斂的過程變成收斂。但如果ω值選用不合適，則效果相反,有時甚至會使原來收斂的過程變得不收斂。松弛因子的數值往往要根據經驗選定,但選用較小的松弛因子,一般可以保證迭代過程的收斂。2.12.82.72.52.32.2總目錄2.92.62.42.4松弛迭代法例:用松弛迭代法求解下面非線性方程組,并分析松弛因子對迭代次數及收斂過程的影響。已知迭代初值x和y均為0,收斂精度ε=0.001

解:取以下迭代表達式:

若取松弛因子為1.1則其迭代過程如左表。若改變松弛因子,迭代過程及迭代所需的次數亦將發(fā)生變化，詳見右表。由右表數據可知,當松弛因子小于1時,增大松弛因子,可加速迭代過程,減少迭代次數,但當松弛因子大于1時,迭代次數反而增加,當松弛因子達到1.56

時,迭代過程分散。

2.12.82.72.52.32.2總目錄2.92.62.4此法是一種迭代加速方法。在圖2-3中,從x0開始,曲線y=φ(x)和直線y=x之間的階梯形折線為單變量的直接迭代過程。其迭代步長為每個階梯的水平距離。若利用圖2-3中曲線φ(x)上的兩個點作一直線,通過求該直線和y=x的交點C的橫坐標xW作為新的迭代點,這樣的迭代計算過程就是韋格斯坦法。即韋格斯坦法需要在已知兩點的前提下迭代計算第三點。一般第一點為人為設定,第二點利用直接迭代計算而得,則第三點就可以用韋格斯坦法迭代求解。已知A點的坐標為(x1,y1),B點的坐標為(x2,y2)。則過A、B兩點的直線方程為：求該直線和y=x的交點可得:

韋格斯坦法的一般計算通式為

2.5韋格斯坦法

A

y1

y0

φ(x)Bx0x1x2x*

xw圖2.-3韋格斯坦法迭代示意圖2.12.82.72.52.32.2總目錄2.92.62.42.5韋格斯坦法由上述公式可知,韋格斯坦法也是一種松弛法,其松弛因子為

一般情況下,當1>k>0時,迭代過程為單調收斂過程。當-1<k<0時,迭代過程為振蕩收斂過程,但當k=1時,收斂將發(fā)散,故在編程計算時應注意當k=1時則取k=0進行計算。韋格斯坦法求解方程x3-2x2+5x-4=0根的QB程序見課本2.12.82.72.52.32.2總目錄2.92.62.42.6牛頓迭代法

2.6.1牛頓法的理論推導2.6.2牛頓法的幾何意義2.12.82.72.52.32.2總目錄2.92.62.4對方程f(x)=0可構造多種迭代格式，牛頓迭代法是借助于對函數f(x)=0的泰勒展開而得到的一種迭代格式。將f(x)=0在初始值x0做泰勒展開得：取展開式的線性部分作為的近似值，則有：設則令類似地，再將f(x)=0在x1作泰勒展開并取其線性部分得到：一直做下去得到牛頓法的迭代格式：牛頓迭代格式對應于f(x)=0的等價方程為：若b是f(x)的單根時，，則有，只要初值x0充分接近b，牛頓迭代都收斂。牛頓迭代是二階迭代方法?？梢宰C明，b為f(x)的a重根時，迭代也收斂，但這是一階迭代，收斂因子為，若這時取下面迭代格式，它仍是二階方法：

2.6.1牛頓法的理論推導

2.12.82.72.52.32.2總目錄2.92.62.4以為斜率作過(x0，f(x0))點的直線，即作f(x)在x0的切線方程：令y=0，則在x1處的切線與x軸的交點x1，即：再作f(x)在x1處的切線，得交點x2，逐步逼近方程的根b。如圖2-4所示。在區(qū)域[x0,x0+h]的局部“以直代曲”是處理非線性問題的常用手法。在泰勒展開中，截取函數展開的線性部分替代f(x)。

yxx2x1x02.6.2牛頓法的幾何意義

圖2-4牛頓切線法示意圖

2.12.82.72.52.32.2總目錄2.92.62.42.6.2牛頓法的幾何意義例：用牛頓迭代法求方程f(x)=x3-7.7x2+19.2x-15.3，在x0=1附近的零點。解：

計算結果列于表2-4中。

比較表2-1和表2-4的數值，可以看到牛頓迭代法的收斂速度明顯快于對分法。牛頓迭代法也有局限性。在牛頓迭代法中，選取適當迭代初始值x0是求解的前提，當迭代的初始值x0在某根的附近時迭代才能收斂到這個根，有時會發(fā)生從一個根附近跳向另一個根附近的情況，尤其在導數數值很小時，如圖2-5所示。x2yx1xx0x2x0x3x1如果f(x)=0沒有實根，初始值x0是實數，則迭代序列不收斂。圖2-6給出迭代函數f(x)=2+x2，初始值x0=2的發(fā)散的迭代序列。圖2-5圖2-62.12.82.72.52.32.2總目錄2.92.62.42.7割線法

在牛頓迭代格式中：

，用差商代導數，并給定初始值x0和x1

，那么迭代格式可寫成如下形式：上式稱為割線法。用割線法迭代求根，每次只需計算一次函數值，而用牛頓迭代法每次要計算一次函數值和一次導數值。但割線收斂速度稍慢于牛頓迭代法，割線法為1.618階迭代方法。做過兩點(x0,f(x0))和(x1,f(x1))的一條直線（弦），該直線與x軸交點就是生成的迭代點x2，再做過(x1,f(x1))和(x2,f(x2))的一條直線，x3是該直線與x軸的交點，繼續(xù)做下去得到方程的根f(a)=0，如圖2.4所示。例2.5：用割線法求方程的根，取x0=1.5，x1=4.0。解：

計算結果列于表

VB程序2.12.82.72.52.32.2總目錄2.92.62.42.8非線性方程組的牛頓方法

為了敘述的簡單，我們以解二階非線性方程組為例演示解題的方法和步驟，類似地可以得到解更高階非線性方程組的方法和步驟。設二階方程組

其中x,y為自變量。為了方便起見，將方程組寫成向量形式：將在（x0,y0）附近進行二元泰勒展開，并取其線性部分，得到下面方程組：令則有

2.12.82.72.52.32.2總目錄2.92.62.42.8非線性方程組的牛頓方法如果再將原方程組在u1處進行二元泰勒展開，并取其線性部分，得到下面方程組：解出得出

繼續(xù)做下去，每一次迭代都是一個方程組

2.12.82.72.52.32.2總目錄2.92.62.4即為止。2.8非線性方程組的牛頓方法例2.6：求解下面非線性方程組取初始值解：解方程得繼續(xù)做下去，直到時停止。

2.12.82.72.52.32.2總目錄2.92.62.42.9化工生產中非線性方程組求解應用實例

在化工生產中,為了求解反應前后各物料的濃度,常常要聯立求解一些非線性方程組,這些方程組難以用常規(guī)的解析方法求解,一般只能利用數值求解的方法加以求解。下面是在合成氨生產中利用非線性方程組求解方法求解烴類蒸氣轉化反應前后各物料濃度的實例。例2.7

在合成氨生產中,