BFGS优化算法的C++实现

头文件&＃xff1a;

/** Copyright (c) 2008-2011 Zhang Ming (M. Zhang), zmjerry&＃64;163.com** This program is free software; you can redistribute it and/or modify it* under the terms of the GNU General Public License as published by the* Free Software Foundation, either version 2 or any later version.** Redistribution and use in source and binary forms, with or without* modification, are permitted provided that the following conditions are met:** 1. Redistributions of source code must retain the above copyright notice,* this list of conditions and the following disclaimer.** 2. Redistributions in binary form must reproduce the above copyright* notice, this list of conditions and the following disclaimer in the* documentation and/or other materials provided with the distribution.** This program is distributed in the hope that it will be useful, but WITHOUT* ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or* FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License for* more details. A copy of the GNU General Public License is available at:* http://www.fsf.org/licensing/licenses*//****************************************************************************** bfgs.h** BFGS quasi-Newton method.** This class is designed for finding the minimum value of objective function* in one dimension or multidimension. Inexact line search algorithm is used* to get the step size in each iteration. BFGS (Broyden-Fletcher-Goldfarb* -Shanno) modifier formula is used to compute the inverse of Hesse matrix.** Zhang Ming, 2010-03, Xi&＃39;an Jiaotong University.*****************************************************************************/#ifndef BFGS_H
#define BFGS_H#include
#include namespace splab
{template class BFGS : public LineSearch{public:BFGS();~BFGS();void optimize( Ftype &func, Vector &x0, Dtype tol&＃61;Dtype(1.0e-6),int maxItr&＃61;100 );Vector getOptValue() const;Vector getGradNorm() const;Dtype getFuncMin() const;int getItrNum() const;private:// iteration numberint itrNum;// minimum value of objective functionDtype fMin;// optimal solutionVector xOpt;// gradient norm for each iterationVector gradNorm;};// class BFGS#include }
// namespace splab#endif
// BFGS_H

实现文件&＃xff1a;

/** Copyright (c) 2008-2011 Zhang Ming (M. Zhang), zmjerry&＃64;163.com** This program is free software; you can redistribute it and/or modify it* under the terms of the GNU General Public License as published by the* Free Software Foundation, either version 2 or any later version.** Redistribution and use in source and binary forms, with or without* modification, are permitted provided that the following conditions are met:** 1. Redistributions of source code must retain the above copyright notice,* this list of conditions and the following disclaimer.** 2. Redistributions in binary form must reproduce the above copyright* notice, this list of conditions and the following disclaimer in the* documentation and/or other materials provided with the distribution.** This program is distributed in the hope that it will be useful, but WITHOUT* ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or* FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License for* more details. A copy of the GNU General Public License is available at:* http://www.fsf.org/licensing/licenses*//****************************************************************************** bfgs-impl.h** Implementation for BFGS class.** Zhang Ming, 2010-03, Xi&＃39;an Jiaotong University.*****************************************************************************//*** constructors and destructor*/
template
BFGS::BFGS() : LineSearch()
{
}template
BFGS::~BFGS()
{
}/*** Finding the optimal solution. The default tolerance error and maximum* iteratin number are "tol&＃61;1.0e-6" and "maxItr&＃61;100", respectively.*/
template
void BFGS::optimize( Ftype &func, Vector &x0,Dtype tol, int maxItr )
{// initialize parameters.int k &＃61; 0,cnt &＃61; 0,N &＃61; x0.dim();Dtype ys,yHy,alpha;Vector d(N),s(N),y(N),v(N),Hy(N),gPrev(N);Matrix H &＃61; eye( N, Dtype(1.0) );Vector x(x0);Dtype fx &＃61; func(x);this->funcNum&＃43;&＃43;;Vector gnorm(maxItr);Vector g &＃61; func.grad(x);gnorm[k&＃43;&＃43;]&＃61; norm(g);while( ( gnorm(k) > tol ) && ( k getStep( func, x, d );// check flag for restartif( !this->success )if( norm(H-eye(N,Dtype(1.0))) funcNum&＃43;&＃43;;gPrev &＃61; g;g &＃61; func.grad(x);y &＃61; g - gPrev;Hy &＃61; H * y;ys &＃61; dotProd( y, s );yHy &＃61; dotProd( y, Hy );if( (ys tol )this->success &＃61; false;
}/*** Get the optimum point.*/
template
inline Vector BFGS::getOptValue() const
{return xOpt;
}/*** Get the norm of gradient in each iteration.*/
template
inline Vector BFGS::getGradNorm() const
{return gradNorm;
}/*** Get the minimum value of objective function.*/
template
inline Dtype BFGS::getFuncMin() const
{return fMin;
}/*** Get the iteration number.*/
template
inline int BFGS::getItrNum() const
{return gradNorm.dim()-1;
}

运行结果&＃xff1a;

The iterative number is: 7The number of function calculation is: 16The optimal value of x is: size: 2 by 1
-0.7071
0.0000The minimum value of f(x) is: -0.4289The gradient&＃39;s norm at x is: 0.0000Process returned 0 (0x0) execution time : 0.078 s
Press any key to continue.