Separate implementation of secant solver from bracketing solvers.

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

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

@@ -1,71 +0,0 @@
-/*
- * Licensed to the Apache Software Foundation (ASF) under one or more
- * contributor license agreements.  See the NOTICE file distributed with
- * this work for additional information regarding copyright ownership.
- * The ASF licenses this file to You under the Apache License, Version 2.0
- * (the "License"); you may not use this file except in compliance with
- * the License.  You may obtain a copy of the License at
- *
- *      http://www.apache.org/licenses/LICENSE-2.0
- *
- * Unless required by applicable law or agreed to in writing, software
- * distributed under the License is distributed on an "AS IS" BASIS,
- * WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied.
- * See the License for the specific language governing permissions and
- * limitations under the License.
- */
-
-package org.apache.commons.math.analysis.solvers;
-
-/**
- * Base class for <em>Secant</em> methods that guarantee convergence
- * by maintaining a {@link BracketedUnivariateRealSolver bracketed solution}.
- *
- * @since 3.0
- * @version $Id$
- */
-public class BaseBracketedSecantSolver extends BaseSecantSolver
-    implements BracketedUnivariateRealSolver {
-    /**
-     * Construct a solver with default accuracy (1e-6).
-     *
-     * @param method Method.
-     */
-    protected BaseBracketedSecantSolver(Method method) {
-        super(DEFAULT_ABSOLUTE_ACCURACY, method);
-    }
-
-    /**
-     * Construct a solver.
-     *
-     * @param absoluteAccuracy absolute accuracy
-     * @param method Method.
-     */
-    protected BaseBracketedSecantSolver(final double absoluteAccuracy,
-                                        Method method) {
-        super(absoluteAccuracy, method);
-    }
-
-    /**
-     * Construct a solver.
-     *
-     * @param relativeAccuracy relative accuracy
-     * @param absoluteAccuracy absolute accuracy
-     * @param method Method.
-     */
-    protected BaseBracketedSecantSolver(final double relativeAccuracy,
-                                        final double absoluteAccuracy,
-                                        Method method) {
-        super(relativeAccuracy, absoluteAccuracy, method);
-    }
-
-    /** {@inheritDoc} */
-    public AllowedSolutions getAllowedSolutions() {
-        return allowedSolutions;
-    }
-
-    /** {@inheritDoc} */
-    public void setAllowedSolutions(final AllowedSolutions allowedSolutions) {
-        this.allowedSolutions = allowedSolutions;
-    }
-}

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

@@ -18,14 +18,14 @@
 package org.apache.commons.math.analysis.solvers;
 
 import org.apache.commons.math.util.FastMath;
+import org.apache.commons.math.analysis.UnivariateRealFunction;
 import org.apache.commons.math.exception.MathInternalError;
 
 /**
- * Base class for all <em>Secant</em>-based methods for root-finding
+ * Base class for all bracketing <em>Secant</em>-based methods for root-finding
  * (approximating a zero of a univariate real function).
  *
- * <p>Implementation of the {@link SecantSolver <em>Secant</em>},
+ * <p>Implementation of the {@link RegulaFalsiSolver <em>Regula Falsi</em>}, and
- * {@link RegulaFalsiSolver <em>Regula Falsi</em>}, and
  * {@link IllinoisSolver <em>Illinois</em>} methods is based on the
  * following article: M. Dowell and P. Jarratt,
  * <em>A modified regula falsi method for computing the root of an
@@ -38,14 +38,23 @@ import org.apache.commons.math.exception.MathInternalError;
  * BIT Numerical Mathematics, volume 12, number 4, pages 503-508, Springer,
  * 1972.</p>
  *
+ * <p>The  {@link SecantSolver <em>secant<em>} method is <em>not</emp> a
+ * bracketing method so it is not implemented here. It has a separate
+ * implementation.</p>
+ *
  * @since 3.0
  * @version $Id$
  */
-public abstract class BaseSecantSolver extends AbstractUnivariateRealSolver {
+public abstract class BaseSecantSolver
+    extends AbstractUnivariateRealSolver
+    implements BracketedUnivariateRealSolver<UnivariateRealFunction> {
+
     /** Default absolute accuracy. */
     protected static final double DEFAULT_ABSOLUTE_ACCURACY = 1e-6;
+
     /** The kinds of solutions that the algorithm may accept. */
-    protected AllowedSolutions allowedSolutions = AllowedSolutions.EITHER_SIDE;
+    private AllowedSolutions allowedSolutions;
+
     /** The <em>Secant</em>-based root-finding method to use. */
     private final Method method;
 
@@ -57,6 +66,7 @@ public abstract class BaseSecantSolver extends AbstractUnivariateRealSolver {
      */
     protected BaseSecantSolver(final double absoluteAccuracy, final Method method) {
         super(absoluteAccuracy);
+        this.allowedSolutions = AllowedSolutions.ANY_SIDE;
         this.method = method;
     }
 
@@ -71,11 +81,33 @@ public abstract class BaseSecantSolver extends AbstractUnivariateRealSolver {
                                final double absoluteAccuracy,
                                final Method method) {
         super(relativeAccuracy, absoluteAccuracy);
+        this.allowedSolutions = AllowedSolutions.ANY_SIDE;
         this.method = method;
     }
 
+    /** {@inheritDoc} */
+    public double solve(final int maxEval, final UnivariateRealFunction f,
+                        final double min, final double max,
+                        final AllowedSolutions allowedSolutions) {
+        return solve(maxEval, f, min, max, min + 0.5 * (max - min), allowedSolutions);
+    }
+
+    /** {@inheritDoc} */
+    public double solve(final int maxEval, final UnivariateRealFunction f,
+                        final double min, final double max, final double startValue,
+                        final AllowedSolutions allowedSolutions) {
+        this.allowedSolutions = allowedSolutions;
+        return super.solve(maxEval, f, min, max, startValue);
+    }
+
     /** {@inheritDoc} */
     @Override
+    public double solve(final int maxEval, final UnivariateRealFunction f,
+                        final double min, final double max, final double startValue) {
+        return solve(maxEval, f, min, max, startValue, AllowedSolutions.ANY_SIDE);
+    }
+
+    /** {@inheritDoc} */
     protected final double doSolve() {
         // Get initial solution
         double x0 = getMin();
@@ -102,8 +134,7 @@ public abstract class BaseSecantSolver extends AbstractUnivariateRealSolver {
         final double rtol = getRelativeAccuracy();
 
         // Keep track of inverted intervals, meaning that the left bound is
-        // larger than the right bound. Not used for the original Secant
+        // larger than the right bound.
-        // method.
         boolean inverted = false;
 
         // Keep finding better approximations.
@@ -120,12 +151,7 @@ public abstract class BaseSecantSolver extends AbstractUnivariateRealSolver {
             }
 
             // Update the bounds with the new approximation.
-            if (method == Method.SECANT) {
+            if (f1 * fx < 0) {
-                x0 = x1;
-                f0 = f1;
-                x1 = x;
-                f1 = fx;
-            } else if (f1 * fx < 0) {
                 // We had [x0..x1]. We update it to [x1, x]. Note that the
                 // value of x1 has switched to the other bound, thus inverting
                 // the interval.
@@ -151,7 +177,7 @@ public abstract class BaseSecantSolver extends AbstractUnivariateRealSolver {
             // the root than we already are.
             if (FastMath.abs(f1) <= ftol) {
                 switch (allowedSolutions) {
-                case EITHER_SIDE:
+                case ANY_SIDE:
                     return x1;
                 case LEFT_SIDE:
                     if (inverted) {
@@ -183,7 +209,7 @@ public abstract class BaseSecantSolver extends AbstractUnivariateRealSolver {
             if (FastMath.abs(x1 - x0) < FastMath.max(rtol * FastMath.abs(x1),
                                                      atol)) {
                 switch (allowedSolutions) {
-                case EITHER_SIDE:
+                case ANY_SIDE:
                     return x1;
                 case LEFT_SIDE:
                     return inverted ? x1 : x0;
@@ -202,8 +228,6 @@ public abstract class BaseSecantSolver extends AbstractUnivariateRealSolver {
 
     /** <em>Secant</em>-based root-finding methods. */
     protected enum Method {
-        /** The original {@link SecantSolver <em>Secant</em>} method. */
-        SECANT,
 
         /**
          * The {@link RegulaFalsiSolver <em>Regula Falsi</em>} or
@@ -215,6 +239,7 @@ public abstract class BaseSecantSolver extends AbstractUnivariateRealSolver {
         ILLINOIS,
 
         /** The {@link PegasusSolver <em>Pegasus</em>} method. */
-        PEGASUS,
+        PEGASUS;
+
     }
 }

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

@@ -17,6 +17,8 @@
 
 package org.apache.commons.math.analysis.solvers;
 
+import org.apache.commons.math.util.FastMath;
+
 /**
  * Implements the <em>Secant</em> method for root-finding (approximating a
  * zero of a univariate real function). The solution that is maintained is
@@ -36,10 +38,14 @@ package org.apache.commons.math.analysis.solvers;
  *
  * @version $Id$
  */
-public class SecantSolver extends BaseSecantSolver {
+public class SecantSolver extends AbstractUnivariateRealSolver {
+
+    /** Default absolute accuracy. */
+    protected static final double DEFAULT_ABSOLUTE_ACCURACY = 1e-6;
+
     /** Construct a solver with default accuracy (1e-6). */
     public SecantSolver() {
-        super(DEFAULT_ABSOLUTE_ACCURACY, Method.SECANT);
+        super(DEFAULT_ABSOLUTE_ACCURACY);
     }
 
     /**
@@ -48,7 +54,7 @@ public class SecantSolver extends BaseSecantSolver {
      * @param absoluteAccuracy absolute accuracy
      */
     public SecantSolver(final double absoluteAccuracy) {
-        super(absoluteAccuracy, Method.SECANT);
+        super(absoluteAccuracy);
     }
 
     /**
@@ -59,6 +65,69 @@ public class SecantSolver extends BaseSecantSolver {
      */
     public SecantSolver(final double relativeAccuracy,
                         final double absoluteAccuracy) {
-        super(relativeAccuracy, absoluteAccuracy, Method.SECANT);
+        super(relativeAccuracy, absoluteAccuracy);
     }
+
+    /** {@inheritDoc} */
+    @Override
+    protected final double doSolve() {
+        // Get initial solution
+        double x0 = getMin();
+        double x1 = getMax();
+        double f0 = computeObjectiveValue(x0);
+        double f1 = computeObjectiveValue(x1);
+
+        // If one of the bounds is the exact root, return it. Since these are
+        // not under-approximations or over-approximations, we can return them
+        // regardless of the allowed solutions.
+        if (f0 == 0.0) {
+            return x0;
+        }
+        if (f1 == 0.0) {
+            return x1;
+        }
+
+        // Verify bracketing of initial solution.
+        verifyBracketing(x0, x1);
+
+        // Get accuracies.
+        final double ftol = getFunctionValueAccuracy();
+        final double atol = getAbsoluteAccuracy();
+        final double rtol = getRelativeAccuracy();
+
+        // Keep finding better approximations.
+        while (true) {
+
+            // Calculate the next approximation.
+            final double x = x1 - ((f1 * (x1 - x0)) / (f1 - f0));
+            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,
+            // we can return it regardless of the allowed solutions.
+            if (fx == 0.0) {
+                return x;
+            }
+
+            // Update the bounds with the new approximation.
+            x0 = x1;
+            f0 = f1;
+            x1 = x;
+            f1 = fx;
+
+            // If the function value of the last approximation is too small,
+            // given the function value accuracy, then we can't get closer to
+            // the root than we already are.
+            if (FastMath.abs(f1) <= ftol) {
+                return x1;
+            }
+
+            // If the current interval is within the given accuracies, we
+            // are satisfied with the current approximation.
+            if (FastMath.abs(x1 - x0) < FastMath.max(rtol * FastMath.abs(x1), atol)) {
+                return x1;
+            }
+        }
+    }
+
 }