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

Posted 2021-11-04 白水baishui

tags:

篇首语：本文由小常识网(cha138.com)小编为大家整理，主要介绍了无约束优化问题中牛顿法与拟牛顿法四种迭代方法的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}$