- Direct Known Subclasses:
RootFinder.Brent
,RootFinder.Halley
,RootFinder.Newton
Solves the equation y = f(x,p) for x, where x and y have a type of double and p has a type set by a type parameter and represents a set of parameters that do not vary. There are Three implementations: RootFinder.Brent, RootFinder.Newton and RootFinder.Halley that implement Brent's, Newton's and Halley's method respectively. If the type parameter is not used, one should not be provided. Otherwise it is a user-defined type used in computing the values of f and its derivatives with respect to x.
Once a RootFinder rf is constructed, calling rf.findRoot(guess) finds a root by solving 0 = f(x, p). To solve the equation y = f(x, p), call solve(y, guess) for the Newton and Halley case, or solve(y, lower, upper) for the Brent case, and the value of x will be returned. The argument "guess" is an initial guess and the arguments "lower" and "upper" are values of x bracketing the value giving the solution, where f(lower, p) - y and f(upper, p) - y have opposite signs.
Parameters allow a single root finder to use a number of related functions without the need to allocate a new root finder for each case. If parameters are used, one will define a class to hold the parameters. For example,
class Parameters {
double p1;
double p2;
public Parameters(double x1, double y1) {
p1 = x1;
p2 = x2;
}
}
Then create the RootFinder:
Parameters parameters = new Parameters(1.0, 2.0);
RootFinder<Parameters> rf = new RootFinder.Newton<Parameters>(parameters)
{
public double function(double x) {
Parameters p = getParameters();
double p1 = p.x1;
double p2 = p.x2;
double result;
...
return result;
}
public double firstDerivative(double x) {
Parameters p = getParameters();
double p1 = p.x1;
double p2 = p.x2;
double result;
...
return result;
}
};
To change the parameters, either create a new instance parameters2 and call
rf.setParameters(parameters2) or call rf.getParameters() and manipulate the
values (e.g., by setting parameters.p1 and parameters.p2).
The criteria for convergence can be changed by calling setEpsilon. Limits on the number of iteration steps and some control over those can be set by calling setLimits.
-
Nested Class Summary
Nested ClassesModifier and TypeClassDescriptionstatic class
RootFinder class using Brent's method.static class
RootFinder exception class.static class
RootFinder class using Halley's method.static class
RootFinder class using Newton's method. -
Field Summary
Fields -
Constructor Summary
ConstructorsModifierConstructorDescriptionprotected
Constructor.protected
RootFinder
(P parameters) Constructor given parameters. -
Method Summary
Modifier and TypeMethodDescriptionabstract double
ferror
(double x) The error in the value returned byfunction
(x).final double
findRoot
(double... initialArgs) Find a root.abstract double
function
(double x) The function f in the equation y = f(x, p).double
Get the current value of epsilon The value of epsilon determines when the algorithm terminates normally.Get a RootFinder's parameters.static double
refineSolution
(DoubleUnaryOperator f, DoubleUnaryOperator df, double y, double guess) Use Newton's method to refine a solution to the equation f(x) - y = 0.void
setEpsilon
(double epsilon) Set the current value of epsilon.void
setEpsilon
(double epsilon, boolean relative) Set the current value of epsilon, indicating if it is a relative or absolute value.void
setLimit
(int maxLimit) Set the iterations limit.void
setParameters
(P parameters) Set a RootFinder's parameters.abstract double
solve
(double y, double... initialArgs) Solve an equation.static int
solveBezier
(double[] eqn, int offset, int n, double[] res) Find the zeros of a sum of Bernstein polynomial of degree n.static int
solveCubic
(double[] eqn) Solves the cubic equation ax3 + bx2 + cx + d and places the real roots into its argument array, which is used for both input and output.static int
solveCubic
(double[] eqn, double[] res) Solves the cubic equation ax3 + bx2 + cx + d = 0.static int
solvePolynomial
(double[] eqn, int n, double[] res) Find the zeros of a polynomial of degree n.static int
solvePolynomial
(BezierPolynomial p, double[] res) Find the zeros of aBezierPolynomial
.static int
solvePolynomial
(Polynomial p, double[] res) Find the zeros of aPolynomial
.static int
solveQuadratic
(double[] eqn) Solves the quadratic equation ax2 + bx + c and places the real roots into its argument array, which is used for both input and output.static int
solveQuadratic
(double[] eqn, double[] res) Solves the quadratic equation ax2 + bx + c.static int
solveQuartic
(double[] eqn) Solve the quartic equation ax4 + bx3 + cx2 + dx + e = 0 and store the results in the argument array.
-
Field Details
-
maxLimit
protected int maxLimitThe iteration limit.
-
-
Constructor Details
-
RootFinder
protected RootFinder()Constructor. -
RootFinder
Constructor given parameters.- Parameters:
parameters
- the parameters
-
-
Method Details
-
setParameters
Set a RootFinder's parameters. Parameters are used to provide additional values for computing. root-finder's function and their derivatives.- Parameters:
parameters
- the parameters
-
getParameters
Get a RootFinder's parameters. Parameters are used to provide additional values for computing. root-finder's function and their derivatives.- Returns:
- an instance of the class representing a root finder's parameters (this will be the same instance passed to setParameters)
-
findRoot
Find a root. Starting from an initial guess, findRoot returns the value of x that satisfies f(x, p) = 0, where p represents the current parameters. The implementation simply calls solve with the parameter y set to 0, and with the initial arguments that solve requires (which depend on the root-finding algorithm being used).- Parameters:
initialArgs
- the initial arguments required by a subclass- Returns:
- the root
- Throws:
RootFinder.ConvergenceException
- the method failed to convergeMathException
- an error occurred calling the function f or one of its derivatives.
-
solve
public abstract double solve(double y, double... initialArgs) throws RootFinder.ConvergenceException Solve an equation. Starting from an initial guess, findRoot returns the value of x that satisfies f(x, p) = y, where p represents the current parameters.- Parameters:
y
- the desired value of f(x, p)initialArgs
- the initial arguments required by a subclass- Returns:
- the value of x that satisfies f(x, p) = y
- Throws:
RootFinder.ConvergenceException
- the method failedMathException
- an error occurred calling the function f or one of its derivatives.
-
setLimit
public void setLimit(int maxLimit) Set the iterations limit. The limit is a bound on the maximum number of iterations.- Parameters:
maxLimit
- the maximum-iteration limit
-
getEpsilon
public double getEpsilon()Get the current value of epsilon The value of epsilon determines when the algorithm terminates normally. If the current value of x is in the range [-1.0, 1.0] the requirement for convergence is that the change in the value of x from the previous iteration is less than epsilon and that the corresponding value of y is within epsilon of its desired value. Otherwise the requirement for convergence is that the change in the value of x, with that change divided by x, is less than epsilon and that the corresponding value of y is within epsilon of its desired value.- Returns:
- the value of epsilon
-
setEpsilon
Set the current value of epsilon. The value of epsilon determines when the algorithm terminates normally. The error associated with evaluating the function will be set to epsilon.- Parameters:
epsilon
- the value of epsilon- Throws:
IllegalArgumentException
- the argument was not a positive real number
-
setEpsilon
Set the current value of epsilon, indicating if it is a relative or absolute value. When absolute, the error in evaluating the function will be epsilon when the default implementation offerror(double)
is used. When relative, the error will be the value of the function after that value is multiplied by epsilon.- Parameters:
epsilon
- the value of epsilonrelative
- true if the value is relative; false if it is absolute.- Throws:
IllegalArgumentException
- the first argument was not a positive real number
-
function
public abstract double function(double x) The function f in the equation y = f(x, p). The function f(x, p) is used with x varied and with the parameters p constant.- Parameters:
x
- the argument of the function- Returns:
- the value of the function f
- See Also:
-
ferror
public abstract double ferror(double x) The error in the value returned byfunction
(x). The default implementation returns the same value asgetEpsilon()
. Callers are encouraged to override this method.- Parameters:
x
- the argument passed tofunction(double)
- Returns:
- the error for the value
function
(x)
-
solveQuadratic
public static int solveQuadratic(double[] eqn) Solves the quadratic equation ax2 + bx + c and places the real roots into its argument array, which is used for both input and output.This method duplicates the behavior of the method
QuadCurve2D.solveQuadratic(double[])
. It is provided to break a dependency on the jdk.desktop module.- Parameters:
eqn
- the input/output array containing the roots (starting at index 0) for output and the coefficients ordered so that eqn[0] = c, eqn[1] = b, and eqn[2] = a.- Returns:
- the number of roots; -1 if the coefficients a and b are zero
-
solveQuadratic
public static int solveQuadratic(double[] eqn, double[] res) Solves the quadratic equation ax2 + bx + c. The length of the array res should be large enough to hold all of the roots (a length of at least 2 should be used).If multiple roots have the same value, only one is shown. If roots are sufficiently close that they cannot be distinguished due to floating-point limitations, only one root may be shown. This method duplicates the calling convention used by the method
QuadCurve2D.solveQuadratic(double[],double[])
. It was provided to break a dependency on the jdk.desktop module.- Parameters:
eqn
- the input/output array containing the roots (starting at index 0) for output and the coefficients ordered so that eqn[0] = c, eqn[1] = b, and eqn[2] = a.res
- the output array containing the roots (starting at index 0)- Returns:
- the number of roots; -1 if the coefficients a and b are zero
-
solveCubic
public static int solveCubic(double[] eqn) Solves the cubic equation ax3 + bx2 + cx + d and places the real roots into its argument array, which is used for both input and output.If multiple roots have the same value, only one of those roots is shown. If roots are sufficiently close that they cannot be distinguished due to floating-point limitations, only one root may be shown even though multiple distinct roots actually exist. This method duplicates the calling convention used by the method
CubicCurve2D.solveCubic(double[])
. It was provided to break a dependency on the jdk.desktop module, but experience with it indicates that is is more accurate than the implementation ofCubicCurve2D.solveCubic(double[])
provided in Java 11.0.1.A cubic function f(x) = ax3 + bx2 + cx + d has an inflection point at -b/(3a). The number of critical points is dependent on the value Δ0 = b2 - 3ac:
- if Δ is negative, there are no critical points and the cubic function is strictly monotonic. The implementation uses Newton's method, with the inflection point as an initial guess.
- if Δ is zero, there is one critical point, which is also an inflection point. There is only one real root, which is a triple root if the inflection point is a root, or a single root if the inflection point is not a root. In the later case, the root's value is (-b - (27ad - b3)1/3 / (3a).
- if Δ is positive, there are two critical points
with the inflection point between these two. The
critical points are $(-b-\sqrt{\Delta_0})/(3a)$
and $(-b + sqrt{\Delta_0})/(3a)$.
If one of these
is a root r (both cannot be roots), it is a double root, and
the other root has the value -(2ar + b)/a. Otherwise,
let x1 be the lower of these two critical points
and let x2 be the higher of these two values.
- If f(x1) and f″(x1) have opposite signs, there is a real root below x1. One can find the root using Newton's method with an initial guess $x_g = x_1 - (\frac32)^m\sqrt{-2f(x_1)/f\prime(x_1)}$, , where m is the smallest non-negative integer such that f(xg) and f(x1 have opposite signs.
- If f(x1) and f(x1) have
opposite signs, there is a root between these two
values. One can find the root using Newton's method
with a lower and upper bound and with the inflection
point as an initial guess. When upper and lower
bounds are provided, the
RootFinder.Newton
class will temporarily switch to Brent's method when Newton's method's next value is out of range, and will tighten the upper and lower bounds as new values are tried. - If f(x2) and f″(x2) have opposite signs, there is a real root above x2. One can find the root using Newton's method with an initial guess $x_g = x_2 - (\frac32)^m\sqrt{-2f(x_2)/f\prime(x_2)}$, where m is the smallest non-negative integer such that f(xg) and f(x2 have opposite signs.
- Parameters:
eqn
- the input/output array containing the roots (starting at index 0) for output and the coefficients ordered so that eqn[0] = d, eqn[1] = c, eqn[2] = b, and eqn[3] = a.- Returns:
- the number of roots; -1 if the coefficients a, b, and c are zero
-
refineSolution
public static double refineSolution(DoubleUnaryOperator f, DoubleUnaryOperator df, double y, double guess) Use Newton's method to refine a solution to the equation f(x) - y = 0. Newton's method uses a guess that is improved by approximating f(x) by x0 + f'(x0)(x - x0) and solving the linear equation y = x0 + f'(x0)(x - x0), iterating until a convergence criteria is satisfied or until the implementation detects that it will not converge. This implementation assumes the convergence criteria has been satisfied and iterates as long as the solution can be improved.The algorithm starts by setting x1 to the value of a guess that is provided by one of this methods arguments.
A sequence of values {x1, ... xn, ...} is generated as follows: For xn, if f(xn) - y = 0, then xn+1 = x
. Otherwise if f(xn) - y != 0, Newton's method is used to determine a value w = xn = (f(xn)-y)/f'(x1)) unless f'(x1) = 0, in which case w = xn. If the sign if f(xn) - y and f(w) - y differ, a point z = x1 + (u - x1) (f(x1 - y) / (f(x1) - f(w)) is also considered. Note that z is an interpolation point between x1 and w. The value of |f(x) - y| is compared at x1, w, and (if calculated) z and the value of xn+1 is set to the one which provides the minimal value of |f(x)-y|. This sequence converges when xn+1 = xn, and that must occur because there are a finite number of double-precision numbers between f(x1) and 0.
The final step in the algorithm is to compute f(x) at adjacent floating-point numbers and pick the one that provides the minimum value for |f(x) - y|. The initial guess will be used in this final step if the argument df is null.
- Parameters:
f
- the functiondf
- the function's derivative; null if only the final step is usedy
- the value for f(x)guess
- an initial guess (assumed to be very close to the solution)- Returns:
- a value x such that f(x) = y to within floating-point errors
-
solveCubic
public static int solveCubic(double[] eqn, double[] res) Solves the cubic equation ax3 + bx2 + cx + d = 0. The length of the array res should be large enough to hold all of the real roots (a length of at least 3 should be used unless one knows a priori that there are less than 3 real roots).If multiple real roots have the same value, only one of those roots is shown. If roots are sufficiently close that they cannot be distinguished due to floating-point limitations, only one root may be shown even though multiple distinct roots actually exist. This method duplicates the calling convention used by the method
CubicCurve2D.solveCubic(double[],double[])
. It was provided to break a dependency on the jdk.desktop module, but experience with it indicates that is is more accurate than the implementation ofCubicCurve2D.solveCubic(double[],double[])
provided in Java 11.0.1.A cubic function f(x) = ax3 + bx2 + cx + d has an inflection point at -b/(3a). The number of critical points is dependent on the value Δ0 = b2 - 3ac:
- if Δ is negative, there are no critical points and the cubic function is strictly monotonic. The implementation uses Newton's method, with the inflection point as an initial guess.
- if Δ is zero, there is one critical point, which is also an inflection point. There is only one real root, which is a triple root if the inflection point is a root, or a single root if the inflection point is not a root. In the later case, the root's value is (-b - (27ad - b3)1/3 / (3a).
- if Δ is positive, there are two critical points
with the inflection point between these two. The
critical points are
$(-b - \sqrt{\Delta_0})/(3a)$
and
$(-b + \sqrt{\Delta_0})/(3a)$.
If one of these is a root r (both cannot be roots), it
is a double root, and the other root has the value
-(2ar + b)/a. Otherwise, let x1 be the lower of
these two critical points and let x2 be the
higher of these two values.
- If f(x1) and f″(x1) have opposite signs, there is a real root below x1. One can find the root using Newton's method with an initial guess $x_g = x_1 - (\frac32)^m\sqrt{-2f(x_1)/f\prime(x_1)}$, where m is the smallest non-negative integer such that f(xg) and f(x1 have opposite signs.
- If f(x1) and f(x1) have
opposite signs, there is a root between these two
values. One can find the root using Newton's method
with a lower and upper bound and with the inflection
point as an initial guess. When upper and lower
bounds are provided, the
RootFinder.Newton
class will temporarily switch to Brent's method when Newton's method's next value is out of range, and will tighten the upper and lower bounds as new values are tried. - If f(x2) and f″(x2) have opposite signs, there is a real root above x2. One can find the root using Newton's method with an initial guess $x_g = x_2 - (\frac32)^m\sqrt{-2f(x_2)/f\prime(x_2)}$, xg = x2 - (3/2)msqrt(-2f″(x2)/f(x2), where m is the smallest non-negative integer such that f(xg) and f(x2 have opposite signs.
- Parameters:
eqn
- the input/output array containing the roots (starting at index 0) for output and the coefficients ordered so that eqn[0] = d, eqn[1] = c, eqn[2] = b, and eqn[3] = a.res
- the output array containing the roots (starting at index 0) in numerical order- Returns:
- the number of roots; -1 if the coefficients a, b, and c are zero
-
solvePolynomial
Find the zeros of aPolynomial
.- Parameters:
p
- the polynomialres
- an array to hold the results- Returns:
- the number of roots; -1 if the degree of the polynomial is zero.
- Throws:
IllegalArgumentException
- if the array argument is too short
-
solvePolynomial
public static int solvePolynomial(double[] eqn, int n, double[] res) Find the zeros of a polynomial of degree n.- Parameters:
eqn
- the coefficients of the polynomial for input, ordered so that eqn[i] provides the coefficient for xi where x is the polynomial's indeterminaten
- the degree of the polynomialres
- the roots that were found- Returns:
- the number of roots; -1 if the degree of the polynomial is zero.
- Throws:
IllegalArgumentException
- n is negative or the array arguments are too short
-
solveQuartic
public static int solveQuartic(double[] eqn) Solve the quartic equation ax4 + bx3 + cx2 + dx + e = 0 and store the results in the argument array. The algorithm is based on the description in the Wikipedia article on the quartic function.A number of special cases are handled explicitly and for coefficients that are integer values, a computation of the discriminant is performed using two's complement arithmetic with long integers. This is equivalent to modular arithmetic (modulo 264) and allows one to reliably determine if the discriminant is zero, and to improve the accuracy of determining if the discriminant is positive, negative, or zero. Also, the real-valued form of the discriminant uses Kahan's additional algorithm to improve the accuracy of the sum.The same procedure is used for some other quantities as well. These are combined to determine the number of roots based on whether these quantities are positive, zero, or negative. Cases that are not covered should not occur unless as the result of floating point errors, in which case a numerical method is used instead. To further reduce floating-point errors, Newton's method and some floating-point adjustments are made to improve the accuracy of the solutions. In addition, for cases where there are known to be multiple roots, the critical points are computed by solving a cubic equation and those solutions are used when they are also solutions. If floating-point errors split a single multiple root of the polynomial into multiple zeros, the critical points are used to eliminate these spurious solutions.
- Parameters:
eqn
- an array containing the coefficients of the polynomial, where eqn[0] = e, eqn[1] = d, eqn[2] = c, eqn[3] = d, and eqn[4] = e (the array length must be at least 5), also used as the output array- Returns:
- the number of real-valued solutions; -1 if a, b, c, and d are 0
-
solvePolynomial
Find the zeros of aBezierPolynomial
.- Parameters:
p
- the polynomialres
- an array to hold the results- Returns:
- the number of roots; -1 if the degree of the polynomial is zero.
- Throws:
IllegalArgumentException
- if the array argument is too short
-
solveBezier
public static int solveBezier(double[] eqn, int offset, int n, double[] res) throws IllegalArgumentException Find the zeros of a sum of Bernstein polynomial of degree n. The name of the method, solveBezier, reflects the use of Bernstein polynomials for Bézier curves. The roots returned are in the range [0.0, 1.0], the range appropriate for Bézier curves.- Parameters:
eqn
- the coefficients of the polynomial for input, where eqn[i] is the coefficient for the Bernstein polynomial Bi,n(t)offset
- the offset into the array eqn at which the coefficients startn
- the degree of the polynomialres
- the roots that were found- Returns:
- the number of roots; -1 if the degree of the polynomial is zero.
- Throws:
IllegalArgumentException
- n is negative or the array arguments are too short
-