In this post I develop the DFP / Davidon Fletcher Powell class which is a quasi Newton numerical method that can be used for unconstrained optimization problems.
Most of the functionality of this class resides in the previously created quasi Newton base class, Quasi Newton Base, from which the DFP class derives. The primary thing that needs to be done for this class after inheriting from Quasi Newton Base is to override the Update Inverse Hessian virtual function and use it to compute the approximate inverse Hessian using the DFP update formula.
The DFP inverse Hessian update formula is shown below.
This equation is implemented in the overridden virtual function Update Inverse Hessian.
The code for the DFP class is given below.
public class DFP : QuasiNewtonBase
{
public DFP(Func<
double> function, Parameter[] optimizationParameters,
double inverseHessianResetMultiple,
int numberOfDerivativePoints) :
base(function, optimizationParameters, inverseHessianResetMultiple, numberOfDerivativePoints)
{
}
public DFP(Func<double> function, Parameter[] optimizationParameters)
: this(function, optimizationParameters, 2.0, 3)
{
}
protected override void UpdateInverseHessian()
{
int numberOfParameters = _solvedForParameters.Length;
double[] positionDelta = new double[numberOfParameters];
double[] gradientDelta = new double[numberOfParameters];
for (int i = 0; i < numberOfParameters; i++)
{
gradientDelta[i] = _currentGradient[i] - _lastGradient[i];
positionDelta[i] = _currentPosition[i] - _lastPosition[i];
}
Matrix XXt = ComputeQNProduct(positionDelta, positionDelta);
double YtX = ComputeInnerProduct(gradientDelta, positionDelta);
Matrix T1 = XXt / YtX;
Matrix YYt = ComputeQNProduct(gradientDelta, gradientDelta);
Matrix AYYtAt = _inverseHessian * YYt * _inverseHessian.Transpose();
double YtAY = ComputeInnerProduct(gradientDelta, _inverseHessian * gradientDelta);
Matrix T2 = AYYtAt / YtAY;
Matrix matInverseHessianNew = _inverseHessian + T1 - T2;
_inverseHessian = matInverseHessianNew;
}
}
Subscribe