牛頓迭代法解非線性方程組（MATLAB版）

牛頓迭代法，又名切線法，這里不詳細介紹，簡單說明每一次牛頓迭代的運算：首先將各個方程式在一個根的估計值處線性化（泰勒展開式忽略高階余項），然后求解線性化后的方程組，最后再更新根的估計值。下面以求解最簡單的非線性二元方程組為例（平面二維定位最基本原理），貼出源代碼：

1、新建函數fun.m，定義方程組

 

function f=fun(x);
%定義非線性方程組如下
%變量x1 x2
%函數f1 f2
syms x1 x2
f1 = sqrt((x1-4)^2 + x2^2)-sqrt(17);
f2 = sqrt(x1^2 + (x2-4)^2)-5;
f=[f1 f2];

 

2、新建dfun.m，求出一階微分方程

 

function df=dfun(x);
f=fun(x);
df=[diff(f,'x1');diff(f,'x2')]; %雅克比矩陣

 

3、建立newton.m，執行牛頓迭代過程

 

clear;clc
format;
x0=[0 0];   % 迭代初始值
eps = 0.00001;  % 定位精度要求
for i = 1:10
    f = double(subs(fun(x0),{'x1' 'x2'},{x0(1) x0(2)}));
    df = double(subs(dfun(x0),{'x1' 'x2'},{x0(1) x0(2)}));  % 得到雅克比矩陣
    x = x0 - f/df;
    if(abs(x-x0) < eps)
        break;
    end
    x0 = x; % 更新迭代結果
end
disp('定位坐標：');
x
disp('迭代次數：');
i

 

結果如下：

定位坐標：

x =


    0.0000   -1.0000

迭代次數：

i =


     4