拟牛顿法-DFP算法举例与matlab代码实现(转载＋整理)

算法来自于[1],如下：
值得一提的是,[1]中的python代码实现了对Rosenbrock函数的求极值测试.
例子来自于[2]:

----------------------------------------------------
用DFP算法求解：
minf(x)=x12+2x22−2x1x2−4x1minf(x)=x_1^2+2x_2^2-2x_1x_2-4x_1minf(x)=x12+2x22−2x1x2−4x1
取x0=(1,1)T,H0=x_0=(1,1)^T,H_0=x0=(1,1)T,H0=
[1001]\left[ \begin{matrix} 1 & 0 \\ 0&1 \\ \end{matrix} \right] [1001]

----------------------------------------------------

解答:
g(x)=(2x1−2x2−4,−2x1+4x2)Tg(x)=(2x_1-2x_2-4,-2x_1+4x_2)^Tg(x)=(2x1−2x2−4,−2x1+4x2)T
g0=(−4,2)Tg_0=(-4,2)^Tg0=(−4,2)T
p0=−H0g0=(4,−2)Tp_0=-H_0g_0=(4,-2)^Tp0=−H0g0=(4,−2)T,
(i)求迭代点x1x_1x1,令
ϕ0(α)=f(x0+αp0)=40α2−20α−3\phi_0(\alpha)=f(x_0+\alpha p_0)=40\alpha^2-20\alpha-3ϕ0(α)=f(x0+αp0)=40α2−20α−3,
得到ϕ(α)\phi(\alpha)ϕ(α)的极小值点为α0=14\alpha_0=\frac{1}{4}α0=41,
所以得：
x1=x0+α0p0=(2,0.5)T,g1=(−1,−2)T,x_1=x_0+\alpha_0p_0=(2,0.5)^T,g_1=(-1,-2)^T,x1=x0+α0p0=(2,0.5)T,g1=(−1,−2)T,
s0=x1−x0=(1,−0.5)T,y0=g1−g0=(3,−4)Ts_0=x_1-x_0=(1,-0.5)^T,y_0=g_1-g_0=(3,-4)^Ts0=x1−x0=(1,−0.5)T,y0=g1−g0=(3,−4)T
这里的 s0s_0s0是因为需要满足一个拟Newton条件，可以参考[4]

于是,
由DFP修正公式有H1=H0−H0y0y0TH0y0TH0y0+s0s0Ty0Ts0=1100H_1=H_0-\frac{H_0 y_0 y_0^TH_0}{y_0^TH_0y_0}+\frac{s_0s_0^T}{y_0^Ts_0}=\frac{1}{100}H1=H0−y0TH0y0H0y0y0TH0+y0Ts0s0s0T=1001[84383841]\left[ \begin{matrix} 84 & 38 \\ 38&41 \\ \end{matrix} \right][84383841]

所以下一个搜索方向为p1=−H1g1=15(8,6)Tp_1=-H_1g_1=\frac{1}{5}(8,6)^Tp1=−H1g1=51(8,6)T

(2)求迭代点x2
令
ϕ1(α)=f(x1+αp1)=85α2−4α−5.5\phi_1(\alpha)=f(x_1+\alpha p_1)=\frac{8}{5}\alpha^2-4\alpha -5.5ϕ1(α)=f(x1+αp1)=58α2−4α−5.5,

得到ϕ(α)\phi(\alpha)ϕ(α)的极小值点α1=54\alpha_1=\frac{5}{4}α1=45
于是得:
x2=x1+α1p1=(4,2)T,g2=(0,0)T,所以:x∗=x2=(4,2)T,f∗=−8x_2=x_1+\alpha_1p_1=(4,2)^T,g_2=(0,0)^T,所以:x^*=x_2=(4,2)^T,f^*=-8x2=x1+α1p1=(4,2)T,g2=(0,0)T,所以:x∗=x2=(4,2)T,f∗=−8
因为Hessian矩阵G(x)=G=
[2−2−24]T\left[ \begin{matrix} 2 & -2 \\ -2&4 \\ \end{matrix} \right]^T[2−2−24]T
为正定矩阵，f(x)f(x)f(x)为严格凸函数，所以x∗x*x∗为整体极小点

[3]提供了matlab代码,建立一个文件DFP.m(文件名必须和代码中的函数名保持一致)，代码如下：

function [best_x,best_fx,count]=DFP(x0,ess)
colormap Jet
% ###########################
syms x1 x2 t;
f=x1*x1+2*x2*x2-2*x1*x2-4*x1;
fx=diff(f,x1);%求表达式f对x1的一阶求导
fy=diff(f,x2);%求表达式f对x2的一阶求导
fi=[fx fy];%构造函数f的梯度函数
%初始点的梯度和函数值
g0=subs(fi,[x1 x2],x0);
f0=subs(f,[x1 x2],x0);
H0=eye(2); %输出x0,f0,g0
x0
f0
g0
xk=x0;
fk=f0;
gk=g0;
Hk=H0;
k=1;
while(norm(gk)>ess)%迭代终止条件||gk||<=ess   disp('************************************************************')     disp(['第' num2str(k) '次寻优'])
%确定搜索方向   pk=-Hk*gk';
%由步长找到下一点x(k+1)    xk=xk+t*pk';     f_t=subs(f,[x1 x2],xk); %构造一元搜索的一元函数φ(t) %由一维搜索找到最优步长    df_t=diff(f_t,t);    tk=solve(df_t);
if tk~=0         tk=double(tk);
elsebreak;
end
%计算下一点的函数值和梯度xk = subs(xk,t,tk)    fk=subs(f,[x1 x2],xk)    gk0=gk;     gk=subs(fi,[x1 x2],xk)
%DPF校正公式，找到修正矩阵    yk=gk-gk0;    sk=tk*pk';Hk=Hk-(Hk*yk'*yk*Hk)/(yk*Hk*yk')+sk'*sk/(yk*sk')%修正公式    k=k+1;
enddisp('结果如下：')
best_x=xk;%最优点
best_fx=fk;%最优值
count=k-1;
end

matlab终端运行方法如下:
>> x0=[1 1];
>> ess=1e-6
>> [best_x,best_fx,count]=DFP(x0,ess)

输出如下：

x0 =1     1f0 =-3g0 =[ -4, 2]************************************************************
第1次寻优xk =[ 2, 1/2]fk =-11/2gk =[ -1, -2]Hk =[ 21/25,  19/50]
[ 19/50, 41/100]************************************************************
第2次寻优xk =[ 4, 2]fk =-8gk =[ 0, 0]Hk =[   1, 1/2]
[ 1/2, 1/2]结果如下：best_x =[ 4, 2]best_fx =-8count =2>> x0=[1 1];
>> ess=1e-6ess =1.0000e-06

Reference:

[1]优化算法——拟牛顿法之DFP算法
[2]拟牛顿法－最优化方法－百度文库
[3]DFP算法及Matlab程序
[4]DFP算法

拟牛顿法-DFP算法举例与matlab代码实现(转载＋整理)相关推荐

粒子群(pso)算法详解matlab代码,粒子群(pso)算法详解matlab代码
粒子群(pso)算法详解matlab代码 (1)---- 一.粒子群算法的历史粒子群算法源于复杂适应系统(Complex Adaptive System,CAS).CAS理论于1994年正式提出,C ...
智能优化算法：海鸥算法原理及Matlab代码
智能优化算法:海鸥算法原理及Matlab代码文章导航 1. 算法原理 2. SOA算法流程 2.1 迁徙行为(exploration ability) 2.2 攻击行为(exploitation a ...
music的matlab程序,DOA经典算法MUSIC的MATLAB代码(作者：Nikhil Shetty).pdf
DOA经典算法MUSIC的MATLAB代码(作者:Nikhil Shetty) Direction of arrival estimation algorithms 作者:Nikhil Shetty ...
MUSIC算法原理及MATLAB代码阵列信号处理
MUSIC算法原理及MATLAB代码阵列信号处理 MUSIC(multiple signal classification algorithm)算法是一种基于矩阵特征空间分解的方法.从几何角度讲,信 ...
果蝇（FOA）优化算法(附完整Matlab代码，可直接复制）
果蝇优化算法的核心是利用果蝇搜索实物的机制来对问题进行寻优.果蝇根据气味来确定食物的位置,食物腐烂程度越高,气味越大,果蝇对其越敏感.果蝇的觅食行为如下图: Fig1. 果蝇觅食行为示意图 1.果蝇位 ...
COI实验室技能：常见的图像增强算法（含MATLAB代码）
COI实验室技能:常见的图像增强算法(含MATLAB代码) 成像系统采集到的图片有时不一定具有很好的图像质量,往往需要进行一些图像增强的操作.本文主要整理了几种常见的图像增强算法,包括:同态滤波. ...
回声消除AEC算法（含Matlab代码）
基于自适应滤波器的回声消除AEC算法(含Matlab代码) 摘要自适应滤波器声学回波抵消 AEC算法解析 LMS算法 NLMS算法 VSNLMS算法 APLMS算法 LMS-Newton算法 PF ...
孔雀优化算法(POA)——(含MATLAB代码)
先要声明的是,这个算法的争议较大,评论是两边倒的,一些学者觉得是一锅乱炖,我个人不发表观点,这里只是把这个算法介绍一下.我不推荐它,也不反对大家去了解,只做单纯的讲解.觉得这个算法不行的就不看吧~ 孔 ...
遗传算法及其matlab实现,现代智能化算法—遗传算法及其matlab代码实现
遗传算法(Genetic Algorithm)是模拟达尔文生物进化论的自然选择和遗传学机理的生物进化过程的计算模型,是一种通过模拟自然进化过程搜索最优解的方法,它最初由美国Michigan大学J.Ho ...

拟牛顿法-DFP算法举例与matlab代码实现(转载＋整理)

拟牛顿法-DFP算法举例与matlab代码实现(转载＋整理)相关推荐

最新文章

热门文章