1 牛頓法簡(jiǎn)介

牛頓迭代法（Newton’s method）又稱為牛頓-拉夫遜（拉弗森）方法（Newton-Raphson method），它是牛頓在17世紀(jì)提出的一種在實(shí)數(shù)域和復(fù)數(shù)域上近似求解方程的方法。

多數(shù)方程不存在求根公式，因此求精確根非常困難，甚至不可解，從而尋找方程的近似根就顯得特別重要。方法使用函數(shù) $f(x)$ 的泰勒級(jí)數(shù)的前面幾項(xiàng)來尋找方程 $f(x)=0$ 的根。牛頓迭代法是求方程根的重要方法之一，其最大優(yōu)點(diǎn)是在方程 $f(x)=0$ 的單根附近具有平方收斂，而且該法還可以用來求方程的重根、復(fù)根，此時(shí)線性收斂，但是可通過一些方法變成超線 性收?。另外該方法廣泛用于計(jì)算機(jī)編程中。

2 牛頓法原理

設(shè) $r$ 是 $f(x)=0$ 的根，選取 $x_0$ 作為 $r$ 的初始近似值，過點(diǎn) $\left(x_0,f\left(x_0\right)\right)$ 做曲線 $y=f(x)$ 的切線 $L$ ，
$L:y=f\left(x_0\right)+f^{\prime }\left(x_0\right)\left(x?x_0\right)$ ，則 $L$ 與 $x$ 軸交點(diǎn)的橫坐標(biāo) $x_1=x_0?\frac{f\left(x_0\right)}{f^{\prime }\left(x_0\right)}$ ，稱 $x_1$ 為 $r$ 的一次近似值。過點(diǎn) $\left(x_1,f\left(x_1\right)\right)$ 做曲線 $y=f(x)$ 的切線，并求該切線與 $\times$ 軸交點(diǎn)的橫坐標(biāo) $x_2=x_1?\frac{f\left(x_1\right)}{f^{\prime }\left(x_1\right)}$ ，稱 $x_2$ 為 $\mathrm{r}$ 的二次近似值。重曷 以上過程，得 $r$ 的近似值序列，其中， $x_{n+1}=x_n?\frac{f\left(x_n\right)}{f^{\prime }\left(x_n\right)}$ 稱為 $r$ 的 $n+1$ 次近似值，上式稱為牛頓迭代公式。

用牛頓迭代法解非線性方程，是把非線性方程 $f(x)=0$ 線性化的一種近似方法。把 $f(x)$ 在點(diǎn) $x_0$ 的桌鄰域內(nèi)展開成泰勒 級(jí)數(shù) $f(x)=f\left(x_0\right)+f^{\prime }\left(x_0\right)\left(x?x_0\right)+\frac{f^{\prime \prime }\left(x_0\right){\left(x?x_0\right)}^2}{2!}+?+\frac{f^{(n)}\left(x_0\right){\left(x?x_0\right)}^n}{n!}+R_n(x)$ ，取其線性部分 (即泰勒展開的前兩項(xiàng))，并令其等于 0 ，即 $f\left(x_0\right)+f^{\prime }\left(x_0\right)\left(x?x_0\right)=0$ ，以此作為非線性方程 $f(x)=0$ 的近似方程， 若 $f^{\prime }\left(x_0\right)\ne 0$ ，則其解為 $x_1=x_0?\frac{f\left(x_0\right)}{f^{\prime }\left(x_0\right)}$ ，這樣，得到牛頓迭代法的一個(gè)朱代關(guān)系式: $x_{n+1}=x_n?\frac{f\left(x_n\right)}{f^{\prime }\left(x_n\right)}$ 。

已經(jīng)證明，如果是連續(xù)的，并且待求的零點(diǎn)是孤立的，那么在零點(diǎn)周圍存在一個(gè)區(qū)域，只要初始值位于這個(gè)鄰近區(qū)域內(nèi)，那 么牛頓法必定收斂。并且，如果不為 0 , 那么牛頓法將具有平方收斂的性能. 粗略的說，這意味著每造代一次，牛頓法結(jié)果的有效 數(shù)字將增加一倍。

3 牛頓法推導(dǎo)

4 Matlab代碼實(shí)現(xiàn)

下面用Matlab代碼求解上面的示例。

clear;clc;

% 定義原函數(shù)
syms xx yy
fy(xx,yy) = 0.5 * xx^2 + 2 * yy^2;

% 確定迭代次數(shù)
n = 10
% 確定初始點(diǎn)
x0 = 1
y0 = 1
% 求初始點(diǎn)函數(shù)值
fy(x0,y0)
% 求函數(shù)梯度
xf = -5:0.2:5;
yf = xf';
ff = 0.5 * xf.^2 + 2 * yf.^2;
surf(xf,yf,ff)
xlabel('x')
ylabel('y')
zlabel('z')
view([119.1 40.8])
[fx,fy] = gradient(ff,0.2);

% 提取點(diǎn)初始點(diǎn)處的梯度值
t = (xf == x0) & (yf == y0);
indt = find(t);
f_grad = [fx(indt) fy(indt)]
% 求海森矩陣
syms x y
f(x,y) = 0.5 * x^2 + 2 * y^2;
H = hessian(f,[x,y])
% 迭代
for i=1:n

    % 判斷是否可以跳出(如果梯度向量都接近0，就跳出)
    b = 0;
    for j = 1:length(f_grad)
        if f_grad(j) > 0.000001
            b = 1;
            break
        end
    end
    if b==0
        break
    end

    % 確定下降方向
    d = -inv(H)*(f_grad)';
    dk = d(x0,y0);
    
    % 確定步長(zhǎng)，牛頓法步長(zhǎng)為1
    a = 1;

    % 獲取下一狀態(tài)的點(diǎn)
    newX = [x0,y0] + dk' .* a
    x0 = newX(1);
    y0 = newX(2);

    % 更新梯度信息
    t = (xf == x0) & (yf == y0);
    indt = find(t);
    f_grad = [fx(indt) fy(indt)];

end

5 低版本Matlab報(bào)錯(cuò)

最近有朋友向我反應(yīng)代碼運(yùn)行會(huì)報(bào)錯(cuò)，具體報(bào)錯(cuò)內(nèi)容如下：

他使用的matlab版本是2016a，推測(cè)可能是低版本不支持(151)的矩陣和(511)的矩陣直接做運(yùn)算，如果大家有遇到這樣的報(bào)錯(cuò)的話，可以試一下將原代碼的16、17行刪去，換成以下代碼應(yīng)該就可以了：文章來源地址http://www.zghlxwxcb.cn/news/detail-406565.html

yf = xf;
s = size(xf,2);
ff = zeros(s,s);
for i = 1 : s
    for j = 1 : s
        obj = 0.5 * xf(i)^2 + 2*yf(j)^2;
        ff(i,j) = obj;
    end
end

到了這里，關(guān)于【最優(yōu)化理論】牛頓法+Matlab代碼實(shí)現(xiàn)的文章就介紹完了。如果您還想了解更多內(nèi)容，請(qǐng)?jiān)谟疑辖撬阉鱐OY模板網(wǎng)以前的文章或繼續(xù)瀏覽下面的相關(guān)文章，希望大家以后多多支持TOY模板網(wǎng)！