Author: tn
Date: Sun Sep 15 19:42:24 2013
New Revision: 1523495
URL: http://svn.apache.org/r1523495
Log:
[MATH-842] Fix Blands rule by applying it also to pivot column selection, add
getter/setter for cycle prevention, reduce default cut-off threshold to 1e-10.
Modified:
commons/proper/math/trunk/src/main/java/org/apache/commons/math3/optimization/linear/SimplexSolver.java
commons/proper/math/trunk/src/main/java/org/apache/commons/math3/optimization/linear/SimplexTableau.java
commons/proper/math/trunk/src/test/java/org/apache/commons/math3/optimization/linear/SimplexSolverTest.java
Modified:
commons/proper/math/trunk/src/main/java/org/apache/commons/math3/optimization/linear/SimplexSolver.java
URL:
http://svn.apache.org/viewvc/commons/proper/math/trunk/src/main/java/org/apache/commons/math3/optimization/linear/SimplexSolver.java?rev=1523495&r1=1523494&r2=1523495&view=diff
==============================================================================
---
commons/proper/math/trunk/src/main/java/org/apache/commons/math3/optimization/linear/SimplexSolver.java
(original)
+++
commons/proper/math/trunk/src/main/java/org/apache/commons/math3/optimization/linear/SimplexSolver.java
Sun Sep 15 19:42:24 2013
@@ -47,6 +47,9 @@ public class SimplexSolver extends Abstr
/** Amount of error to accept in floating point comparisons (as ulps). */
private final int maxUlps;
+ /** Prevent cycles via Bland's rule. */
+ private boolean cyclePrevention;
+
/**
* Build a simplex solver with default settings.
*/
@@ -62,6 +65,33 @@ public class SimplexSolver extends Abstr
public SimplexSolver(final double epsilon, final int maxUlps) {
this.epsilon = epsilon;
this.maxUlps = maxUlps;
+ this.cyclePrevention = true;
+ }
+
+ /**
+ * Enables/disables cycle prevention using Bland's rule.
+ * <p>
+ * In case of a degeneracy the simplex algorithm might cycle. This can be
prevented
+ * by using Bland's pivot selection rule. The drawback of this rule is
that it may result
+ * in pivot choices that do not significantly improve the objective
function value,
+ * thus increasing the number of required iterations.
+ * <p>
+ * By default, the cycle prevention is enabled.
+ *
+ * @param cyclePrevention if cycle prevention shall be used
+ *
+ * @see <a
href="http://www.stanford.edu/class/msande310/blandrule.pdf";>Bland's rule</a>
+ */
+ public void setCyclePrevention(boolean cyclePrevention) {
+ this.cyclePrevention = cyclePrevention;
+ }
+
+ /**
+ * Returns whether cycle prevention (using Bland's rule) is enabled.
+ * @return {@code true} if cycle prevention is enable, {@code false}
otherwise
+ */
+ public boolean getCyclePrevention() {
+ return cyclePrevention;
}
/**
@@ -79,14 +109,42 @@ public class SimplexSolver extends Abstr
if (entry < minValue) {
minValue = entry;
minPos = i;
+
+ // Bland's rule: chose the entering column with the lowest
index
+ if (cyclePrevention && isValidPivotColumn(tableau, i)) {
+ break;
+ }
}
}
return minPos;
}
/**
+ * Checks whether the given column is valid pivot column, i.e. will result
+ * in a valid pivot row.
+ * <p>
+ * When applying Bland's rule to select the pivot column, it may happen
that
+ * there is no corresponding pivot row. This method will check if the
selected
+ * pivot column will return a valid pivot row.
+ *
+ * @param tableau simplex tableau for the problem
+ * @param col the column to test
+ * @return {@code true} if the pivot column is valid, {@code false}
otherwise
+ */
+ private boolean isValidPivotColumn(SimplexTableau tableau, int col) {
+ for (int i = tableau.getNumObjectiveFunctions(); i <
tableau.getHeight(); i++) {
+ final double entry = tableau.getEntry(i, col);
+
+ if (Precision.compareTo(entry, 0d, maxUlps) > 0) {
+ return true;
+ }
+ }
+ return false;
+ }
+
+ /**
* Returns the row with the minimum ratio as given by the minimum ratio
test (MRT).
- * @param tableau simple tableau for the problem
+ * @param tableau simplex tableau for the problem
* @param col the column to test the ratio of. See {@link
#getPivotColumn(SimplexTableau)}
* @return row with the minimum ratio
*/
@@ -136,11 +194,7 @@ public class SimplexSolver extends Abstr
//
// see http://www.stanford.edu/class/msande310/blandrule.pdf
// see http://en.wikipedia.org/wiki/Bland%27s_rule (not equivalent
to the above paper)
- //
- // Additional heuristic: if we did not get a solution after half
of maxIterations
- // revert to the simple case of just
returning the top-most row
- // This heuristic is based on empirical data gathered while
investigating MATH-828.
- if (getIterations() < getMaxIterations() / 2) {
+ if (cyclePrevention) {
Integer minRow = null;
int minIndex = tableau.getWidth();
final int varStart = tableau.getNumObjectiveFunctions();
Modified:
commons/proper/math/trunk/src/main/java/org/apache/commons/math3/optimization/linear/SimplexTableau.java
URL:
http://svn.apache.org/viewvc/commons/proper/math/trunk/src/main/java/org/apache/commons/math3/optimization/linear/SimplexTableau.java?rev=1523495&r1=1523494&r2=1523495&view=diff
==============================================================================
---
commons/proper/math/trunk/src/main/java/org/apache/commons/math3/optimization/linear/SimplexTableau.java
(original)
+++
commons/proper/math/trunk/src/main/java/org/apache/commons/math3/optimization/linear/SimplexTableau.java
Sun Sep 15 19:42:24 2013
@@ -73,7 +73,7 @@ class SimplexTableau implements Serializ
private static final int DEFAULT_ULPS = 10;
/** The cut-off threshold to zero-out entries. */
- private static final double CUTOFF_THRESHOLD = 1e-12;
+ private static final double CUTOFF_THRESHOLD = 1e-10;
/** Serializable version identifier. */
private static final long serialVersionUID = -1369660067587938365L;
Modified:
commons/proper/math/trunk/src/test/java/org/apache/commons/math3/optimization/linear/SimplexSolverTest.java
URL:
http://svn.apache.org/viewvc/commons/proper/math/trunk/src/test/java/org/apache/commons/math3/optimization/linear/SimplexSolverTest.java?rev=1523495&r1=1523494&r2=1523495&view=diff
==============================================================================
---
commons/proper/math/trunk/src/test/java/org/apache/commons/math3/optimization/linear/SimplexSolverTest.java
(original)
+++
commons/proper/math/trunk/src/test/java/org/apache/commons/math3/optimization/linear/SimplexSolverTest.java
Sun Sep 15 19:42:24 2013
@@ -31,6 +31,31 @@ import org.junit.Test;
public class SimplexSolverTest {
@Test
+ public void testMath842Cycle() {
+ // from http://www.math.toronto.edu/mpugh/Teaching/APM236_04/bland
+ // maximize 10 x1 - 57 x2 - 9 x3 - 24 x4
+ // subject to
+ // 1/2 x1 - 11/2 x2 - 5/2 x3 + 9 x4 <= 0
+ // 1/2 x1 - 3/2 x2 - 1/2 x3 + x4 <= 0
+ // x1 <= 1
+ // x1,x2,x3,x4 >= 0
+
+ LinearObjectiveFunction f = new LinearObjectiveFunction(new double[] {
10, -57, -9, -24}, 0);
+
+ ArrayList <LinearConstraint>constraints = new
ArrayList<LinearConstraint>();
+
+ constraints.add(new LinearConstraint(new double[] {0.5, -5.5, -2.5,
9}, Relationship.LEQ, 0));
+ constraints.add(new LinearConstraint(new double[] {0.5, -1.5, -0.5,
1}, Relationship.LEQ, 0));
+ constraints.add(new LinearConstraint(new double[] { 1, 0, 0,
0}, Relationship.LEQ, 1));
+
+ double epsilon = 1e-6;
+ SimplexSolver solver = new SimplexSolver();
+ PointValuePair solution = solver.optimize(f, constraints,
GoalType.MAXIMIZE, true);
+ Assert.assertEquals(1.0d, solution.getValue(), epsilon);
+ Assert.assertTrue(validSolution(solution, constraints, epsilon));
+ }
+
+ @Test
public void testMath828() {
LinearObjectiveFunction f = new LinearObjectiveFunction(
new double[] { 1.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0,
0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0}, 0.0);
