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
018package org.apache.commons.math3.analysis.solvers;
019
020import org.apache.commons.math3.util.FastMath;
021import org.apache.commons.math3.exception.NoBracketingException;
022import org.apache.commons.math3.exception.TooManyEvaluationsException;
023
024/**
025 * Implements the <em>Secant</em> method for root-finding (approximating a
026 * zero of a univariate real function). The solution that is maintained is
027 * not bracketed, and as such convergence is not guaranteed.
028 *
029 * <p>Implementation based on the following article: M. Dowell and P. Jarratt,
030 * <em>A modified regula falsi method for computing the root of an
031 * equation</em>, BIT Numerical Mathematics, volume 11, number 2,
032 * pages 168-174, Springer, 1971.</p>
033 *
034 * <p>Note that since release 3.0 this class implements the actual
035 * <em>Secant</em> algorithm, and not a modified one. As such, the 3.0 version
036 * is not backwards compatible with previous versions. To use an algorithm
037 * similar to the pre-3.0 releases, use the
038 * {@link IllinoisSolver <em>Illinois</em>} algorithm or the
039 * {@link PegasusSolver <em>Pegasus</em>} algorithm.</p>
040 *
041 * @version $Id: SecantSolver.java 1379560 2012-08-31 19:40:30Z erans $
042 */
043public class SecantSolver extends AbstractUnivariateSolver {
044
045    /** Default absolute accuracy. */
046    protected static final double DEFAULT_ABSOLUTE_ACCURACY = 1e-6;
047
048    /** Construct a solver with default accuracy (1e-6). */
049    public SecantSolver() {
050        super(DEFAULT_ABSOLUTE_ACCURACY);
051    }
052
053    /**
054     * Construct a solver.
055     *
056     * @param absoluteAccuracy absolute accuracy
057     */
058    public SecantSolver(final double absoluteAccuracy) {
059        super(absoluteAccuracy);
060    }
061
062    /**
063     * Construct a solver.
064     *
065     * @param relativeAccuracy relative accuracy
066     * @param absoluteAccuracy absolute accuracy
067     */
068    public SecantSolver(final double relativeAccuracy,
069                        final double absoluteAccuracy) {
070        super(relativeAccuracy, absoluteAccuracy);
071    }
072
073    /** {@inheritDoc} */
074    @Override
075    protected final double doSolve()
076        throws TooManyEvaluationsException,
077               NoBracketingException {
078        // Get initial solution
079        double x0 = getMin();
080        double x1 = getMax();
081        double f0 = computeObjectiveValue(x0);
082        double f1 = computeObjectiveValue(x1);
083
084        // If one of the bounds is the exact root, return it. Since these are
085        // not under-approximations or over-approximations, we can return them
086        // regardless of the allowed solutions.
087        if (f0 == 0.0) {
088            return x0;
089        }
090        if (f1 == 0.0) {
091            return x1;
092        }
093
094        // Verify bracketing of initial solution.
095        verifyBracketing(x0, x1);
096
097        // Get accuracies.
098        final double ftol = getFunctionValueAccuracy();
099        final double atol = getAbsoluteAccuracy();
100        final double rtol = getRelativeAccuracy();
101
102        // Keep finding better approximations.
103        while (true) {
104            // Calculate the next approximation.
105            final double x = x1 - ((f1 * (x1 - x0)) / (f1 - f0));
106            final double fx = computeObjectiveValue(x);
107
108            // If the new approximation is the exact root, return it. Since
109            // this is not an under-approximation or an over-approximation,
110            // we can return it regardless of the allowed solutions.
111            if (fx == 0.0) {
112                return x;
113            }
114
115            // Update the bounds with the new approximation.
116            x0 = x1;
117            f0 = f1;
118            x1 = x;
119            f1 = fx;
120
121            // If the function value of the last approximation is too small,
122            // given the function value accuracy, then we can't get closer to
123            // the root than we already are.
124            if (FastMath.abs(f1) <= ftol) {
125                return x1;
126            }
127
128            // If the current interval is within the given accuracies, we
129            // are satisfied with the current approximation.
130            if (FastMath.abs(x1 - x0) < FastMath.max(rtol * FastMath.abs(x1), atol)) {
131                return x1;
132            }
133        }
134    }
135
136}