matlab中牛顿法程序

Posted 2023-05-15

tags:

篇首语：本文由小常识网(cha138.com)小编为大家整理，主要介绍了matlab中牛顿法程序相关的知识，希望对你有一定的参考价值。

我手里有两个牛顿迭代法的程序，但是两种程序计算出来的数值有误差，而我又对数值的要求很高，我是初学牛顿法，所以不知道应该采用哪个。请大家帮我看看！多谢！最好能解释一下为什么采用这个，多谢！
第一种：
function x=nanewton(fname,dfname,x0,e,N)
if nargin<5,N=500;
end
if nargin<4,e=1e-4;
end
x=x0;
x0=x+2*e;
k=0;
while abs(x0-x)>e&k<N,
k=k+1;x0=x;x=x0-feval(fname,x0)/feval(dfname,x0);
disp(x);
end
if k==N,warning('已达到迭代次数上限');
end

第二种：
function x=nanewton1(fname,dfname,x0,e,N)
if nargin<5
N=500;
end
if nargin<4
e=1e-4;
end
x=x0+2*e;
k=0;
while abs(x0-x)>e&(k<N)
k=k+1;
x0=x;
x=x0-feval(fname,x0)/feval(dfname,x0);
disp(x)
if k==N
warning('已达迭代次数上限');
end
end

采用第一个。
首先你的两个代码的计算过程和方法以及步骤是一致的。
只不过第二个将k==N放在循环内部判断是没有必要的。
放在while外面，可以节省点计算量。

如果你要求结果精度高一些的话，你调用：
x=nanewton1(fname,dfname,x0,e,N)
时e要小一些，比如说取1e-6这样。
另外：
if nargin<4
e=1e-4; %这个值也下调几个量级，作为缺省的精度。
end 参考技术A x0,tol,N)%改了这里
%root是系统保留字。不推荐使用。可以用exist('root')命令测试
% Newton Method
% The first parameter f is a external function with respect to viable x.
% The second parameter df is the first order diffential function of fx.
% x0 is initial iteration point.
% tol is the tolerance of the loop.
% N is the maximum number of iterations.
x=x0;
f0=eval(f);df0=eval(df);
n=0;
disp(' [ n xn xn+1 fn+1 ]');
while n<=N
x1=x0-f0/df0;
x=x1;
f1=eval(f);
X=[n,x0,x1,f1];
disp(X);
if abs(x0-x1)<tol
fprintf('The procedure was successful.\n')%改了这里
kk=X;%改了这里
return
else
n=n+1;
x0=x1;f0=f1;
end
end
if n==N+1
fprintf('the method failed after N iterations. '),
kk=0;%改了这里
end

无约束优化问题中牛顿法与拟牛顿法四种迭代方法的matlab实现

文章目录

1. 无约束优化问题的解法

在无约束优化问题中，有四种经典的迭代优化方法：Newton’s method(牛顿法)、Levenberg-Marquardt’s method(非线性最小二乘法,LM)、Broyden-Fletcher-Goldfarb-Shanno’s method(BFGS)、Davidon-Fletcher-Powell’s method(DFP)。

2. Matlab实现

假设有如下无约束优化问题： $f(x_1,x_2)=(x_1-3)^4+(x_1-3x_2)^2$ 为了方便四种算法的比较，我们统一设置初始迭代点为 $x_0=[0,0]^T$ ，则初始海森矩阵为： $H_0=\\begin{pmatrix} 110 & -6 &\\\\ -6 & 18 &\\\\ \\end{pmatrix}$

上述优化问题用matlab代码表示为：

syms x1; % 变量x1
syms x2; % 变量x2
f = (x1 - 3).^4 + (x1 - 3*x2).^2; $ 函数表达
x0=[0 0]'; % 初始迭代点
H0=[110 -6;-6 18]; % 初始海森矩阵
m=2; % 变量个数
k=30; % 迭代次数

该函数的最优 $X=(x_1,x_2)$ 以及 $f(x_1,x_2)$ 应该为： $x_1=3, x_2=1$ $f(x_1,x_2)=0$

下面展示一下如何用Matlab实现对函数的优化：

2.1. Newton’s method(牛顿法)

牛顿法的迭代公式为： $x_{k+1}=x_k-H^{-1}(x_k)\\nabla f(x_k)$ matlab实现函数实现为：

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%文件名:Newton.m
%
%f表示函数表达式
%H0表示初始的海森矩阵
%x0表示初始的迭代点 为列向量
%m表示变量的个数
%k表示迭代次数
%X存储每次迭代的x,F为函数值，G为每次的梯度，H为海森阵，HN为海森矩阵的逆
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

function[X, H, F, G,HN] = Newton(f,H0,x0,m,k)
    x1 = sym('x',[1,m]); % [x1, x2]
    %f = (x1(1) - 3)^4 + (x1(1) - 3*x1(2))^2;
    c = num2cell(x1); % c=变量[x1, x2]
    g = sym('x',[m,1]); % [x1, x2]^T
    
    X = zeros(m, k+1); % x1、x2的迭代值
    H = zeros( m, m, k+1); % hessian的迭代值
    F = zeros(1, k+1); % function的迭代值
    G = zeros(m, k+1); % function‘的迭代值
    HN = zeros( m, m, k+1); % hessian的逆阵的迭代值
    
    H(:,:,1) = H0; % hessian初始化
    HN(:,:,1) = inv(H0); % hessian逆初始化
    X(:,1) = x0; % X(x1, x2)初始化
    F(1,1) = subs(f, c, {X(:,1)'}); % 初始X值赋予F
    h = hessian(f,x1);%求海森矩阵
    
    for n = 1:m % f对x1、x2分别求偏导
        g(n) = diff(f,x1(n));
    end
     G(:,1) = subs(g,c,{X(:,1)'}); % 初始X导赋予G
     
     % 迭代
    for n = 1:k
        X(:,n+1) = X(:,n) - (H(:,:,n))\\G(:,n);
        F(1,n+1) = subs(f,c,{X(:,n+1)'}); 
        G(:,n+1) = subs(g,c,{X(:,n+1)'});
        H(:,:,n+1) = subs(h,c,{X(:,n+1)'});
        HN(:,:,n+1) = inv(H(:,:,n+1));
    end
end

执行matlab代码：

[X, H, F, G, HN] = Newton(f,H0,x0,m,k);

即可得到优化结果，下表是迭代次数 k 分别为：0、1、2、3时的输出值：

$k$	$x_k$	$f(x_k)$	$\\nabla f(x_k)$	$H(x_k)$
0	$(0, 0)$	81	$(- 108, 0)$	$\\begin{pmatrix}110 & -6 &\\\\-6 & 18 &\\\\\\end{pmatrix}$
1	$(1, 0.3333)$	16	$(- 32, 0)$	$\\begin{pmatrix}50 & -6 &\\\\-6 & 18 &\\\\\\end{pmatrix}$