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
18 package org.apache.commons.math4.legacy.analysis.solvers;
19
20 import org.apache.commons.math4.legacy.exception.NoBracketingException;
21 import org.apache.commons.math4.legacy.exception.TooManyEvaluationsException;
22 import org.apache.commons.math4.core.jdkmath.JdkMath;
23
24 /**
25 * Implements the <em>Secant</em> method for root-finding (approximating a
26 * zero of a univariate real function). The solution that is maintained is
27 * not bracketed, and as such convergence is not guaranteed.
28 *
29 * <p>Implementation based on the following article: M. Dowell and P. Jarratt,
30 * <em>A modified regula falsi method for computing the root of an
31 * equation</em>, BIT Numerical Mathematics, volume 11, number 2,
32 * pages 168-174, Springer, 1971.</p>
33 *
34 * <p>Note that since release 3.0 this class implements the actual
35 * <em>Secant</em> algorithm, and not a modified one. As such, the 3.0 version
36 * is not backwards compatible with previous versions. To use an algorithm
37 * similar to the pre-3.0 releases, use the
38 * {@link IllinoisSolver <em>Illinois</em>} algorithm or the
39 * {@link PegasusSolver <em>Pegasus</em>} algorithm.</p>
40 *
41 */
42 public class SecantSolver extends AbstractUnivariateSolver {
43
44 /** Default absolute accuracy. */
45 protected static final double DEFAULT_ABSOLUTE_ACCURACY = 1e-6;
46
47 /** Construct a solver with default accuracy (1e-6). */
48 public SecantSolver() {
49 super(DEFAULT_ABSOLUTE_ACCURACY);
50 }
51
52 /**
53 * Construct a solver.
54 *
55 * @param absoluteAccuracy absolute accuracy
56 */
57 public SecantSolver(final double absoluteAccuracy) {
58 super(absoluteAccuracy);
59 }
60
61 /**
62 * Construct a solver.
63 *
64 * @param relativeAccuracy relative accuracy
65 * @param absoluteAccuracy absolute accuracy
66 */
67 public SecantSolver(final double relativeAccuracy,
68 final double absoluteAccuracy) {
69 super(relativeAccuracy, absoluteAccuracy);
70 }
71
72 /** {@inheritDoc} */
73 @Override
74 protected final double doSolve()
75 throws TooManyEvaluationsException,
76 NoBracketingException {
77 // Get initial solution
78 double x0 = getMin();
79 double x1 = getMax();
80 double f0 = computeObjectiveValue(x0);
81 double f1 = computeObjectiveValue(x1);
82
83 // If one of the bounds is the exact root, return it. Since these are
84 // not under-approximations or over-approximations, we can return them
85 // regardless of the allowed solutions.
86 if (f0 == 0.0) {
87 return x0;
88 }
89 if (f1 == 0.0) {
90 return x1;
91 }
92
93 // Verify bracketing of initial solution.
94 verifyBracketing(x0, x1);
95
96 // Get accuracies.
97 final double ftol = getFunctionValueAccuracy();
98 final double atol = getAbsoluteAccuracy();
99 final double rtol = getRelativeAccuracy();
100
101 // Keep finding better approximations.
102 while (true) {
103 // Calculate the next approximation.
104 final double x = x1 - ((f1 * (x1 - x0)) / (f1 - f0));
105 final double fx = computeObjectiveValue(x);
106
107 // If the new approximation is the exact root, return it. Since
108 // this is not an under-approximation or an over-approximation,
109 // we can return it regardless of the allowed solutions.
110 if (fx == 0.0) {
111 return x;
112 }
113
114 // Update the bounds with the new approximation.
115 x0 = x1;
116 f0 = f1;
117 x1 = x;
118 f1 = fx;
119
120 // If the function value of the last approximation is too small,
121 // given the function value accuracy, then we can't get closer to
122 // the root than we already are.
123 if (JdkMath.abs(f1) <= ftol) {
124 return x1;
125 }
126
127 // If the current interval is within the given accuracies, we
128 // are satisfied with the current approximation.
129 if (JdkMath.abs(x1 - x0) < JdkMath.max(rtol * JdkMath.abs(x1), atol)) {
130 return x1;
131 }
132 }
133 }
134 }