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.solvers;
18
19 import org.apache.commons.math.FunctionEvaluationException;
20 import org.apache.commons.math.ConvergenceException;
21 import org.apache.commons.math.MathRuntimeException;
22 import org.apache.commons.math.analysis.UnivariateRealFunction;
23
24 /**
25 * Utility routines for {@link UnivariateRealSolver} objects.
26 *
27 * @version $Revision: 799857 $ $Date: 2009-08-01 09:07:12 -0400 (Sat, 01 Aug 2009) $
28 */
29 public class UnivariateRealSolverUtils {
30 /**
31 * Default constructor.
32 */
33 private UnivariateRealSolverUtils() {
34 super();
35 }
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 LazyHolder.FACTORY.newDefaultSolver().solve(f, 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 = LazyHolder.FACTORY.newDefaultSolver();
79 solver.setAbsoluteAccuracy(absoluteAccuracy);
80 return solver.solve(f, 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></p>
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.</p>
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></p>
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 MathRuntimeException.createIllegalArgumentException("function is null");
172 }
173 if (maximumIterations <= 0) {
174 throw MathRuntimeException.createIllegalArgumentException(
175 "bad value for maximum iterations number: {0}", maximumIterations);
176 }
177 if (initial < lowerBound || initial > upperBound || lowerBound >= upperBound) {
178 throw MathRuntimeException.createIllegalArgumentException(
179 "invalid bracketing parameters: lower bound={0}, initial={1}, upper bound={2}",
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 = 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={0}, maximum iterations={1}, " +
201 "initial={2}, lower bound={3}, upper bound={4}, final a value={5}, " +
202 "final b value={6}, f(a)={7}, f(b)={8}",
203 numIterations, maximumIterations, initial,
204 lowerBound, upperBound, a, b, fa, 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 * @param f input function
224 * @throws IllegalArgumentException if f is null
225 */
226 private static void setup(UnivariateRealFunction f) {
227 if (f == null) {
228 throw MathRuntimeException.createIllegalArgumentException("function is null");
229 }
230 }
231
232 /** Holder for the factory.
233 * <p>We use here the Initialization On Demand Holder Idiom.</p>
234 */
235 private static class LazyHolder {
236 /** Cached solver factory */
237 private static final UnivariateRealSolverFactory FACTORY =
238 UnivariateRealSolverFactory.newInstance();
239 }
240
241 }