任務三　求解一元非線性方程

知識目標

掌握一元三次以上非線性方程求根的原理，熟悉求解一元三次以上非線性方程的方法。

能力目標

能用Newton迭代法手工求解一元三次以上非線性方程的根，能用Excel的單變量求解法求一元三次以上非線性方程的根。

化學化工中的許多問題常常歸結為解函數方程f（x）=0。若f（x）是一元線性方程式或一元二次方程式，就可用代數方法求解析解。若f（x）是一元三次或高次的代數方程式，求解析解是不可能的，只能用數值方法求近似解。例如，用于處理真實氣體的狀態方程Redlich-Kwong（RK）方程（以1mol氣體為基準）為

其展開式為

顯然，式（1-21）為體積V的一元三次非線性方程。如果要計算在一定溫度T和壓力p時氣體的摩爾體積，則需求解f（V）這樣一個一元三次非線性方程。

當然，采用數值計算法能處理上述問題，但計算過程往往較繁雜且需花費大量時間，但若以計算機為工具，利用數值計算法處理此問題，就很容易得到極為準確的數值解。本任務主要介紹處理一元非線性代數方程的常用數值計算方法及借助于計算機的具體處理步驟。

一、一元非線性方程求根方法

（一）逐步掃描法

逐步掃描法是一元非線性方程求根初始近似值常用的方法，若f（x）是n次多項式，對應的方程為n次代數方程，這時方程的根也稱為多項式的根。根有實根和復根之分，這里僅介紹實根的求法。用數值方法求方程的根可分為兩步，先找出根的某個粗略近似值，又稱為“初始近似值”，然后再將初始近似值逐步加工成滿足精度要求的結果。

設待解方程為

f（x）=0

在直角坐標系中給出相應于y=f（x）的曲線，如圖1-28所示。顯然，此曲線與x軸的交點就是方程f（x）=0的根。若函數f（x）在區間（a，b）連續，且f（a）與f（b）異號，則區間（a，b）內必定至少有一個實根。若函數f（x）在區間（a，b）連續并單調（單調上升或單調下降），則在區間（a，b）必定只有一個實根（圖1-29）。這時，選定一個步長h，并計算函數f（a）、f（a+h）、f（a）與f（a+h）兩函數值乘積，若兩函數值乘積大于零，說明該區間內無實根。這時，把區間的終點作為下一次邁步的始點。再計算f（a+2h）的值，如此反復向前邁步，直至兩函數值的乘積小于或等于零，即直至相鄰兩個函數值異號，此時可把兩函數乘積小于或等于零的區間的始點作為方程式根的近似值。這個方法叫邁步法或逐步掃描法。

顯然，步長縮小，精度提高。當精度要求較高時，步長要求很小，計算機循環計算的次數將會很大，所耗機時較大。因此，通常不用該方法求根的精確值，只用于求根的初始近似值。

用邁步法求得非線性代數方程根的初始近似值后，可用以下介紹的二分法、牛頓法和弦截法等尋求根的精確值。

☆思考：何謂逐步掃描法？

（二）二分法

設求解方程f（x）=0通過逐步掃描法已知有根區間為（x1，x2），如圖1-30所示。現取x1與x2的中點x0，即

x0=（x1+x2）/2

從而將區間（x1，x2）分為相等的兩部分。然后檢查f（x0）與f（x1）的符號是否相同，如為同號，根必在x0與x2之間，如圖1-30（a）。令

x1=x0

如為異號，則根必在x1與x0之間，如圖1-30（b）。令

x2=x0

然后再取新區間（x1，x2）的中點，重復以上步驟，直至x1與x2之間的距離小于某指定值ε為止。取前一區間的中點x0作為方程式根的精確值。這個方法稱為“二分法”或“對分法”。只要區間（x1，x2）內有根，用二分法一定能夠求出結果，因而它是最保險的一種方法，但其最大的缺點是收斂速度較慢。

需要說明的是，用上述方法求根，只能得到方程式的一個實根。實際上，若方程式有數個實根時，只需增加一個終值B，當求出第一個實根之后，檢查有根區間的終點x2是否小于終值B。若小于終值B，則再施行邁步法求下一個實根的近似值，進而求出精確值。如此反復進行下去，直至達到或超過終值B為止。

☆思考：何謂二分法？

（三）牛頓法

牛頓法亦稱牛頓-拉福森（Newton-Raphson）法。由于這個方法的計算結果頗佳，而計算過程亦較簡單，故被普遍地使用。

如圖1-31所示，假設方程f（x）=0有一個實根x*。取一初值x0，過x0作x軸垂線交于曲線f（x）于點P0，過P0點作曲線f（x）的切線并于x軸相交于x1點，顯然x1點較x0點更接近于根x*，如果|x1-x0|<ε，則方程根x*=x1，否則按上述同樣方法過x1作x軸垂線交于曲線f（x）于點P1，過P1點作曲線f（x）的切線并于x軸相交于x2……，直到|xk+1-xk|<ε為止，方程的根為

x*=xk+1

由圖1-31可得

則

于是可寫出牛頓迭代法的通式

式中　xk，xk+1——第k、k+1次求得的方程近似根。

只要函數y=f（x）的導數f'（x）易得，初值適當，用牛頓迭代法求非線性方程f（x）=0的根既快又正確。但是必須注意，起始值若不在根的附近，牛頓迭代公式不一定收斂。所以，在實際使用中，牛頓迭代法最好與逐步掃描法結合起來，先通過逐步掃描法求出根的近似值，然后再用牛頓迭代公式求其精確值，以發揮牛頓法收斂速度快的優點。

☆思考：何謂牛頓迭代法？

（四）弦截法

雖然牛頓迭代法收斂速度快，但它要求計算f'（x）的值。在科學與工程計算中，常會遇到f'（x）不易計算或算式復雜而不便計算的情況。下面介紹一種既具有牛頓迭代法收斂速度快的優點、又不用求導數f'（x）的弦截法。

弦截法的基本思想與牛頓迭代法相似，即將非線性函數f（x）線性化后求解。兩者的差別在于弦截法實現函數線性化的手段采用的是兩點間的弦線，而不是某點的切線。

如圖1-32所示，若已知非線性函數f（x）的根區間（x0，x1），過x0、x1作垂線交函數f（x）于P0、P1點，連接P0、P1交x軸于x2點。

若f（x2）=0，則方程的解為x*=x2；若f（x2）f（x1）>0，如圖1-32（a），則用x2代替x1；若f（x2）f（x1）<0，如圖1-32（b），則用x2代替x0；再用新得到的兩個點用以上方法繼續迭代，直到相鄰兩次值滿足|xk+1-xk|<ε為止，不難可得弦截法的迭代格式為

式中　xk-1，xk，xk+1——第k-1、k、k+1次求得的方程近似根。

與牛頓法只需給出一個初值不同，弦截法需要給出兩個迭代初值x0和x1。如果與逐步掃描法結合起來，則最后搜索的區間的兩個端點值常可作為初值x0和x1。

☆思考：何謂弦截法？

弦截法雖比牛頓法收斂速度稍慢，但在每次迭代中只需計算一次函數值，又不必求函數的導數，且對初值x0和x1要求不甚苛刻，因此是工程計算中常用的有效計算方法之一。

二、手工求解一元非線性方程

上面講解了幾種求解一元非線性方程的方法，下面通過實例介紹手工求解一元非線性方程的具體步驟。

【例1-8】　某一常壓氣相反應體系，當反應達平衡時其中某一組分的平衡分壓p（單位為atm）符合以下方程：

4p3-1.640p2+1.640p-0.410=0　　（1-24）

試求其分壓的數值？

解：以Newton迭代法為例

步驟1：令　　y=4p3-1.640p2+1.640p-0.410　　（1）

則　　y'=12p2-3.280p+1.640　　（2）

步驟2：確定p的初值，由題意知為常壓反應系統，則其分壓可取為p0=0.3atm。

步驟3：將p0=0.3atm代入式（1）、式（2）計算得：y0=0.04240、。

步驟4：將y0、代入式（1-22）計算得：

步驟5：計算p0、p1之間的相對誤差ε1

步驟6：判斷誤差ε1是否符合計算精度（一人為規定的兩次計算值之間的相對誤差值，手工計算時可取10-3，計算機編程計算可取10-4，甚至更小些）要求，若滿足精度要求，計算結束，此次計算值即為解，否則將本次計算值作為初值，重復步驟3～6，直到滿足計算精度要求。

本次計算的誤差ε1>10-3，需重復計算步驟3～6，具體計算結果如表1-21所示，從表中數據可知，通過三次迭代計算，其迭代精度已達10-6數量級，所以，該組分的平衡分壓p為：

p=0.2749atm

☆思考：如何利用Excel求解一元非線性方程？

三、采用Excel求解一元非線性方程

用Excel求解一元非線性方程可分為三種方法：①基于Excel單元格的Newton迭代法；②Excel自帶的單變量求解法；③規劃求解法。以下仍以例1-8為例，詳細說明采用以上三種方法如何求解非線性方程的過程。

（一）基于Excel單元格的Newton法

步驟1：將原方程、一階導數式的各系數輸入單元格，如圖1-33中單元格D3、E3、F3、G3、E4、F4、G4。

步驟2：第一次迭代取初值為0.3，如圖1-33中單元格D6，迭代精度取1×10-4，如單元格F6。建計算過程表格，如單元格A7～G7，在單元格A8～G8中對應分別輸入：“1”、“0.3”、“=$D$3*B8 3+$E$3*B8 2+$F$3*B8+$G$3”、“=$E$4*B8 2+$F$4*B8+$G$4”、“=B8-C8/D8”、“=（E8-B8）/B8”、“=IF（ABS（F8）<$F$6”，“滿足精度”，“重新計算”）。

在單元格C8中計算得函數值0.0424，單元格D8中得一階導數值1.736，單元格E8中輸入的是牛頓迭代公式（1-22）計算得到第一次迭代結果值0.275576，單元格F8中計算的是第一次迭代結果與假設的初值之間的相對誤差值-8.141×10-2，單元格G8用來判斷計算是否可以結果，采用的是Excel自帶的IF函數，其中ABS為絕對值函數。當單元格G8中判斷需“重新計算”時，需進行第二次迭代計算。

步驟3：在單元格A9中輸入“2”，B9中輸入“=E8”，即將第一次計算的結果作為第二次計算的初值，選中單元格C8～G8，利用Excel的自動填充功能向下拉，讓Excel自動填充單元格C9～G9，由圖1-33可見第二次迭代的計算結果。當G9中仍為：“重新計算”時，則需進行第三次迭代計算。

步驟4：由于第三次迭代計算方法與前兩次完全一樣，可以讓Excel自動進行。選中單元格B9～G9，利用Excel自動填充得B10～G10的結果，由圖1-33可見，第三次迭代計算結果“0.274901”與第二次計算的結果“0.274902”相對誤差已達10-6數量級，滿足計算精度要求，計算到此結束，計算結果與手工計算完全一致。

若第三次迭代計算精度仍不能滿足要求，則需進行第四、五……，方法與步驟4完全相同，直到計算滿足精度為止。

注意：此例中為便于使用Excel單元格的自動填充功能，單元中的輸入格式有兩種：絕對引用格式和相對引用格式，這兩種格式的詳細區別請參閱Excel相關章節的內容。

（二）Excel自帶的單變量求解法

步驟1：打開Excel，選擇某一單元格，如圖1-34中的單元格B14，輸入：

“=4*A14 3-1.64*A14 2+1.64*A14-0.41”。

此時由于單元A14中沒有任何信息，Excel將其視為“0”，輸入完畢確定，在單元格B14中將得到原函數常數項數值，即“-0.41”。

單元格A14相當于原表達式中的分壓p，若在此單元格中輸入某一數據可使單元格B14中的值為0，則此數據為該函數的一個解，單變量求解就據此原理。

步驟2：選中單元格B14，如圖1-34，單擊“數據”?“假設分析”?“單變量求解”，出現“單變量求解”對話框（見圖1-35），在各空格中依次填入以下信息。

目標單元格：B14；

目標值：0；

可變單元格：$A$14或A14。

步驟3：在圖1-35中按“確定”按鈕后?圖1-36“單變量求解狀態”對話框。單擊確定按鈕，答案出現在單元格A14中（如圖1-35所示）：0.27490124，很顯然，此值與前面計算的結果完全一致。

需要說明的是單元格B14中的值為0時，說明此根為一精確解。但有時單元格B14中是一個很小的數值，Excel可能顯示的值為“0”而并非真為0。若要減小誤差使解更精確，可采取以下方法操作：單擊左上角“Excel按鈕”?“Excel選項”?“公式”?出現“Excel選項”卡，如圖1-37所示，可將“計算選項”中的“最多迭代次數”和“最大誤差”值進行修改，如將最多迭代次數由“100”改為“10000”，此值最大值為32767。最小誤差由“0.001”改為“0.0001”等，這樣可使計算精度得到提高。

注意，此處的最小誤差不能改得過小，否則有可能導致超出最大循環次數而無法計算。

☆思考：如何利用ExcelVBA編寫一段求解一元三次非線性方程的小程序？

（三）規劃求解法

步驟1：打開Excel，選擇某一單元格如圖1-38中的單元格B20，按式（1-24）輸入前三項：“=4*A20 3-1.64*A20 2+1.64*A20”。此時由于單元A20中沒有任何信息，Excel將其視為“0”，輸入完畢確定，在單元格B20中將顯示“0”。

步驟2：選中單元格B20，單擊“數據”，選擇“規劃求解”，調出“規劃求解參數”對話框，如圖1-39所示。選擇“值為（ <文字顏色>V</文字顏色> ）”，輸入方程的常數項值“0.41”。可變單元格輸入框中輸入：$A$20。單擊“求解”按鈕?圖1-40的“規劃求解結果”對話框。

步驟3：單擊圖1-40中的“確定”按鈕?規劃求解結果，如圖中單元格A20中的數據“0.27490143”，可見計算結果與前面兩種方法完全一致。

注意：在Office軟件安裝時，規劃求解往往不會被安裝，可以通過加載實現。步驟是：單擊軟件左上角“Office按鈕”?“Excel選項”?“加載項”?選擇“規劃求解”加載項?“加載宏”對話框?在“規劃求解”項打勾即可。

四、技能拓展——ExcelVBA自定義函數求解一元三次非線性方程

應用牛頓迭代法計算原理，采用VBA自編迭代函數求解，可使求解過程大為簡化，以下仍以例1-8為例說明具體步驟。

步驟1：打開Excel?“開發工具”?“VisualBasic”編輯器?插入?模塊?“添加過程”對話框，如處理的是一元三次非線性方程，則可輸入函數名“Newton3”，如圖1-41所示。

步驟2：選擇好相關類型后，按“確定”，出現VBA代碼編輯窗口，在編輯窗口編寫VBA代碼，如圖1-42所示。

步驟3：在Excel表格中分別輸入方程初值、迭代精度及方程中各系數之值，如圖1-43中單元格D6、F6、D3、E3、F3、G3。

步驟4：在圖1-43單元格A17中輸入：=Newton3（D6，F6，D3，E3，F3，G3），按“Enter”即得方程的解為：0.27490124。

顯然，采用自編VBA函數求解結果與單變量求解結果完全相同。