数学软件 之 基于MATLAB的DFP算法

　　DFP算法是本科数学系中最优化方法的知识，也是无约束最优化方法中非常重要的两个拟Newton算法之一，上一周写了一周的数学软件课程论文，姑且将DFP算法的实现细节贴出来分享给学弟学妹参考吧，由于博客不支持数学公式，所以就不累述算法原理及推倒公式了。


DFP算法流程图

先给出DFP算法迭代流程图，总体上是拟Newton方法的通用迭代步骤，唯独在校正公式的地方有所区别。

MATLAB实现DFP

　　基于此图便可以设计DFP算法的MATLAB程序：

　　对分法及加步探索法的实现

　　首先由于DFP算法中需要利用一维搜索得到最优步长，因此需要先设计一个一维搜索函数，博主选用的是简单的对分法（二分法）：

　　

%本函数利用二分法求解X = ls(Xk,Pk)问题
%目标函数：f
%符号参数：var
%终止限：eps
function x = dichotomy(f,var,eps)
    g = diff(f,var);
    [a, b] = search(f,var);
    x = (a + b)/2;      %防止eps过大导致x无值
    while b - a > eps
        x = (a+b)/2;
        gx = subs(g, var, x);
        if gx > 0
            b = x;
        elseif gx < 0
            a = x;
        else break;
        end
    end
    %加步搜索法-确定搜索区间
    function [a, b] = search(g,var)
        gt = matlabFunction(g);
        X = 0;  tmp = X;
        h = 1;  k = 0;
        while 1
            Xk = X + h; k = k+1;
            Y = subs(gt,var,X);
            Yk = subs(gt,var,Xk);
            if Y > Yk   %加大步长搜索
                h = 2 * h;
                tmp = X;
                X = Xk;
            elseif Y == Yk  %缩小步长搜索
                h = h/2;
            elseif k == 1
                h = -h;     %反向搜索
            else break;
            end
        end
        a = min(tmp, Xk);
        b = max(tmp, Xk);
    end
end

　　DFP算法的实现

　　有了一维搜索函数，那么实现DFP算法也就能依照算法流程图来设计了：

　　

%DFP算法主程序
%目标函数：f
%初始点：X0
%参数：var
%终止限：eps
function DFP(f, X0, var, eps)
%初始化符号函数，梯度，维数等
syms var t;
g = jacobian(f)‘;   %Jacobian转置->Grad
fx = matlabFunction(f); %符号函数->函数句柄(R2009以上支持)
gx = matlabFunction(g);
n = length(var);    %维数
X = X0; Xk = X0;
while 1
    fx0 = fx(X(1),X(2));   gx0 = gx(X(1),X(2));
    Hk = eye(2);    Pk = -gx0;  %初始方向
    k = 0;  %迭代次数
    while 1
        Y = Xk + t*Pk;
        y = fx(Y(1),Y(2));
        tk = dichotomy(y, t, eps);  %一维搜索
        Xk = Xk + tk*Pk;
        fx1 = fx(Xk(1),Xk(2));
        gx1 = gx(Xk(1),Xk(2));
        if norm(gx1) < eps || k == n
            X = Xk; fx0 = fx1;
            break;
        end
        Sk = Xk - X;    Yk = gx1 - gx0;
        Hk = Hk + Sk*Sk‘/(Sk‘*Yk) - Hk*(Yk)*Yk‘*Hk/(Yk‘*Hk*Yk);
        Pk = -Hk*gx1;   %校正方向
        k = k+1;
    end
    if norm(gx1) < eps
        disp(‘X(k+1) = ‘);  disp(Xk);
        disp(‘F(K+1) = ‘);  disp(fx0);
        break;
    end
end

实例验证

　　有了DFP算法的实现函数，那么应用于实例也就不难了。

　　可以在命令文件下输入如下代码就能得到目标函数极值点及极值

clear; clc; format long;
syms x1 x2;
f = 4*(x1-5)^2 + (x2-6)^2;
tic;    %初始时间
DFP(f, [8;9], [x1, x2], 0.00000001);
toc;    %结束时间

Normal
0

7.8 磅
0
2

false
false
false

EN-US
ZH-CN
X-NONE

/* Style Definitions */
table.MsoNormalTable
{mso-style-name:普通表格;
mso-tstyle-rowband-size:0;
mso-tstyle-colband-size:0;
mso-style-noshow:yes;
mso-style-priority:99;
mso-style-parent:"";
mso-padding-alt:0cm 5.4pt 0cm 5.4pt;
mso-para-margin:0cm;
mso-para-margin-bottom:.0001pt;
mso-pagination:widow-orphan;
font-size:10.0pt;
font-family:"Times New Roman",serif;}

输出结果如下：

X(k+1) =

4.999995811278565

5.999767686222325

F(K+1) =

5.403987284687523e-08

Elapsed time is 8.229108 seconds.

　　算法时间度分析：

Normal
0

7.8 磅
0
2

false
false
false

EN-US
ZH-CN
X-NONE

/* Style Definitions */
table.MsoNormalTable
{mso-style-name:普通表格;
mso-tstyle-rowband-size:0;
mso-tstyle-colband-size:0;
mso-style-noshow:yes;
mso-style-priority:99;
mso-style-parent:"";
mso-padding-alt:0cm 5.4pt 0cm 5.4pt;
mso-para-margin:0cm;
mso-para-margin-bottom:.0001pt;
mso-pagination:widow-orphan;
font-size:10.0pt;
font-family:"Times New Roman",serif;}

由此可知，函数 在[8,9]附近的点[5.00,6.00]处取得局部最小值，其中局部极值点约为5.40e-8.

此算法运行时间约为8.23s，并且我们在降低终止限eps后，针对本题，算法运行时间增长较快，例如若eps = 1e-3，耗时11.6s，若eps = 1e-5，耗时22.94s，而eps = 1e-7，耗时甚至超过15分钟.这说明DFP算法在求解高精度运算时的运行效率表现得并不是那么好，甚至有可能无法得出最优解.

　　

　　实例搜索图

　　基于该实例，对算法的迭代过程进行绘图，得到如下搜索图

　　

Normal
0

7.8 磅
0
2

false
false
false

EN-US
ZH-CN
X-NONE

/* Style Definitions */
table.MsoNormalTable
{mso-style-name:普通表格;
mso-tstyle-rowband-size:0;
mso-tstyle-colband-size:0;
mso-style-noshow:yes;
mso-style-priority:99;
mso-style-parent:"";
mso-padding-alt:0cm 5.4pt 0cm 5.4pt;
mso-para-margin:0cm;
mso-para-margin-bottom:.0001pt;
mso-pagination:widow-orphan;
font-size:10.0pt;
font-family:"Times New Roman",serif;}

可以由以上两个搜索图像得出一个结论：DFP算法的实质是在每一次迭代过程中调整自己的搜索方向，以使得该方向能够尽可能接近极值点，这也正是几乎所有拟Newton算法中校正矩阵的作用.

目   录

引言1

1 
DFP算法原理1

1.1  DFP算法基本原理1

1.2  DFP算法迭代步骤2

1.3  DFP算法流程图4

2 
DFP算法实现4

2.1  对分法及加步探索法的实现4

2.2  DFP算法的实现6

3 
实例分析7

3.1  实例求解7

3.2  DFP的实例搜索过程图8

结论9

参考文献10

Normal
0

7.8 磅
0
2

false
false
false

EN-US
ZH-CN
X-NONE

/* Style Definitions */
table.MsoNormalTable
{mso-style-name:普通表格;
mso-tstyle-rowband-size:0;
mso-tstyle-colband-size:0;
mso-style-noshow:yes;
mso-style-priority:99;
mso-style-parent:"";
mso-padding-alt:0cm 5.4pt 0cm 5.4pt;
mso-para-margin:0cm;
mso-para-margin-bottom:.0001pt;
mso-pagination:widow-orphan;
font-size:10.0pt;
font-family:"Times New Roman",serif;}

时间： 2024-10-31 09:47:50

数学软件 之 基于MATLAB的DFP算法的相关文章

Mathematica 和 MATLAB、Maple 并称为三大数学软件

Mathematica是一款科学计算软件,很好地结合了数值和符号计算引擎.图形系统.编程语言.文本系统.和与其他应用程序的高级连接.很多功能在相应领域内处于世界领先地位,它也是使用最广泛的数学软件之一.Mathematica的发布标志着现代科技计算的开始.Mathematica是世界上通用计算系统中最强大的系统.自从1988发布以来,它已经对如何在科技和其它领域运用计算机产生了深刻的影响. Mathematica 和 MATLAB.Maple 并称为三大数学软件. Mathematica的功能包

基于matlab mex的平面点集按重心逆时针排序算法

基于matlab mex的平面点集按重心逆时针排序算法,可用于求凸集,代码如下: #include <mex.h> #include <mat.h> #include <stdlib.h> #include <vector> #include <algorithm> using namespace std; typedef struct PointF { int id; double x, y; PointF() { id = -1; } Poi

基于MATLAB的腐蚀膨胀算法实现

本篇文章要分享的是基于MATLAB的腐蚀膨胀算法实现,腐蚀膨胀是形态学图像处理的基础,腐蚀在二值图像的基础上做“收缩”或“细化”操作,膨胀在二值图像的基础上做“加长”或“变粗”的操作. 什么是二值图像呢?把一幅图片看做成一个二维的数组,那么二值图像是一个只有0和1的逻辑数组,我们前面Sobel边缘检测后的图像输出边缘效果,设置个阈值,大于阈值输出为1,小于阈值输出为0,最后输出就是一幅二维图像了. 腐蚀 腐蚀是一种消除边界点,使边界向内部收缩的过程.可以用来消除小且无意义的物体.用3X3的结构元

基于MATLAB的多项式数据拟合方法研究-毕业论文

摘要:本论文先介绍了多项式数据拟合的相关背景,以及对整个课题做了一个完整的认识.接下来对拟合模型,多项式数学原理进行了详细的讲解,通过对文献的阅读以及自己的知识积累对原理有了一个系统的认识.介绍多项式曲线拟合的基本理论,对多项式数据拟合原理进行了全方面的理论阐述,同时也阐述了曲线拟合的基本原理及多项式曲线拟合模型的建立.具体记录了多项式曲线拟合的具体步骤,在建立理论的基础上具体实现多项式曲线的MATLAB实现方法的研究,采用MATLAB R2016a的平台对测量的数据进行多项式数据拟合,介绍了M

基于乐理模型的算法作曲研究

Abstraction Musical composition is a sophisticated movement of thought that is unique for humans. More concretely, composition reflects the basic technology and theory of music, i.e. harmonics, polyphony, orchestration, structure and so on, which are

《基于改进随机抽取算法的信息论题库和智能组卷系统的设计与实现》

一:基本信息 1标题:<基于改进随机抽取算法的信息论题库和智能组卷系统的设计与实现> 2时间:2015 3来源:中国地质大学硕士学位论文 4关键词:随机抽取,题库,智能组卷. 二:内容: 研究意义:为了真正实现教学与考试分离,推进无纸化教学改革的浪潮,使命题组卷更加科学合理.公正公平,应用先进的计算机技术完成试题数据库系统建设,当需要考卷时,使用系统的组卷功能就可以很方便地从试题库中抽取符合要求的试题直接组成试卷.这样,不仅节省了教师宝贵的时间和精力,提高教学考试工作效率,更加客观.科学.全面

ArcGIS水分分析工具的流向分析是基于D8单流向算法

ArcGIS水分分析工具的流向分析是基于D8单流向算法,如果分析使用的DEM存在凹陷点,就会产生汇,导致径流断流从而影响了分析结果.在前面章节<ArcGIS水文分析实战教程(2)ArcGIS水文分析工具的基本原理>中又介绍过D8算法,而<ArcGIS水文分析实战教程(4)地形预处理>章节中笔者也较少过如何创建无凹陷点得DEM数据,在使用流向分析工具之前可以先行阅读. 首先流向分析要使用填洼过的数据,确保DEM数据没有凹陷点.如果数据准备妥当,直接使用水文分析工具箱中的[流向]工具进

基于直方图的图像增强算法（HE、CLAHE、Retinex）之（二）

作为图像增强算法系列的第二篇文章,下面我们将要介绍功能强大.用途广泛.影响深远的对比度有限的自适应直方图均衡(CLAHE,Contrast Limited Adaptive Histogram Equalization)算法.尽管最初它仅仅是被当作一种图像增强算法被提出,但是现今在图像去雾.低照度图像增强,水下图像效果调节.以及数码照片改善等方面都有应用.这个算法的算法原理看似简单,但是实现起来却并不那么容易.我们将结合相应的Matlab代码来对其进行解释.希望你在阅读本文之前对朴素的直方图均衡

牛顿法与拟牛顿法学习笔记（三）DFP 算法

机器学习算法中经常碰到非线性优化问题,如 Sparse Filtering 算法,其主要工作在于求解一个非线性极小化问题.在具体实现中,大多调用的是成熟的软件包做支撑,其中最常用的一个算法是 L-BFGS.为了解这个算法的数学机理,这几天做了一些调研,现把学习过程中理解的一些东西整理出来. 目录链接 (1) 牛顿法 (2) 拟牛顿条件 (3) DFP 算法 (4) BFGS 算法 (5) L-BFGS 算法 作者: peghoty 出处: http://blog.csdn.net/itplus/