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 }