View Javadoc

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