Browse Source

adds code from port, passing first test

pull/728/head
rfellers 6 years ago
parent
commit
915c0b2c82
  1. 2
      src/Numerics.Tests/OptimizationTests/BrentMinimizerTests.cs
  2. 2
      src/Numerics.Tests/OptimizationTests/GoldenSectionMinimizerTests.cs
  3. 128
      src/Numerics/Optimization/BrentMinimizer.cs

2
src/Numerics.Tests/OptimizationTests/BrentMinimizerTests.cs

@ -52,7 +52,7 @@ namespace MathNet.Numerics.UnitTests.OptimizationTests
var algorithm = new BrentMinimizer(1e-5, 1000);
var f1 = new Func<double, double>(x => (x - 3) * (x - 3));
var obj = ObjectiveFunction.ScalarValue(f1);
var r1 = algorithm.FindMinimum(obj, -5, 5);
var r1 = algorithm.FindMinimum(obj, -2, 2);
Assert.That(Math.Abs(r1.MinimizingPoint - 3.0), Is.LessThan(1e-4));
}

2
src/Numerics.Tests/OptimizationTests/GoldenSectionMinimizerTests.cs

@ -53,7 +53,7 @@ namespace MathNet.Numerics.UnitTests.OptimizationTests
var algorithm = new GoldenSectionMinimizer(1e-5, 1000);
var f1 = new Func<double, double>(x => (x - 3)*(x - 3));
var obj = ObjectiveFunction.ScalarValue(f1);
var r1 = algorithm.FindMinimum(obj, -5, 5);
var r1 = algorithm.FindMinimum(obj, -1, 1);
Assert.That(Math.Abs(r1.MinimizingPoint - 3.0), Is.LessThan(1e-4));
}

128
src/Numerics/Optimization/BrentMinimizer.cs

@ -27,6 +27,7 @@
// OTHER DEALINGS IN THE SOFTWARE.
// </copyright>
using MathNet.Numerics.Optimization.ObjectiveFunctions;
using System;
namespace MathNet.Numerics.Optimization
@ -53,9 +54,132 @@ namespace MathNet.Numerics.Optimization
return Minimum(objective, lowerBound, upperBound, XTolerance, MaximumIterations, MaximumExpansionSteps, LowerExpansionFactor, UpperExpansionFactor);
}
public static ScalarMinimizationResult Minimum(IScalarObjectiveFunction objective, double lowerBound, double upperBound, double xTolerance = 1e-5, int maxIterations = 1000, int maxExpansionSteps = 10, double lowerExpansionFactor = 2.0, double upperExpansionFactor = 2.0)
public static ScalarMinimizationResult Minimum(IScalarObjectiveFunction objective, double lowerBound, double upperBound, double xTolerance = 1e-5,
int maxIterations = 1000, int maxExpansionSteps = 10, double lowerExpansionFactor = 2.0, double upperExpansionFactor = 2.0)
{
return null;
int maxfun = maxIterations;
if (lowerBound > upperBound)
throw new OptimizationException("Lower bound must be lower than upper bound.");
double sqrt_eps = Math.Sqrt(2.2e-16);
// This is not the golden_mean, but golden angle. Not sure why.
// https://en.wikipedia.org/wiki/Golden_angle
double golden_angle = 0.5 * (3.0 - Math.Sqrt(5.0));
double a = lowerBound;
double b = upperBound;
double fulc = a + golden_angle * (b - a);
double nfc = fulc, xf = fulc;
double rat = 0.0, e = 0.0;
double x = xf;
var evaluation = objective.Evaluate(x);
double fx = evaluation.Value;
int num = 1;
double fu = double.PositiveInfinity;
double ffulc = fx, fnfc = fx;
double xm = 0.5 * (a + b);
double tol1 = sqrt_eps * Math.Abs(xf) + xTolerance / 3.0;
double tol2 = 2.0 * tol1;
while (Math.Abs(xf - xm) > (tol2 - 0.5 * (b - a)))
{
bool golden = true;
// Check for parabolic fit
if (Math.Abs(e) > tol1)
{
golden = false;
double r = (xf - nfc) * (fx - ffulc);
double q = (xf - fulc) * (fx - fnfc);
double p = (xf - fulc) * q - (xf - nfc) * r;
q = 2.0 * (q - r);
if (q > 0.0)
p = -p;
q = Math.Abs(q);
r = e;
e = rat;
// Check for acceptability of parabola
if ((Math.Abs(p) < Math.Abs(0.5 * q * r)) && (p > q * (a - xf)) && (p < q * (b - xf)))
{
rat = (p + 0.0) / q;
x = xf + rat;
if (((x - a) < tol2) || ((b - x) < tol2))
{
int si_2 = Math.Sign(xm - xf) + ((xm - xf) == 0 ? 1 : 0);
rat = tol1 * si_2;
}
}
else // do a golden-section step
golden = true;
}
if (golden) // do a golden-section step
{
if (xf >= xm)
e = a - xf;
else
e = b - xf;
rat = golden_angle * e;
}
int si = Math.Sign(rat) + (rat == 0 ? 1 : 0);
x = xf + si * Math.Max(Math.Abs(rat), tol1);
evaluation = objective.Evaluate(x);
fu = evaluation.Value;
num += 1;
if (fu <= fx)
{
if (x >= xf)
a = xf;
else
b = xf;
fulc = nfc; ffulc = fnfc;
nfc = xf; fnfc = fx;
xf = x; fx = fu;
}
else
{
if (x < xf)
a = x;
else
b = x;
if ((fu <= fnfc) || (nfc == xf))
{
fulc = nfc; ffulc = fnfc;
nfc = x; fnfc = fu;
}
else if ((fu <= ffulc) || (fulc == xf) || (fulc == nfc))
{
fulc = x; ffulc = fu;
}
}
xm = 0.5 * (a + b);
tol1 = sqrt_eps * Math.Abs(xf) + xTolerance / 3.0;
tol2 = 2.0 * tol1;
if (num >= maxfun)
break;
}
var exitCondition = ExitCondition.BoundTolerance;
if (num >= maxfun)
exitCondition = ExitCondition.ExceedIterations;
else if (double.IsNaN(xf) || double.IsNaN(fx) || double.IsNaN(fu))
exitCondition = ExitCondition.InvalidValues;
return new ScalarMinimizationResult(new ScalarValueObjectiveFunctionEvaluation(xf, fx), num, exitCondition);
}
static void ValueChecker(double value, double point)

Loading…
Cancel
Save