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

桑海風(fēng)，萬保成

（1.吉林大學(xué) 數(shù)學(xué)學(xué)院，長春130012；2.北華大學(xué) 數(shù)學(xué)與統(tǒng)計學(xué)院，吉林 吉林132013；3.吉林農(nóng)業(yè)大學(xué) 信息技術(shù)學(xué)院，長春130118）

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

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

1 預(yù)備知識

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

（H1）f在x＊的鄰域內(nèi)滿足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階矩陣

基于上述符號，定義邊界系統(tǒng)

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

可得：

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

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

進(jìn)而

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

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

2 主要結(jié)果

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

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

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

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

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

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

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

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

成立，則稱F為函數(shù)f的區(qū)間擴展.

定義2 設(shè)F：I（?n）→I（?），X，Y∈I（?n）且滿足X?Y，若有F（X）?F（Y），則稱區(qū)間值函數(shù)F具有包含單調(diào)性.

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

為Krawczyk算子.

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

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

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

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

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

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

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

如果

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

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

由式（4）知

利用區(qū)間Newton算子的性質(zhì)，如果

基于上述理論，設(shè)計算法如下.

算法1

2）初始化：利用區(qū)間轉(zhuǎn)換函數(shù)intval，得到初始區(qū)間向量

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

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

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

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

3 數(shù)值算例

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

輸出：

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

例2 考慮非線性系統(tǒng)

輸出：

例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.