| [section:roots2 Root Finding Without Derivatives] |
| |
| [h4 Synopsis] |
| |
| `` |
| #include <boost/math/tools/roots.hpp> |
| `` |
| |
| namespace boost{ namespace math{ namespace tools{ |
| |
| template <class F, class T, class Tol> |
| std::pair<T, T> |
| bisect( |
| F f, |
| T min, |
| T max, |
| Tol tol, |
| boost::uintmax_t& max_iter); |
| |
| template <class F, class T, class Tol> |
| std::pair<T, T> |
| bisect( |
| F f, |
| T min, |
| T max, |
| Tol tol); |
| |
| template <class F, class T, class Tol, class ``__Policy``> |
| std::pair<T, T> |
| bisect( |
| F f, |
| T min, |
| T max, |
| Tol tol, |
| boost::uintmax_t& max_iter, |
| const ``__Policy``&); |
| |
| template <class F, class T, class Tol> |
| std::pair<T, T> |
| bracket_and_solve_root( |
| F f, |
| const T& guess, |
| const T& factor, |
| bool rising, |
| Tol tol, |
| boost::uintmax_t& max_iter); |
| |
| template <class F, class T, class Tol, class ``__Policy``> |
| std::pair<T, T> |
| bracket_and_solve_root( |
| F f, |
| const T& guess, |
| const T& factor, |
| bool rising, |
| Tol tol, |
| boost::uintmax_t& max_iter, |
| const ``__Policy``&); |
| |
| template <class F, class T, class Tol> |
| std::pair<T, T> |
| toms748_solve( |
| F f, |
| const T& a, |
| const T& b, |
| Tol tol, |
| boost::uintmax_t& max_iter); |
| |
| template <class F, class T, class Tol, class ``__Policy``> |
| std::pair<T, T> |
| toms748_solve( |
| F f, |
| const T& a, |
| const T& b, |
| Tol tol, |
| boost::uintmax_t& max_iter, |
| const ``__Policy``&); |
| |
| template <class F, class T, class Tol> |
| std::pair<T, T> |
| toms748_solve( |
| F f, |
| const T& a, |
| const T& b, |
| const T& fa, |
| const T& fb, |
| Tol tol, |
| boost::uintmax_t& max_iter); |
| |
| template <class F, class T, class Tol, class ``__Policy``> |
| std::pair<T, T> |
| toms748_solve( |
| F f, |
| const T& a, |
| const T& b, |
| const T& fa, |
| const T& fb, |
| Tol tol, |
| boost::uintmax_t& max_iter, |
| const ``__Policy``&); |
| |
| // Termination conditions: |
| template <class T> |
| struct eps_tolerance; |
| |
| struct equal_floor; |
| struct equal_ceil; |
| struct equal_nearest_integer; |
| |
| }}} // namespaces |
| |
| [h4 Description] |
| |
| These functions solve the root of some function /f(x)/ without the |
| need for the derivatives of /f(x)/. The functions here that use TOMS |
| Algorithm 748 are asymptotically the most efficient known, and have |
| been shown to be optimal for a certain classes of smooth functions. |
| |
| Alternatively, there is a simple bisection routine which can be useful |
| in its own right in some situations, or alternatively for narrowing |
| down the range containing the root, prior to calling a more advanced |
| algorithm. |
| |
| All the algorithms in this section reduce the diameter of the enclosing |
| interval with the same asymptotic efficiency with which they locate the |
| root. This is in contrast to the derivative based methods which may /never/ |
| significantly reduce the enclosing interval, even though they rapidly approach |
| the root. This is also in contrast to some other derivative-free methods |
| (for example the methods of [@http://en.wikipedia.org/wiki/Brent%27s_method Brent or Dekker)] |
| which only reduce the enclosing interval on the final step. |
| Therefore these methods return a std::pair containing the enclosing interval found, |
| and accept a function object specifying the termination condition. |
| Three function objects are provided for ready-made termination conditions: |
| /eps_tolerance/ causes termination when the relative error in the enclosing |
| interval is below a certain threshold, while /equal_floor/ and /equal_ceil/ are |
| useful for certain statistical applications where the result is known to be |
| an integer. Other user-defined termination conditions are likely to be used |
| only rarely, but may be useful in some specific circumstances. |
| |
| template <class F, class T, class Tol> |
| std::pair<T, T> |
| bisect( |
| F f, |
| T min, |
| T max, |
| Tol tol, |
| boost::uintmax_t& max_iter); |
| |
| template <class F, class T, class Tol> |
| std::pair<T, T> |
| bisect( |
| F f, |
| T min, |
| T max, |
| Tol tol); |
| |
| template <class F, class T, class Tol, class ``__Policy``> |
| std::pair<T, T> |
| bisect( |
| F f, |
| T min, |
| T max, |
| Tol tol, |
| boost::uintmax_t& max_iter, |
| const ``__Policy``&); |
| |
| These functions locate the root using bisection, function arguments are: |
| |
| [variablelist |
| [[f] [A unary functor which is the function whose root is to be found.]] |
| [[min] [The left bracket of the interval known to contain the root.]] |
| [[max] [The right bracket of the interval known to contain the root. |
| It is a precondition that /min < max/ and /f(min)*f(max) <= 0/, |
| the function signals evaluation error if these preconditions are violated. |
| The action taken is controlled by the evaluation error policy. |
| A best guess may be returned, perhaps significantly wrong.]] |
| [[tol] [A binary functor that specifies the termination condition: the function |
| will return the current brackets enclosing the root when /tol(min,max)/ becomes true.]] |
| [[max_iter][The maximum number of invocations of /f(x)/ to make while searching for the root.]] |
| ] |
| |
| [optional_policy] |
| |
| Returns: a pair of values /r/ that bracket the root so that: |
| |
| f(r.first) * f(r.second) <= 0 |
| |
| and either |
| |
| tol(r.first, r.second) == true |
| |
| or |
| |
| max_iter >= m |
| |
| where /m/ is the initial value of /max_iter/ passed to the function. |
| |
| In other words, it's up to the caller to verify whether termination occurred |
| as a result of exceeding /max_iter/ function invocations (easily done by |
| checking the value of /max_iter/ when the function returns), rather than |
| because the termination condition /tol/ was satisfied. |
| |
| template <class F, class T, class Tol> |
| std::pair<T, T> |
| bracket_and_solve_root( |
| F f, |
| const T& guess, |
| const T& factor, |
| bool rising, |
| Tol tol, |
| boost::uintmax_t& max_iter); |
| |
| template <class F, class T, class Tol, class ``__Policy``> |
| std::pair<T, T> |
| bracket_and_solve_root( |
| F f, |
| const T& guess, |
| const T& factor, |
| bool rising, |
| Tol tol, |
| boost::uintmax_t& max_iter, |
| const ``__Policy``&); |
| |
| This is a convenience function that calls /toms748_solve/ internally |
| to find the root of /f(x)/. It's usable only when /f(x)/ is a monotonic |
| function, and the location of the root is known approximately, and in |
| particular it is known whether the root is occurs for positive or negative |
| /x/. The parameters are: |
| |
| [variablelist |
| [[f][A unary functor that is the function whose root is to be solved. |
| f(x) must be uniformly increasing or decreasing on /x/.]] |
| [[guess][An initial approximation to the root]] |
| [[factor][A scaling factor that is used to bracket the root: the value |
| /guess/ is multiplied (or divided as appropriate) by /factor/ |
| until two values are found that bracket the root. A value |
| such as 2 is a typical choice for /factor/.]] |
| [[rising][Set to /true/ if /f(x)/ is rising on /x/ and /false/ if /f(x)/ |
| is falling on /x/. This value is used along with the result |
| of /f(guess)/ to determine if /guess/ is |
| above or below the root.]] |
| [[tol] [A binary functor that determines the termination condition for the search |
| for the root. /tol/ is passed the current brackets at each step, |
| when it returns true then the current brackets are returned as the result.]] |
| [[max_iter] [The maximum number of function invocations to perform in the search |
| for the root.]] |
| ] |
| |
| [optional_policy] |
| |
| Returns: a pair of values /r/ that bracket the root so that: |
| |
| f(r.first) * f(r.second) <= 0 |
| |
| and either |
| |
| tol(r.first, r.second) == true |
| |
| or |
| |
| max_iter >= m |
| |
| where /m/ is the initial value of /max_iter/ passed to the function. |
| |
| In other words, it's up to the caller to verify whether termination occurred |
| as a result of exceeding /max_iter/ function invocations (easily done by |
| checking the value of /max_iter/ when the function returns), rather than |
| because the termination condition /tol/ was satisfied. |
| |
| template <class F, class T, class Tol> |
| std::pair<T, T> |
| toms748_solve( |
| F f, |
| const T& a, |
| const T& b, |
| Tol tol, |
| boost::uintmax_t& max_iter); |
| |
| template <class F, class T, class Tol, class ``__Policy``> |
| std::pair<T, T> |
| toms748_solve( |
| F f, |
| const T& a, |
| const T& b, |
| Tol tol, |
| boost::uintmax_t& max_iter, |
| const ``__Policy``&); |
| |
| template <class F, class T, class Tol> |
| std::pair<T, T> |
| toms748_solve( |
| F f, |
| const T& a, |
| const T& b, |
| const T& fa, |
| const T& fb, |
| Tol tol, |
| boost::uintmax_t& max_iter); |
| |
| template <class F, class T, class Tol, class ``__Policy``> |
| std::pair<T, T> |
| toms748_solve( |
| F f, |
| const T& a, |
| const T& b, |
| const T& fa, |
| const T& fb, |
| Tol tol, |
| boost::uintmax_t& max_iter, |
| const ``__Policy``&); |
| |
| These two functions implement TOMS Algorithm 748: it uses a mixture of |
| cubic, quadratic and linear (secant) interpolation to locate the root of |
| /f(x)/. The two functions differ only by whether values for /f(a)/ and |
| /f(b)/ are already available. The parameters are: |
| |
| [variablelist |
| [[f] [A unary functor that is the function whose root is to be solved. |
| f(x) need not be uniformly increasing or decreasing on /x/ and |
| may have multiple roots.]] |
| [[a] [ The lower bound for the initial bracket of the root.]] |
| [[b] [The upper bound for the initial bracket of the root. |
| It is a precondition that /a < b/ and that /a/ and /b/ |
| bracket the root to find so that /f(a)*f(b) < 0/.]] |
| [[fa] [Optional: the value of /f(a)/.]] |
| [[fb] [Optional: the value of /f(b)/.]] |
| [[tol] [A binary functor that determines the termination condition for the search |
| for the root. /tol/ is passed the current brackets at each step, |
| when it returns true then the current brackets are returned as the result.]] |
| [[max_iter] [The maximum number of function invocations to perform in the search |
| for the root. On exit /max_iter/ is set to actual number of function |
| invocations used.]] |
| ] |
| |
| [optional_policy] |
| |
| Returns: a pair of values /r/ that bracket the root so that: |
| |
| f(r.first) * f(r.second) <= 0 |
| |
| and either |
| |
| tol(r.first, r.second) == true |
| |
| or |
| |
| max_iter >= m |
| |
| where /m/ is the initial value of /max_iter/ passed to the function. |
| |
| In other words, it's up to the caller to verify whether termination occurred |
| as a result of exceeding /max_iter/ function invocations (easily done by |
| checking the value of /max_iter/), rather than because the termination |
| condition /tol/ was satisfied. |
| |
| template <class T> |
| struct eps_tolerance |
| { |
| eps_tolerance(int bits); |
| bool operator()(const T& a, const T& b)const; |
| }; |
| |
| This is the usual termination condition used with these root finding functions. |
| Its operator() will return true when the relative distance between /a/ and /b/ |
| is less than twice the machine epsilon for T, or 2[super 1-bits], whichever is |
| the larger. In other words you set /bits/ to the number of bits of precision you |
| want in the result. The minimal tolerance of twice the machine epsilon of T is |
| required to ensure that we get back a bracketing interval: since this must clearly |
| be at least 1 epsilon in size. |
| |
| struct equal_floor |
| { |
| equal_floor(); |
| template <class T> bool operator()(const T& a, const T& b)const; |
| }; |
| |
| This termination condition is used when you want to find an integer result |
| that is the /floor/ of the true root. It will terminate as soon as both ends |
| of the interval have the same /floor/. |
| |
| struct equal_ceil |
| { |
| equal_ceil(); |
| template <class T> bool operator()(const T& a, const T& b)const; |
| }; |
| |
| This termination condition is used when you want to find an integer result |
| that is the /ceil/ of the true root. It will terminate as soon as both ends |
| of the interval have the same /ceil/. |
| |
| struct equal_nearest_integer |
| { |
| equal_nearest_integer(); |
| template <class T> bool operator()(const T& a, const T& b)const; |
| }; |
| |
| This termination condition is used when you want to find an integer result |
| that is the /closest/ to the true root. It will terminate as soon as both ends |
| of the interval round to the same nearest integer. |
| |
| [h4 Implementation] |
| |
| The implementation of the bisection algorithm is extremely straightforward |
| and not detailed here. TOMS algorithm 748 is described in detail in: |
| |
| ['Algorithm 748: Enclosing Zeros of Continuous Functions, |
| G. E. Alefeld, F. A. Potra and Yixun Shi, |
| ACM Transactions on Mathematica1 Software, Vol. 21. No. 3. September 1995. |
| Pages 327-344.] |
| |
| The implementation here is a faithful translation of this paper into C++. |
| |
| [endsect][/section:roots2 Root Finding Without Derivatives] |
| |
| [/ |
| Copyright 2006 John Maddock and Paul A. Bristow. |
| Distributed under the Boost Software License, Version 1.0. |
| (See accompanying file LICENSE_1_0.txt or copy at |
| http://www.boost.org/LICENSE_1_0.txt). |
| ] |
