求解非约束优化问题的拟牛顿方法(BFGS、DFP)

求解非约束优化问题的拟牛顿方法(BFGS、DFP)

拟牛顿法是一种以牛顿法为基础设计的,求解非线性方程组或连续的最优化问题函数的零点或极大、极小值的算法。当牛顿法中所要求计算的雅可比矩阵或Hessian矩阵难以甚至无法计算时,拟牛顿法便可派上用场。

数学原理

考虑模型问题

这里写图片描述

最小化之,可得拟牛顿步长:

这里写图片描述

那么,新的迭代为:

这里写图片描述

问题是,这里的 αk α k Bk B k 如何选取?

我们在割线方程的约束框架下,分别求解子问题:

这里写图片描述

和子问题:

这里写图片描述

可以得到 Hk=B1k H k = B k − 1 的迭代:

BFGS:

这里写图片描述

DFP:

这里写图片描述

当然,还有其他的一些 Hk H k 更新方式:

这里写图片描述

步长 αk α k 如何选取?只要保证Wolf条件,就能保证其收敛性。取精确线搜索步长时,其值为1。

算法简述

基本框架如下:

这里写图片描述

在迭代的过程中,我们不需要显性地去求 Bk B k ,我们用到的只是它的逆,所以,我们可以用其逆 Hk H k 的迭代公式进行搜索迭代。

程序代码

clc
clear

f = @(t)t(1)^2+2*t(2)^2;
x0 = [10,10]';
epsilon = 0.001;
H0 = eye(2);
method = 'BFGS';%DFP
x = quasi_Newton(f,x0,epsilon,H0,method);

function [xk,k] = quasi_Newton(f,x0,epsilon,H0,method)
%使用:quasi_Newton(f,x0,method)

if nargin < 3
    help mfilename;
end
k = 0;
syms t1 t2;
t = [t1,t2]';
fs = f(t);
dfs = gradient(fs);
df = matlabFunction(dfs);
df = @(x) df(x(1),x(2));
df0 = df(x0);
normdf = sqrt(df0'*df0);
H = H0;
xk = x0;
dfk = df0;
while normdf > epsilon
    p = -H*dfk;
    alpha = cal_alpha(H,dfk);
    xk1 = xk + alpha*p;
    dfk1 = df(xk1);
    sk = xk1 - xk;
    yk = dfk1 - dfk;
    eval(['H = ' method '(H,sk,yk);']);
    %H = BFGS(H,sk,yk);
    k = k + 1;
    xk = xk1;
    dfk = dfk1;
    normdf = sqrt(dfk'*dfk);
    xk
end
end
function alpha = cal_alpha(H,dfk)
    %alpha = dfk'*H*dfk/(dfk'*H'*dfk); 
    alpha = 1;
end
function H = BFGS(H,sk,yk)
    gammak = 1/(yk'*sk);
    skykT = sk * yk';
    skskT = sk * sk';
    E = eye(2);
    H = (E-gammak*skykT)*H*(E-gammak*skykT) + gammak*skskT;
end
function H = DFP(H,sk,yk)
    Hyk = H*yk;
    ykTsk =  yk'*sk;
    skskT = sk * sk';
    H = H - (Hyk*yk'*H)/(yk'*Hyk) + skskT/ykTsk;
end

算例

写一个程序用带精确线搜索步长的拟牛顿方法来求解以下问题:

这里写图片描述

初始点 x0=(1,1)T x 0 = ( 1 , 1 ) T 。分别使用BFGS和DFP的更新公式。为两种方法都设置初始矩阵 H0=I H 0 = I

通过简单的8步迭代,我们就达到了一个很高的精度,如下:

这里写图片描述

已标记关键词 清除标记
©️2020 CSDN 皮肤主题: 点我我会动 设计师:上身试试 返回首页