Added BELOW_SIDE and ABOVE_SIDE in the possible allowed solutions for bracketing solvers.

git-svn-id: https://svn.apache.org/repos/asf/commons/proper/math/trunk@1142244 13f79535-47bb-0310-9956-ffa450edef68

src/main/java/org/apache/commons/math/analysis/solvers/AllowedSolutions.java

@@ -18,9 +18,10 @@
 package org.apache.commons.math.analysis.solvers;
 
 
-/** The kinds of solutions that a {@link UnivariateRealSolver (univariate real)
+/** The kinds of solutions that a {@link BracketedUnivariateRealSolver
- * root-finding algorithm} may accept as solutions. This basically controls
+ * (bracketed univariate real) root-finding algorithm} may accept as solutions.
- * whether or not under-approximations and over-approximations are allowed.
+ * This basically controls whether or not under-approximations and
+ * over-approximations are allowed.
  *
  * <p>If all solutions are accepted ({@link #EITHER_SIDE}), then the solution
  * that the root-finding algorithm returns for a given root may be equal to the
@@ -30,8 +31,8 @@ package org.apache.commons.math.analysis.solvers;
  * tolerances. In certain cases however, in particular for
  * {@link org.apache.commons.math.ode.events.EventHandler state events} of
  * {@link org.apache.commons.math.ode.ODEIntegrator ODE solvers}, it
- * may be necessary to guarantee that a solution is returned that does not
+ * may be necessary to guarantee that a solution is returned that lies on a
- * under-approximate the solution.</p>
+ * specific side the solution.</p>
  *
  * @see BracketedUnivariateRealSolver
  * @since 3.0
@@ -40,23 +41,36 @@ package org.apache.commons.math.analysis.solvers;
 public enum AllowedSolutions {
     /** There are no additional side restriction on the solutions for
      * root-finding. That is, both under-approximations and over-approximations
-     * are allowed. So, if a function f(x) has a root at x = y, then the
+     * are allowed. So, if a function f(x) has a root at x = x0, then the
-     * root-finding result s may be smaller than y, equal to y, or greater
+     * root-finding result s may be smaller than x0, equal to x0, or greater
-     * than y.
+     * than x0.
      */
     EITHER_SIDE,
 
     /** Only solutions that are less than or equal to the actual root are
      * acceptable as solutions for root-finding. In other words,
      * over-approximations are not allowed. So, if a function f(x) has a root
-     * at x = y, then the root-finding result s must satisfy s &lt;= y.
+     * at x = x0, then the root-finding result s must satisfy s &lt;= x0.
      */
     LEFT_SIDE,
 
     /** Only solutions that are greater than or equal to the actual root are
      * acceptable as solutions for root-finding. In other words,
      * under-approximations are not allowed. So, if a function f(x) has a root
-     * at x = y, then the root-finding result s must satisfy s &gt;= y.
+     * at x = x0, then the root-finding result s must satisfy s &gt;= x0.
      */
-    RIGHT_SIDE;
+    RIGHT_SIDE,
+
+    /** Only solutions for which values are less than or equal to zero are
+     * acceptable as solutions for root-finding. So, if a function f(x) has
+     * a root at x = x0, then the root-finding result s must satisfy f(s) &lt;= 0.
+     */
+    BELOW_SIDE,
+
+    /** Only solutions for which values are greater than or equal to zero are
+     * acceptable as solutions for root-finding. So, if a function f(x) has
+     * a root at x = x0, then the root-finding result s must satisfy f(s) &gt;= 0.
+     */
+    ABOVE_SIDE;
+
 }

src/main/java/org/apache/commons/math/analysis/solvers/BaseSecantSolver.java

@@ -101,10 +101,6 @@ public abstract class BaseSecantSolver extends AbstractUnivariateRealSolver {
         final double atol = getAbsoluteAccuracy();
         final double rtol = getRelativeAccuracy();
 
-        // Variables to hold new bounds.
-        double x;
-        double fx;
-
         // Keep track of inverted intervals, meaning that the left bound is
         // larger than the right bound. Not used for the original Secant
         // method.
@@ -113,8 +109,8 @@ public abstract class BaseSecantSolver extends AbstractUnivariateRealSolver {
         // Keep finding better approximations.
         while (true) {
             // Calculate the next approximation.
-            x = x1 - ((f1 * (x1 - x0)) / (f1 - f0));
+            final double x = x1 - ((f1 * (x1 - x0)) / (f1 - f0));
-            fx = computeObjectiveValue(x);
+            final double fx = computeObjectiveValue(x);
 
             // If the new approximation is the exact root, return it. Since
             // this is not an under-approximation or an over-approximation,
@@ -151,7 +147,7 @@ public abstract class BaseSecantSolver extends AbstractUnivariateRealSolver {
             }
 
             // If the function value of the last approximation is too small,
-            // given the function value accuracy, then we can't get close to
+            // given the function value accuracy, then we can't get closer to
             // the root than we already are.
             if (FastMath.abs(f1) <= ftol) {
                 switch (allowedSolutions) {
@@ -167,6 +163,16 @@ public abstract class BaseSecantSolver extends AbstractUnivariateRealSolver {
                         return x1;
                     }
                     break;
+                case BELOW_SIDE:
+                    if (f1 <= 0) {
+                        return x1;
+                    }
+                    break;
+                case ABOVE_SIDE:
+                    if (f1 >= 0) {
+                        return x1;
+                    }
+                    break;
                 default:
                     throw new MathInternalError();
                 }
@@ -183,6 +189,10 @@ public abstract class BaseSecantSolver extends AbstractUnivariateRealSolver {
                     return inverted ? x1 : x0;
                 case RIGHT_SIDE:
                     return inverted ? x0 : x1;
+                case BELOW_SIDE:
+                    return (f1 <= 0) ? x1 : x0;
+                case ABOVE_SIDE:
+                    return (f1 >= 0) ? x1 : x0;
                 default:
                     throw new MathInternalError();
                 }

src/site/xdoc/changes.xml

@@ -59,7 +59,8 @@ The <action> type attribute can be add,update,fix,remove.
       </action>
       <action dev="eran" type="add" issue="MATH-599" due-to="Dennis Hendriks">
         Modified "SecantSolver" to comply with the original algorithm. Added several
-        secant-based solvers.
+        secant-based solvers. Added a way to select the side of the root with bracketing
+        solvers.
       </action>
       <action dev="luc" type="fix" issue="MATH-600">
         Fixed javadoc for ODEIntegrator interface

src/test/java/org/apache/commons/math/analysis/solvers/BaseSecantSolverAbstractTest.java

@@ -195,4 +195,42 @@ public abstract class BaseSecantSolverAbstractTest {
             right += 0.3;
         }
     }
+    @Test
+    public void testSolutionBelowSide() {
+        UnivariateRealFunction f = new SinFunction();
+        UnivariateRealSolver solver = getSolver();
+        if (!(solver instanceof BracketedUnivariateRealSolver)) return;
+        ((BracketedUnivariateRealSolver)solver).setAllowedSolutions(AllowedSolutions.BELOW_SIDE);
+        double left = -1.5;
+        double right = 0.05;
+        for(int i = 0; i < 10; i++) {
+            // Test whether the allowed solutions are taken into account.
+            double solution = solver.solve(100, f, left, right);
+            Assert.assertTrue(f.value(solution) <= 0.0);
+
+            // Prepare for next test.
+            left -= 0.1;
+            right += 0.3;
+        }
+    }
+
+    @Test
+    public void testSolutionAboveSide() {
+        UnivariateRealFunction f = new SinFunction();
+        UnivariateRealSolver solver = getSolver();
+        if (!(solver instanceof BracketedUnivariateRealSolver)) return;
+        ((BracketedUnivariateRealSolver)solver).setAllowedSolutions(AllowedSolutions.ABOVE_SIDE);
+        double left = -1.5;
+        double right = 0.05;
+        for(int i = 0; i < 10; i++) {
+            // Test whether the allowed solutions are taken into account.
+            double solution = solver.solve(100, f, left, right);
+            Assert.assertTrue(f.value(solution) >= 0.0);
+
+            // Prepare for next test.
+            left -= 0.1;
+            right += 0.3;
+        }
+    }
+
 }