001 /*
002 * Licensed to the Apache Software Foundation (ASF) under one or more
003 * contributor license agreements. See the NOTICE file distributed with
004 * this work for additional information regarding copyright ownership.
005 * The ASF licenses this file to You under the Apache License, Version 2.0
006 * (the "License"); you may not use this file except in compliance with
007 * the License. You may obtain a copy of the License at
008 *
009 * http://www.apache.org/licenses/LICENSE-2.0
010 *
011 * Unless required by applicable law or agreed to in writing, software
012 * distributed under the License is distributed on an "AS IS" BASIS,
013 * WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied.
014 * See the License for the specific language governing permissions and
015 * limitations under the License.
016 */
017 package org.apache.commons.math.analysis.solvers;
018
019 import org.apache.commons.math.ConvergenceException;
020 import org.apache.commons.math.FunctionEvaluationException;
021 import org.apache.commons.math.MathRuntimeException;
022 import org.apache.commons.math.analysis.UnivariateRealFunction;
023 import org.apache.commons.math.exception.util.LocalizedFormats;
024 import org.apache.commons.math.exception.NullArgumentException;
025 import org.apache.commons.math.util.FastMath;
026
027 /**
028 * Utility routines for {@link UnivariateRealSolver} objects.
029 *
030 * @version $Revision: 1070725 $ $Date: 2011-02-15 02:31:12 +0100 (mar. 15 f??vr. 2011) $
031 */
032 public class UnivariateRealSolverUtils {
033
034 /**
035 * Default constructor.
036 */
037 private UnivariateRealSolverUtils() {
038 super();
039 }
040
041 /**
042 * Convenience method to find a zero of a univariate real function. A default
043 * solver is used.
044 *
045 * @param f the function.
046 * @param x0 the lower bound for the interval.
047 * @param x1 the upper bound for the interval.
048 * @return a value where the function is zero.
049 * @throws ConvergenceException if the iteration count was exceeded
050 * @throws FunctionEvaluationException if an error occurs evaluating the function
051 * @throws IllegalArgumentException if f is null or the endpoints do not
052 * specify a valid interval
053 */
054 public static double solve(UnivariateRealFunction f, double x0, double x1)
055 throws ConvergenceException, FunctionEvaluationException {
056 setup(f);
057 return LazyHolder.FACTORY.newDefaultSolver().solve(f, x0, x1);
058 }
059
060 /**
061 * Convenience method to find a zero of a univariate real function. A default
062 * solver is used.
063 *
064 * @param f the function
065 * @param x0 the lower bound for the interval
066 * @param x1 the upper bound for the interval
067 * @param absoluteAccuracy the accuracy to be used by the solver
068 * @return a value where the function is zero
069 * @throws ConvergenceException if the iteration count is exceeded
070 * @throws FunctionEvaluationException if an error occurs evaluating the function
071 * @throws IllegalArgumentException if f is null, the endpoints do not
072 * specify a valid interval, or the absoluteAccuracy is not valid for the
073 * default solver
074 */
075 public static double solve(UnivariateRealFunction f, double x0, double x1,
076 double absoluteAccuracy) throws ConvergenceException,
077 FunctionEvaluationException {
078
079 setup(f);
080 UnivariateRealSolver solver = LazyHolder.FACTORY.newDefaultSolver();
081 solver.setAbsoluteAccuracy(absoluteAccuracy);
082 return solver.solve(f, x0, x1);
083 }
084
085 /**
086 * This method attempts to find two values a and b satisfying <ul>
087 * <li> <code> lowerBound <= a < initial < b <= upperBound</code> </li>
088 * <li> <code> f(a) * f(b) < 0 </code></li>
089 * </ul>
090 * If f is continuous on <code>[a,b],</code> this means that <code>a</code>
091 * and <code>b</code> bracket a root of f.
092 * <p>
093 * The algorithm starts by setting
094 * <code>a := initial -1; b := initial +1,</code> examines the value of the
095 * function at <code>a</code> and <code>b</code> and keeps moving
096 * the endpoints out by one unit each time through a loop that terminates
097 * when one of the following happens: <ul>
098 * <li> <code> f(a) * f(b) < 0 </code> -- success!</li>
099 * <li> <code> a = lower </code> and <code> b = upper</code>
100 * -- ConvergenceException </li>
101 * <li> <code> Integer.MAX_VALUE</code> iterations elapse
102 * -- ConvergenceException </li>
103 * </ul></p>
104 * <p>
105 * <strong>Note: </strong> this method can take
106 * <code>Integer.MAX_VALUE</code> iterations to throw a
107 * <code>ConvergenceException.</code> Unless you are confident that there
108 * is a root between <code>lowerBound</code> and <code>upperBound</code>
109 * near <code>initial,</code> it is better to use
110 * {@link #bracket(UnivariateRealFunction, double, double, double, int)},
111 * explicitly specifying the maximum number of iterations.</p>
112 *
113 * @param function the function
114 * @param initial initial midpoint of interval being expanded to
115 * bracket a root
116 * @param lowerBound lower bound (a is never lower than this value)
117 * @param upperBound upper bound (b never is greater than this
118 * value)
119 * @return a two element array holding {a, b}
120 * @throws ConvergenceException if a root can not be bracketted
121 * @throws FunctionEvaluationException if an error occurs evaluating the 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 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 NullArgumentException(LocalizedFormats.FUNCTION);
172 }
173 if (maximumIterations <= 0) {
174 throw MathRuntimeException.createIllegalArgumentException(
175 LocalizedFormats.INVALID_MAX_ITERATIONS, maximumIterations);
176 }
177 if (initial < lowerBound || initial > upperBound || lowerBound >= upperBound) {
178 throw MathRuntimeException.createIllegalArgumentException(
179 LocalizedFormats.INVALID_BRACKETING_PARAMETERS,
180 lowerBound, initial, 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 = FastMath.max(a - 1.0, lowerBound);
190 b = FastMath.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 LocalizedFormats.FAILED_BRACKETING,
201 numIterations, maximumIterations, initial,
202 lowerBound, upperBound, a, b, fa, fb);
203 }
204
205 return new double[]{a, b};
206 }
207
208 /**
209 * Compute the midpoint of two values.
210 *
211 * @param a first value.
212 * @param b second value.
213 * @return the midpoint.
214 */
215 public static double midpoint(double a, double b) {
216 return (a + b) * .5;
217 }
218
219 /**
220 * Checks to see if f is null, throwing IllegalArgumentException if so.
221 * @param f input function
222 * @throws IllegalArgumentException if f is null
223 */
224 private static void setup(UnivariateRealFunction f) {
225 if (f == null) {
226 throw new NullArgumentException(LocalizedFormats.FUNCTION);
227 }
228 }
229
230 // CHECKSTYLE: stop HideUtilityClassConstructor
231 /** Holder for the factory.
232 * <p>We use here the Initialization On Demand Holder Idiom.</p>
233 */
234 private static class LazyHolder {
235 /** Cached solver factory */
236 private static final UnivariateRealSolverFactory FACTORY = UnivariateRealSolverFactory.newInstance();
237 }
238 // CHECKSTYLE: resume HideUtilityClassConstructor
239
240 }