Browse Source

More tidying up of minimizers.

pull/173/head
joemoorhouse 13 years ago
parent
commit
5c0466968f
  1. 4
      src/Numerics/Optimization/BrentMinimizer.cs
  2. 57
      src/Numerics/Optimization/PowellMinimizer.cs
  3. 2
      src/UnitTests/OptimizationTests/FunctionMinimizationTests.cs

4
src/Numerics/Optimization/BrentMinimizer.cs

@ -11,7 +11,7 @@ namespace MathNet.Numerics.Optimization
public class BrentOptions public class BrentOptions
{ {
public int MaximumIterations = 1000; public int MaximumIterations = 1000;
public double FunctionTolerance = 1e-4; public double Tolerance = 1e-4;
} }
/// <summary> /// <summary>
@ -77,7 +77,7 @@ namespace MathNet.Numerics.Optimization
double tol1, tol2, xmid; double tol1, tol2, xmid;
double temp1, temp2, p; double temp1, temp2, p;
double u, fu, dx_temp; double u, fu, dx_temp;
tol1 = Options.FunctionTolerance * Math.Abs(x) + minimumTolerance; tol1 = Options.Tolerance * Math.Abs(x) + minimumTolerance;
tol2 = 2.0 * tol1; tol2 = 2.0 * tol1;
xmid = 0.5 * (a + b); xmid = 0.5 * (a + b);
if (Math.Abs(x - xmid) < (tol2 - 0.5 * (b - a))) // check for convergence if (Math.Abs(x - xmid) < (tol2 - 0.5 * (b - a))) // check for convergence

57
src/Numerics/Optimization/PowellMinimizer.cs

@ -12,8 +12,8 @@ namespace MathNet.Numerics.Optimization
{ {
public int? MaximumIterations = null; public int? MaximumIterations = null;
public int? MaximumFunctionCalls = null; public int? MaximumFunctionCalls = null;
public double PointTolerance = 1e-4; public double PointTolerance = 1e-4; // absolute
public double FunctionTolerance = 1e-4; public double FunctionTolerance = 1e-4; // relative
} }
public enum PowellConvergenceType { Success, MaxIterationsExceeded, MaxFunctionCallsExceeded }; public enum PowellConvergenceType { Success, MaxIterationsExceeded, MaxFunctionCallsExceeded };
@ -97,47 +97,46 @@ namespace MathNet.Numerics.Optimization
int maxFunctionCalls = (Options.MaximumFunctionCalls == null) ? n * 1000 : (int)Options.MaximumFunctionCalls; int maxFunctionCalls = (Options.MaximumFunctionCalls == null) ? n * 1000 : (int)Options.MaximumFunctionCalls;
// An array of n directions: // An array of n directions:
double[][] direc = new double[n][]; double[][] directionSet = new double[n][];
for (int i = 0; i < n; ++i) for (int i = 0; i < n; ++i)
{ {
direc[i] = new double[n]; directionSet[i] = new double[n];
direc[i][i] = 1.0; directionSet[i][i] = 1.0;
} }
double[] x = pInitialGuess; double[] x = pInitialGuess;
double[] x1 = (double[])x.Clone(); double[] x1 = (double[])x.Clone();
double[] x2 = new double[n];
brentMinimizer.Options.FunctionTolerance = Options.PointTolerance * 100; double[] direction1 = new double[n];
brentMinimizer.Options.Tolerance = Options.PointTolerance * 100;
fval = function(x); fval = function(x);
double[] x2 = new double[n];
double fx; double fx;
double[] direc1 = new double[n];
double[] xnew;
while (true) while (true)
{ {
fx = fval; fx = fval;
int bigind = 0; int bigIndex = 0;
double delta = 0.0; double delta = 0.0;
double fx2; double fx2;
// loop over all directions:
for (int i = 0; i < n; ++i) for (int i = 0; i < n; ++i)
{ {
direc1 = direc[i];
fx2 = fval; fx2 = fval;
// Do a linesearch with specified starting point and direction. // Do a linesearch along direction
direction = direc1; direction = directionSet[i];
startingPoint = x; startingPoint = x;
double u = brentMinimizer.Minimize(functionAlongLine); double u = brentMinimizer.Minimize(functionAlongLine);
fval = functionAlongLine(u); fval = functionAlongLine(u); // updates point
xnew = point;
for (int j = 0; j < n; ++j) x[j] = xnew[j]; for (int j = 0; j < n; ++j) x[j] = point[j];
if ((fx2 - fval) > delta) if ((fx2 - fval) > delta)
{ {
delta = fx2 - fval; delta = fx2 - fval;
bigind = i; bigIndex = i; // record direction index of the largest decrease seen
} }
} }
iterations++; iterations++;
@ -145,11 +144,10 @@ namespace MathNet.Numerics.Optimization
if (functionCalls >= maxFunctionCalls) break; if (functionCalls >= maxFunctionCalls) break;
if (iterations >= maxIterations) break; if (iterations >= maxIterations) break;
// Construct the extrapolated point // Construct the extrapolated point
direc1 = new double[n];
for (int i = 0; i < n; ++i) for (int i = 0; i < n; ++i)
{ {
direc1[i] = x[i] - x1[i]; direction1[i] = x[i] - x1[i];
x2[i] = 2.0 * x[i] - x1[i]; x2[i] = 2.0 * x[i] - x1[i];
x1[i] = x[i]; x1[i] = x[i];
} }
@ -164,21 +162,18 @@ namespace MathNet.Numerics.Optimization
t -= delta * temp * temp; t -= delta * temp * temp;
if (t < 0.0) if (t < 0.0)
{ {
direction = direc1; // Do a linesearch along direction
direction = direction1;
startingPoint = x; startingPoint = x;
double u = brentMinimizer.Minimize(functionAlongLine); double u = brentMinimizer.Minimize(functionAlongLine);
fval = functionAlongLine(u); fval = functionAlongLine(u); // updates point
xnew = point;
direc1 = new double[n];
for (int i = 0; i < n; ++i) for (int i = 0; i < n; ++i)
{ {
direc1[i] = xnew[i] - x[i]; directionSet[bigIndex][i] = directionSet[n - 1][i];
x[i] = xnew[i]; directionSet[n - 1][i] = point[i] - x[i];
x[i] = point[i];
} }
direc[bigind] = direc[n - 1];
direc[n - 1] = direc1;
} }
} }
} }

2
src/UnitTests/OptimizationTests/FunctionMinimizationTests.cs

@ -50,6 +50,8 @@ namespace MathNet.Numerics.UnitTests.OptimizationTests
var watch = new System.Diagnostics.Stopwatch(); watch.Start(); var watch = new System.Diagnostics.Stopwatch(); watch.Start();
double[] popt = null; double[] popt = null;
minimizer.Options.PointTolerance = 1e-8;
for (int i = 0; i < 1000; ++i) for (int i = 0; i < 1000; ++i)
{ {
popt = minimizer.CurveFit(xin, yin, function, new double[] { 1, 1 }); // 100, 0.75 popt = minimizer.CurveFit(xin, yin, function, new double[] { 1, 1 }); // 100, 0.75

Loading…
Cancel
Save