mdiggory 2003/06/03 19:17:17 Modified: math/src/java/org/apache/commons/math RootFinding.java Log: PR: http://nagoya.apache.org/bugzilla/show_bug.cgi?id=20375 Submitted by: Albert Davidson Chou Revision Changes Path 1.2 +61 -33 jakarta-commons-sandbox/math/src/java/org/apache/commons/math/RootFinding.java Index: RootFinding.java =================================================================== RCS file: /home/cvs/jakarta-commons-sandbox/math/src/java/org/apache/commons/math/RootFinding.java,v retrieving revision 1.1 retrieving revision 1.2 diff -u -r1.1 -r1.2 --- RootFinding.java 29 May 2003 19:34:37 -0000 1.1 +++ RootFinding.java 4 Jun 2003 02:17:17 -0000 1.2 @@ -55,83 +55,111 @@ /** * Utility class comprised of root finding techniques. - * + * * @author Brent Worden */ public class RootFinding { /** Maximum allowed numerical error. */ private static final double EPSILON = 10e-9; - /** - * Default constructor. Prohibit construction. - */ - private RootFinding(){ - super(); - } + /** + * Default constructor. Prohibit construction. + */ + private RootFinding(){ + super(); + } /** * For a function, f, this method returns two values, a and b that bracket * a root of f. That is to say, there exists a value c between a and b * such that f(c) = 0. - * + * * @param function the function * @param initial midpoint of the returned range. * @param lowerBound for numerical safety, a never is less than this value. * @param upperBound for numerical safety, b never is greater than this value. * @return a two element array holding {a, b}. */ - public static double[] bracket(UnivariateFunction function, - double initial, - double lowerBound, + public static double[] bracket(UnivariateFunction function, + double initial, + double lowerBound, double upperBound){ + return bracket( function, initial, lowerBound, upperBound, Integer.MAX_VALUE ) ; + } + + /** + * For a function, f, this method returns two values, a and b that bracket + * a root of f. That is to say, there exists a value c between a and b + * such that f(c) = 0. + * + * @param function the function + * @param initial midpoint of the returned range. + * @param lowerBound for numerical safety, a never is less than this value. + * @param upperBound for numerical safety, b never is greater than this value. + * @param maximumIterations to guard against infinite looping, maximum number of iterations to perform + * @return a two element array holding {a, b}. + */ + public static double[] bracket(UnivariateFunction function, + double initial, + double lowerBound, + double upperBound, + int maximumIterations){ double a = initial; double b = initial; double fa; double fb; - - do { + int numIterations = 0 ; + + do { a = Math.max(a - 1.0, lowerBound); b = Math.min(b + 1.0, upperBound); fa = function.evaluate(a); fb = function.evaluate(b); - } while(fa * fb > 0.0); - - return new double[]{a, b}; + numIterations += 1 ; + } while ( (fa * fb > 0.0) && ( numIterations < maximumIterations ) ); + + return new double[]{a, b}; } - + /** * For a function, f, this method returns a root c that lies between a and * b, and satisfies f(c) = 0. - * + * * @param function the function - * @param a lower bound of a root - * @param b upper bound of a root + * @param a lower (or upper) bound of a root + * @param b upper (or lower) bound of a root * @return a root of f */ - public static double bisection(UnivariateFunction function, - double a, + public static double bisection(UnivariateFunction function, + double a, double b){ double m; double fm; double fa; - double fb; - + + if ( b < a ) + { + double xtemp = a ; + a = b ; + b = xtemp ; + } + + fa = function.evaluate(a); + while(Math.abs(a - b) > EPSILON){ - fa = function.evaluate(a); - fb = function.evaluate(b); - - m = a + (b - a) / 2.0; // midpoint + m = (a + b) * 0.5; // midpoint fm = function.evaluate(m); if(fm * fa > 0.0){ + // b and m bracket the root. + a = m; fa = fm; - a = m; } else { - fb = fm; - b = m; + // a and m bracket the root. + b = m; } } - - return a + (b - a) / 2.0; + + return (a + b) * 0.5; } }
--------------------------------------------------------------------- To unsubscribe, e-mail: [EMAIL PROTECTED] For additional commands, e-mail: [EMAIL PROTECTED]
