ceph/src/boost/libs/math/include/boost/math/distributions/detail/generic_quantile.hpp

   1 //  Copyright John Maddock 2008.
   2 //  Use, modification and distribution are subject to the
   3 //  Boost Software License, Version 1.0. (See accompanying file
   4 //  LICENSE_1_0.txt or copy at http://www.boost.org/LICENSE_1_0.txt)
   5
   6 #ifndef BOOST_MATH_DISTIBUTIONS_DETAIL_GENERIC_QUANTILE_HPP
   7 #define BOOST_MATH_DISTIBUTIONS_DETAIL_GENERIC_QUANTILE_HPP
   8
   9 namespace boost{ namespace math{ namespace detail{
  10
  11 template <class Dist>
  12 struct generic_quantile_finder
  13 {
  14    typedef typename Dist::value_type value_type;
  15    typedef typename Dist::policy_type policy_type;
  16
  17    generic_quantile_finder(const Dist& d, value_type t, bool c)
  18       : dist(d), target(t), comp(c) {}
  19
  20    value_type operator()(const value_type& x)
  21    {
  22       return comp ?
  23          value_type(target - cdf(complement(dist, x)))
  24          : value_type(cdf(dist, x) - target);
  25    }
  26
  27 private:
  28    Dist dist;
  29    value_type target;
  30    bool comp;
  31 };
  32
  33 template <class T, class Policy>
  34 inline T check_range_result(const T& x, const Policy& pol, const char* function)
  35 {
  36    if((x >= 0) && (x < tools::min_value<T>()))
  37       return policies::raise_underflow_error<T>(function, 0, pol);
  38    if(x <= -tools::max_value<T>())
  39       return -policies::raise_overflow_error<T>(function, 0, pol);
  40    if(x >= tools::max_value<T>())
  41       return policies::raise_overflow_error<T>(function, 0, pol);
  42    return x;
  43 }
  44
  45 template <class Dist>
  46 typename Dist::value_type generic_quantile(const Dist& dist, const typename Dist::value_type& p, const typename Dist::value_type& guess, bool comp, const char* function)
  47 {
  48    typedef typename Dist::value_type value_type;
  49    typedef typename Dist::policy_type policy_type;
  50    typedef typename policies::normalise<
  51       policy_type,
  52       policies::promote_float<false>,
  53       policies::promote_double<false>,
  54       policies::discrete_quantile<>,
  55       policies::assert_undefined<> >::type forwarding_policy;
  56
  57    //
  58    // Special cases first:
  59    //
  60    if(p == 0)
  61    {
  62       return comp
  63       ? check_range_result(range(dist).second, forwarding_policy(), function)
  64       : check_range_result(range(dist).first, forwarding_policy(), function);
  65    }
  66    if(p == 1)
  67    {
  68       return !comp
  69       ? check_range_result(range(dist).second, forwarding_policy(), function)
  70       : check_range_result(range(dist).first, forwarding_policy(), function);
  71    }
  72
  73    generic_quantile_finder<Dist> f(dist, p, comp);
  74    tools::eps_tolerance<value_type> tol(policies::digits<value_type, forwarding_policy>() - 3);
  75    boost::uintmax_t max_iter = policies::get_max_root_iterations<forwarding_policy>();
  76    std::pair<value_type, value_type> ir = tools::bracket_and_solve_root(
  77       f, guess, value_type(2), true, tol, max_iter, forwarding_policy());
  78    value_type result = ir.first + (ir.second - ir.first) / 2;
  79    if(max_iter >= policies::get_max_root_iterations<forwarding_policy>())
  80    {
  81       return policies::raise_evaluation_error<value_type>(function, "Unable to locate solution in a reasonable time:"
  82          " either there is no answer to quantile"
  83          " or the answer is infinite.  Current best guess is %1%", result, forwarding_policy());
  84    }
  85    return result;
  86 }
  87
  88 }}} // namespaces
  89
  90 #endif // BOOST_MATH_DISTIBUTIONS_DETAIL_GENERIC_QUANTILE_HPP
  91