In this post I develop the code for the interval halving region elimination class.
Interval halving assumes that the function is unimodal and correctly bracketed. This means that as one moves from the left bracket to the right bracket, the function continuously decreases until it hits a minimum point, and after that point, the function continuously increases up to the right bracket.
An example of a correctly bracketed unimodal function is (x-2)2 between 1 and 3.
An example of a function that is not a correctly bracketed unimodal function would be (x-2)2(x-4)2 between 1 and 5. A graph of this function is shown below. From 1 to 2, the function decreases and from 2 to 3 it increases, however, after 3, which is before the right bracket value of 5, the function starts decreasing again. The decrease after 3 prevents it from being a unimodal function within the bracket of 1 to 5. However, the function is a correctly bracketed unimodal function between 1 and 3, since it only decreases between 1 and 2 and only increases between 2 and 3.
Interval halving region elimination is a fairly simple algorithm. The initial interval is divided up into four equal length sections consisting of five points. These five points are the two boundary points passed in and the points between them at one fourth, two fourths, and three fourths of the distance between the boundary points.
I will refer to the left most point as x0, the point one fourth the distance from the left as x1, the midpoint of the interval as x2, the point one fourth the distance from the right bracket as x3, and the right bracket point as x4.
At each of the five points, the function is evaluated. Three cases can then result.
If the function value at x1 is less than the function value at the midpoint, x2, then the minimum lies between x0 and x2. In that case, the midpoint, x2, becomes the new right bracket, x0 remains the new left bracket, and x1 becomes the midpoint of the new interval. New values for x1 and x3 are then computed.
If the function value at x3 is less than the function value at the midpoint, x2, then the minimum lies between x2 and x4. In that case, the midpoint, x2, becomes the new left bracket, x4 remains the new right bracket, and x3 becomes the midpoint of the new interval. New values for x1 and x3 are then computed.
The third case that could occur involves the midpoint being less than either of the quarter points. In this case, the minimum lies between the two quarter points and x1 becomes the new left bracket, x3 becomes the new right bracket, and x2 remains the interval midpoint. New values for x1 and x3 are then computed.
In each of the above cases, half of the interval is eliminated in each iteration.
The full code for the interval halving region elimination class is shown at the bottom of this post.
In the Loaded event of the Window I added the code below.
IntervalHalvingRegionElimination re =
new IntervalHalvingRegionElimination();
SearchPointCollection col =
new SearchPointCollection();
col.Add(
new SearchPoint(1.0));
col.Add(
new SearchPoint(3.0));
Func<double, double> function = (x) =>
{
return (x - 2) * (x - 2);
};
RegionEliminationResults res = re.RegionEliminate(function, col, 1e-5);
Running this code gives a result of 2.0, which is the result expected.
The interval halving region elimination class is shown below. The code for the region elimination base classes have been posted previously.
public class IntervalHalvingRegionElimination : RegionEliminationBase
{
public IntervalHalvingRegionElimination()
:
base()
{
}
public override RegionEliminationResults RegionEliminate(Func<double, double> function, SearchPointCollection startPoint, double eliminationFraction)
{
Debug.Assert(startPoint.Count >= 2);
Cancel = false;
double result = 0.0;
CalculateMissingValues(startPoint, function);
SearchPointCollection searchPoints = new SearchPointCollection();
searchPoints.MaxSize = 5;
SearchPoint leftBracket = startPoint[0];
SearchPoint rightBracket = startPoint[startPoint.Count - 1];
searchPoints.Add(new SearchPoint(leftBracket.X, leftBracket.Y));
double length = rightBracket.X - leftBracket.X;
double stopRegionSize = length * eliminationFraction;
SearchPoint newPoint;
newPoint = new SearchPoint();
newPoint.X = leftBracket.X + length / 4.0;
newPoint.Y = function(newPoint.X);
searchPoints.Add(newPoint);
newPoint = new SearchPoint();
newPoint.X = leftBracket.X + length / 2.0;
newPoint.Y = function(newPoint.X);
searchPoints.Add(newPoint);
newPoint = new SearchPoint();
newPoint.X = leftBracket.X + (3.0 / 4.0) * length;
newPoint.Y = function(newPoint.X);
searchPoints.Add(newPoint);
searchPoints.Add(new SearchPoint(rightBracket.X, rightBracket.Y));
while (!Cancel)
{
if (searchPoints[1].Y < searchPoints[2].Y)
{
searchPoints[4] = searchPoints[2];
searchPoints[2] = searchPoints[1];
length = searchPoints[4].X - searchPoints[0].X;
searchPoints[1] = new SearchPoint();
searchPoints[1].X = searchPoints[0].X + length / 4.0;
searchPoints[1].Y = function(searchPoints[1].X);
searchPoints[3] = new SearchPoint();
searchPoints[3].X = searchPoints[0].X + (3.0 / 4.0) * length;
searchPoints[3].Y = function(searchPoints[3].X);
}
else if (searchPoints[3].Y < searchPoints[2].Y)
{
searchPoints[0] = searchPoints[2];
searchPoints[2] = searchPoints[3];
length = searchPoints[4].X - searchPoints[0].X;
searchPoints[1] = new SearchPoint();
searchPoints[1].X = searchPoints[0].X + length / 4.0;
searchPoints[1].Y = function(searchPoints[1].X);
searchPoints[3] = new SearchPoint();
searchPoints[3].X = searchPoints[0].X + (3.0 / 4.0) * length;
searchPoints[3].Y = function(searchPoints[3].X);
}
else
{
searchPoints[0] = searchPoints[1];
searchPoints[4] = searchPoints[3];
length = searchPoints[4].X - searchPoints[0].X;
searchPoints[1] = new SearchPoint();
searchPoints[1].X = searchPoints[0].X + length / 4.0;
searchPoints[1].Y = function(searchPoints[1].X);
searchPoints[3] = new SearchPoint();
searchPoints[3].X = searchPoints[0].X + (3.0 / 4.0) * length;
searchPoints[3].Y = function(searchPoints[3].X);
}
if (length < stopRegionSize)
{
Cancel = true;
}
}
result = (searchPoints[0].X + searchPoints[4].X) / 2.0;
return new RegionEliminationResults(result, function(result));
}
}
Subscribe