View Javadoc

    1   /*
    2    * Copyright 2003-2004 The Apache Software Foundation.
    3    *
    4    * Licensed under the Apache License, Version 2.0 (the "License");
    5    * you may not use this file except in compliance with the License.
    6    * You may obtain a copy of the License at
    7    *
    8    *      http://www.apache.org/licenses/LICENSE-2.0
    9    *
    10   * Unless required by applicable law or agreed to in writing, software
    11   * distributed under the License is distributed on an "AS IS" BASIS,
    12   * WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied.
    13   * See the License for the specific language governing permissions and
    14   * limitations under the License.
    15   */
    16  package org.apache.commons.math.analysis;
    17
    18  import org.apache.commons.math.FunctionEvaluationException;
    19  import org.apache.commons.math.ConvergenceException;
    20
    21  /***
    22   * Utility routines for {@link UnivariateRealSolver} objects.
    23   *
    24   * @version $Revision: 155427 $ $Date: 2005-02-26 06:11:52 -0700 (Sat, 26 Feb 2005) $
    25   */
    26  public class UnivariateRealSolverUtils {
    27      /***
    28       * Default constructor.
    29       */
    30      private UnivariateRealSolverUtils() {
    31          super();
    32      }
    33
    34      /*** Cached solver factory */
    35      private static UnivariateRealSolverFactory factory = null;
    36
    37      /***
    38       * Convenience method to find a zero of a univariate real function.  A default
    39       * solver is used.
    40       *
    41       * @param f the function.
    42       * @param x0 the lower bound for the interval.
    43       * @param x1 the upper bound for the interval.
    44       * @return a value where the function is zero.
    45       * @throws ConvergenceException if the iteration count was exceeded
    46       * @throws FunctionEvaluationException if an error occurs evaluating
    47       * the function
    48       * @throws IllegalArgumentException if f is null or the endpoints do not
    49       * specify a valid interval
    50       */
    51      public static double solve(UnivariateRealFunction f, double x0, double x1)
    52      throws ConvergenceException, FunctionEvaluationException {
    53          setup(f);
    54          return factory.newDefaultSolver(f).solve(x0, x1);
    55      }
    56
    57      /***
    58       * Convenience method to find a zero of a univariate real function.  A default
    59       * solver is used.
    60       *
    61       * @param f the function
    62       * @param x0 the lower bound for the interval
    63       * @param x1 the upper bound for the interval
    64       * @param absoluteAccuracy the accuracy to be used by the solver
    65       * @return a value where the function is zero
    66       * @throws ConvergenceException if the iteration count is exceeded
    67       * @throws FunctionEvaluationException if an error occurs evaluating the
    68       * function
    69       * @throws IllegalArgumentException if f is null, the endpoints do not
    70       * specify a valid interval, or the absoluteAccuracy is not valid for the
    71       * default solver
    72       */
    73      public static double solve(UnivariateRealFunction f, double x0, double x1,
    74              double absoluteAccuracy) throws ConvergenceException,
    75              FunctionEvaluationException {
    76
    77          setup(f);
    78          UnivariateRealSolver solver = factory.newDefaultSolver(f);
    79          solver.setAbsoluteAccuracy(absoluteAccuracy);
    80          return solver.solve(x0, x1);
    81      }
    82
    83      /***
    84       * This method attempts to find two values a and b satisfying <ul>
    85      * <li> <code> lowerBound <= a < initial < b <= upperBound</code> </li>
    86       * <li> <code> f(a) * f(b) < 0 </code></li>
    87       * </ul>
    88       * If f is continuous on <code>[a,b],</code> this means that <code>a</code>
    89       * and <code>b</code> bracket a root of f.
    90       * <p>
    91       * The algorithm starts by setting
    92       * <code>a := initial -1; b := initial +1,</code> examines the value of the
    93       * function at <code>a</code> and <code>b</code> and keeps moving
    94       * the endpoints out by one unit each time through a loop that terminates
    95       * when one of the following happens: <ul>
    96       * <li> <code> f(a) * f(b) < 0 </code> --  success!</li>
    97       * <li> <code> a = lower </code> and <code> b = upper</code>
    98       * -- ConvergenceException </li>
    99       * <li> <code> Integer.MAX_VALUE</code> iterations elapse
    100      * -- ConvergenceException </li>
    101      * </ul>
    102      * <p>
    103      * <strong>Note: </strong> this method can take
    104      * <code>Integer.MAX_VALUE</code> iterations to throw a
    105      * <code>ConvergenceException.</code>  Unless you are confident that there
    106      * is a root between <code>lowerBound</code> and <code>upperBound</code>
    107      * near <code>initial,</code> it is better to use
    108      * {@link #bracket(UnivariateRealFunction, double, double, double, int)},
    109      * explicitly specifying the maximum number of iterations.
    110      *
    111      * @param function the function
    112      * @param initial initial midpoint of interval being expanded to
    113      * bracket a root
    114      * @param lowerBound lower bound (a is never lower than this value)
    115      * @param upperBound upper bound (b never is greater than this
    116      * value)
    117      * @return a two element array holding {a, b}
    118      * @throws ConvergenceException if a root can not be bracketted
    119      * @throws FunctionEvaluationException if an error occurs evaluating the
    120      * function
    121      * @throws IllegalArgumentException if function is null, maximumIterations
    122      * is not positive, or initial is not between lowerBound and upperBound
    123      */
    124     public static double[] bracket(UnivariateRealFunction function,
    125             double initial, double lowerBound, double upperBound)
    126     throws ConvergenceException, FunctionEvaluationException {
    127         return bracket( function, initial, lowerBound, upperBound,
    128             Integer.MAX_VALUE ) ;
    129     }
    130
    131      /***
    132      * This method attempts to find two values a and b satisfying <ul>
    133      * <li> <code> lowerBound <= a < initial < b <= upperBound</code> </li>
    134      * <li> <code> f(a) * f(b) < 0 </code> </li>
    135      * </ul>
    136      * If f is continuous on <code>[a,b],</code> this means that <code>a</code>
    137      * and <code>b</code> bracket a root of f.
    138      * <p>
    139      * The algorithm starts by setting
    140      * <code>a := initial -1; b := initial +1,</code> examines the value of the
    141      * function at <code>a</code> and <code>b</code> and keeps moving
    142      * the endpoints out by one unit each time through a loop that terminates
    143      * when one of the following happens: <ul>
    144      * <li> <code> f(a) * f(b) < 0 </code> --  success!</li>
    145      * <li> <code> a = lower </code> and <code> b = upper</code>
    146      * -- ConvergenceException </li>
    147      * <li> <code> maximumIterations</code> iterations elapse
    148      * -- ConvergenceException </li></ul>
    149      *
    150      * @param function the function
    151      * @param initial initial midpoint of interval being expanded to
    152      * bracket a root
    153      * @param lowerBound lower bound (a is never lower than this value)
    154      * @param upperBound upper bound (b never is greater than this
    155      * value)
    156      * @param maximumIterations maximum number of iterations to perform
    157      * @return a two element array holding {a, b}.
    158      * @throws ConvergenceException if the algorithm fails to find a and b
    159      * satisfying the desired conditions
    160      * @throws FunctionEvaluationException if an error occurs evaluating the
    161      * function
    162      * @throws IllegalArgumentException if function is null, maximumIterations
    163      * is not positive, or initial is not between lowerBound and upperBound
    164      */
    165     public static double[] bracket(UnivariateRealFunction function,
    166             double initial, double lowerBound, double upperBound,
    167             int maximumIterations) throws ConvergenceException,
    168             FunctionEvaluationException {
    169
    170         if (function == null) {
    171             throw new IllegalArgumentException ("function is null.");
    172         }
    173         if (maximumIterations <= 0)  {
    174             throw new IllegalArgumentException
    175             ("bad value for maximumIterations: " + maximumIterations);
    176         }
    177         if (initial < lowerBound || initial > upperBound || lowerBound >= upperBound) {
    178             throw new IllegalArgumentException
    179             ("Invalid endpoint parameters:  lowerBound=" + lowerBound +
    180               " initial=" + initial + " upperBound=" + upperBound);
    181         }
    182         double a = initial;
    183         double b = initial;
    184         double fa;
    185         double fb;
    186         int numIterations = 0 ;
    187
    188         do {
    189             a = Math.max(a - 1.0, lowerBound);
    190             b = Math.min(b + 1.0, upperBound);
    191             fa = function.value(a);
    192
    193             fb = function.value(b);
    194             numIterations++ ;
    195         } while ((fa * fb > 0.0) && (numIterations < maximumIterations) &&
    196                 ((a > lowerBound) || (b < upperBound)));
    197
    198         if (fa * fb >= 0.0 ) {
    199             throw new ConvergenceException
    200             ("Number of iterations= " + numIterations +
    201               " maximum iterations= "  + maximumIterations +
    202               " initial= " + initial + " lowerBound=" + lowerBound +
    203               " upperBound=" + upperBound + " final a value=" + a +
    204               " final b value=" + b + " f(a)=" + fa + " f(b)=" + fb);
    205         }
    206
    207         return new double[]{a, b};
    208     }
    209
    210     /***
    211      * Compute the midpoint of two values.
    212      *
    213      * @param a first value.
    214      * @param b second value.
    215      * @return the midpoint.
    216      */
    217     public static double midpoint(double a, double b) {
    218         return (a + b) * .5;
    219     }
    220
    221     /***
    222      * Checks to see if f is null, throwing IllegalArgumentException if so.
    223      * Also initializes factory if factory is null.
    224      *
    225      * @param f  input function
    226      * @throws IllegalArgumentException if f is null
    227      */
    228     private static void setup(UnivariateRealFunction f) {
    229
    230         if (f == null) {
    231             throw new IllegalArgumentException("function can not be null.");
    232         }
    233
    234         if (factory == null) {
    235             factory = UnivariateRealSolverFactory.newInstance();
    236         }
    237     }
    238 }