非線性方程組重根的可信驗證方法

桑海風，萬保成

（1.吉林大學 數學學院，長春130012；2.北華大學 數學與統計學院，吉林 吉林132013；3.吉林農業大學 信息技術學院，長春130118）

在火箭噴口受力、核磁共振機設計和數碼機床的控制等高風險應用領域，計算結果的可靠性至關重要.Rump［1］給出了標準的可信驗證方法，該方法將浮點運算用于嚴格證明，解決了非奇異問題的驗證.如判斷非線性方程組單根的存在性與唯一性問題，利用可信驗證方法可以得到該問題的一個近似解及其相應的誤差界，使得在近似解的誤差界范圍內必存在一個精確解.但驗證非線性方程組多重根的存在性是一個奇異問題，因為對非線性方程組的系數做微小擾動，一個孤立奇異解重數就可能發生改變，而且在驗證非線性方程組是否具有重根的過程中存在舍入誤差，故驗證非常困難.因此，對于非線性方程組多重根的驗證，首先要將該方程組正則化，得到一個新方程組，使原來方程組的多重根成為新方程組單根的一部分，然后再利用標準可信驗證方法驗證新方程組的單根.Rump等［2］給出了Jacobi矩陣秩虧為1的非線性方程組二重根驗證方法，該方法通過將原方程組增加一個光滑變量，并增加一個含n－1個變元、n個方程的方程組，得到一個含2n個變元、2n個方程的新方程組，使得新方程組的Jacobi矩陣非奇異.將驗證原方程組Jacobi矩陣秩虧為1的奇異問題轉化為驗證新方程組Jacobi矩陣滿秩的非奇異問題.Li等［3］給出了Jacobi矩陣秩虧為1的多重根的驗證方法，該方法通過將原方程組增加μ－1個光滑變量，并增加μ－1個方程，得到一個含n＋μ－1個變元及n＋μ－1個方程的新系統，將驗證原方程組Jacobi矩陣秩虧為1的奇異問題轉化為驗證新方程組Jacobi矩陣滿秩的非奇異問題.

本文利用Shen等［4］的數值算法及區間算法，研究Jacobi矩陣秩虧為q的非線性方程組重根的可信性驗證方法.該方法通過將原方程組增加q個光滑變量，并增加q個方程，得到一個Jacobi矩陣非奇異的新方程組，將驗證原方程組多重根的奇異問題轉化為驗證新方程組單根的非奇異問題.基于此方法，提出了一種可信驗證算法，該算法輸出一個近似解及其相應的誤差界，使得在近似解的誤差界范圍內必存在一個精確解.本文推廣了文獻［2］中Jacobi矩陣秩虧為1的二重根驗證方法與文獻［3］中Jacobi矩陣秩虧為1的多重根驗證方法.

1 預備知識

記實數區間集合為I（?），矩陣A 的第i 行為Ai，：＝（Ai，1，Ai，2，…，Ai，n），q 階單位陣為Iq.令f：D??n→?n為非線性系統，x＊∈?n為f（x）＝0的解，Jf（x＊）為f（x）在x＊處的Jacobi矩陣.如果其秩為r＝n－q，則稱Jf（x＊）秩虧為q（1≤q≤n）.本文假設所討論的非線性系統f滿足下列條件：

（H1）f在x＊的鄰域內滿足C2－Lipschitz條件；

（H2）Jf（x＊）秩虧為q；

（H3）存在非零向量μ＊∈Null（（Jf（x＊））T）和Null（Jf（x＊））的一組基η＊＝｛η＊1，…，η＊q｝，使得q×q階矩陣

非奇異.其中q×n階矩陣（η＊）T（（μ＊）TJJf（x＊））的第k列為

對足夠接近x＊的x，令η（x）和h（x）為

的唯一解，μ（x）和g（x）為

的唯一解.其中η（x）和h（x）分別為n×q和q×q階矩陣，μ（x）與g（x）分別為n維和q維向量，α為隨機選取的q維向量.

定義q×q階矩陣

基于上述符號，定義邊界系統

其中λ為q維變元.由F（x，λ）的Jacobi矩陣為

可得：

由上述討論及引理3可得：

證明：由式（1）及條件（H3）（η＊為 Null（Jf（x＊））的一組基），有

進而

由式（5）知g（x＊）＝hT（x＊）·α＝0.于是有

故（x＊，0）為邊界系統F（x，λ）＝0的根.再由引理3知（x＊，0）為邊界系統F（x，λ）＝0的單根.

2 主要結果

問題1 設非線性系統f：D??n→?n滿足條件（H1）～（H3），給定其近似解∈?n，如何確定區間向量X，使得f（x）的精確解x＊∈＋X.

Rump［5］給出了系數陣為一般稠密矩陣的線性方程組求解的可信驗證方法，其中系數陣可以是浮點矩陣，也可以是區間矩陣.

定理2［1］給定矩陣A，T∈?n×n，向量b∈?n×n，區間向量X?I（?n×n）.如果

成立，則矩陣A，T均非奇異，且A－1b∈Tb＋（I－TA）X.

實現線性方程組求解的可信驗證算法函數是INTLAB函數包中的Verifylss函數［1］.對于系數陣為區間矩陣的線性方程組，Verifylss函數輸出區間向量，該區間向量包含此區間線性方程組所有可能的解.

定理3［1］給定區間矩陣∈I（?n×n）和區間向量∈I（?n×n），如果函數Verifylss運行成功，則該函數計算得到的區間向量X?I（?n×n）滿足

區間Newton法利用區間運算，通過在點Newton法的基礎上引進區間變量，構成了Newton法的區間變形，使所得的迭代程序在每次迭代過程中都產生解的界限，從而不僅得到解的近似，同時還可得到相應的誤差.區間Newton法的這個特點，顯然是一般點迭代不具備的.Krawczyk［6］針對區間Newton程序在計算量方面的缺陷，提出一種改進的區間Newton法，建立了不需要計算區間矩陣逆的Krawczyk算子.

定義1 設f：?n→?，若存在區間值映象F：I（?n）→I（?），使得

成立，則稱F為函數f的區間擴展.

定義2 設F：I（?n）→I（?），X，Y∈I（?n）且滿足X?Y，若有F（X）?F（Y），則稱區間值函數F具有包含單調性.

定義3 設函數f：D??n→?n連續可微，X∈I（D），JF是Jf的具包含單調性的區間擴展，Y為任意非奇異矩陣，稱

為Krawczyk算子.

Krawczyk算子的性質：如果X∩K（y，X）＝?，則X中不包含非線性方程組f（x）＝0的解.這個性質給出了解存在與否的檢驗條件.在區間Newton迭代過程中，如果遇到這種情況，則迭代停止.

在Krawczyk算子K（y，X）中，除y可在X中任取外，非奇異矩陣T也是任意的.且實矩陣T的選擇恰當與否直接關系到區間迭代法的收斂速度.取y＝m（X），T＝［m（JF（X））］－1，則相應的K（y，X）為

式（6）稱為Krawczyk算子的Moore形式.由此可構造Krawczyk區間Newton迭代法，取X（0）＝X，

其中k＝0，1，2，….

定理4［7］設函數f：D??n→?n連續可微，JF是Jf的具包含單調性的區間擴展，若初始區間向量X（0）?I（D）是一個n維立方體，K（X（0））滿足K（X（0））?int（X（0）），則非線性方程組f（x）＝0在X（0）中有唯一解，且序列｛X（k）｝至少線性地收斂于該唯一解.

應用可信驗證方法的前提是已知條件能在計算機上經過驗證.在Krawczyk［6］給出了驗證非線性系統解存在性的區間Newton法基礎上，Rump［8］做了進一步的研究，改進區間Newton法使其能更便于實際應用.

定理5［8］設函數f：D??n→?n，其中f＝（f1，…，fn）∈C1.給定向量∈?n，區間向量X∈I（?n），且0∈X，矩陣T∈?n×n，且給定的區間矩陣M∈I（?n×n）滿足條件：

如果

成立，則存在唯一的向量x＊∈＋X，使得f（x＊）＝0，并且每個矩陣∈M都是非奇異的，其中int（X）表示區間X的內部，▽表示梯度.

本文將區間Newton算法應用于文獻［4］中非線性方程組的數值解法，以達到數值計算結果的可信性驗證.令X∈I（?n）且0∈X，（）∈?n＋q為F（x，λ）＝0的近似解，T為JF（）的近似逆矩陣.首先要找到區間矩陣N∈I（?（n＋q）×（n＋q））滿足

由式（4）知

利用區間Newton算子的性質，如果

基于上述理論，設計算法如下.

算法1

2）初始化：利用區間轉換函數intval，得到初始區間向量

返回（X，Λ），算法終止；

由定理5，可得下述命題.

命題1 如果算法1成功返回的包含區間＋X和的包含區間＋Λ，則必存在唯一的向量使得（x＊，0）為F（x，λ）＝0的精確解.進而存在唯一的向量使得x＊為f（x）＝0的精確解.

證明：由數值Newton算法可計算出邊界系統F（x，λ）＝0的數值迭代增量‖Δ（x（k），λ（k））‖＜ε2的近似解如果算法1成功返回的包含區間和的包含區間，則由定理5可知，必存在唯一的向量使得F（x＊，0）＝0.再由邊界系統的定義知

3 數值算例

本文數值實驗基于 Windows7操作系統，軟件分別是 MAPLE 15（Digits∶＝14）和 MATLAB R2011a（INTLAB V6）.在MAPLE和MATLAB中執行可信驗證算法1，可計算出非線性系統的近似解及相應的誤差界，使得在近似解的誤差界范圍內必存在一個精確解.

輸出：

例1的解x＊＝（0，0）為孤立根.

例2 考慮非線性系統

輸出：

例2的解x＊＝（1，1，0，0，0）為孤立根.

［1］Rump S M.Verification Methods：Rigorous Results Using Floating－Point Arithmetic［J］.Acta Numerica，2010，19：287－449.

［2］Rump S M，Graillat S.Verified Error Bounds for Multiple Roots of Systems of Nonlinear Equations ［J］.Numerical Algorithms，2010，54（3）：359－377.

［3］LI Nan，ZHI Lihong.Verified Error Bounds for Isolated Singular Solutions of Polynomial Systems：Case of Breadth One［J］.Theoretical Computer Science，2013，479（1）：163－173.

［4］SHEN Yunqiu，Ypma T J.Newton’s Method for Singular Nonlinear Equations Using Approximate Left and Right Nullspaces of the Jacobian［J］.Applied Numerical Mathematics，2005，54（2）：256－265.

［5］Rump S M.Kleine Fehlerschranken bei Matrixproblemen［D］.Karlsruhe：Universit?t Karlsruhe，1980.

［6］Krawczyk R.Newton－Algorithmen zur Bestimmung von Nullstellen mit Fehlerschranken［J］.Computing，1969，4（3）：187－201.

［7］Moore R E.A Computational Test for Convergence of Iterative Methods for Nonlinear Systems ［J］.SIAM Journal on Numerical Analysis，1978，15（6）：1194－1196.

［8］Rump S M.Solving Algebraic Problems with High Accuracy ［C］／／Proceedings of the Symposium on a New Approach to Scientific Computation.San Diego：Academic Press，1983：51－120.