annotate DEPENDENCIES/generic/include/boost/math/distributions/detail/hypergeometric_cdf.hpp @ 125:34e428693f5d vext

Vext -> Repoint
author Chris Cannam
date Thu, 14 Jun 2018 11:15:39 +0100
parents 2665513ce2d3
children
rev   line source
Chris@16 1 // Copyright 2008 John Maddock
Chris@16 2 //
Chris@16 3 // Use, modification and distribution are subject to the
Chris@16 4 // Boost Software License, Version 1.0.
Chris@16 5 // (See accompanying file LICENSE_1_0.txt
Chris@16 6 // or copy at http://www.boost.org/LICENSE_1_0.txt)
Chris@16 7
Chris@16 8 #ifndef BOOST_MATH_DISTRIBUTIONS_DETAIL_HG_CDF_HPP
Chris@16 9 #define BOOST_MATH_DISTRIBUTIONS_DETAIL_HG_CDF_HPP
Chris@16 10
Chris@16 11 #include <boost/math/policies/error_handling.hpp>
Chris@16 12 #include <boost/math/distributions/detail/hypergeometric_pdf.hpp>
Chris@16 13
Chris@16 14 namespace boost{ namespace math{ namespace detail{
Chris@16 15
Chris@16 16 template <class T, class Policy>
Chris@16 17 T hypergeometric_cdf_imp(unsigned x, unsigned r, unsigned n, unsigned N, bool invert, const Policy& pol)
Chris@16 18 {
Chris@16 19 #ifdef BOOST_MSVC
Chris@16 20 # pragma warning(push)
Chris@16 21 # pragma warning(disable:4267)
Chris@16 22 #endif
Chris@16 23 BOOST_MATH_STD_USING
Chris@16 24 T result = 0;
Chris@16 25 T mode = floor(T(r + 1) * T(n + 1) / (N + 2));
Chris@16 26 if(x < mode)
Chris@16 27 {
Chris@16 28 result = hypergeometric_pdf<T>(x, r, n, N, pol);
Chris@16 29 T diff = result;
Chris@16 30 unsigned lower_limit = static_cast<unsigned>((std::max)(0, (int)(n + r) - (int)(N)));
Chris@16 31 while(diff > (invert ? T(1) : result) * tools::epsilon<T>())
Chris@16 32 {
Chris@16 33 diff = T(x) * T((N + x) - n - r) * diff / (T(1 + n - x) * T(1 + r - x));
Chris@16 34 result += diff;
Chris@16 35 BOOST_MATH_INSTRUMENT_VARIABLE(x);
Chris@16 36 BOOST_MATH_INSTRUMENT_VARIABLE(diff);
Chris@16 37 BOOST_MATH_INSTRUMENT_VARIABLE(result);
Chris@16 38 if(x == lower_limit)
Chris@16 39 break;
Chris@16 40 --x;
Chris@16 41 }
Chris@16 42 }
Chris@16 43 else
Chris@16 44 {
Chris@16 45 invert = !invert;
Chris@16 46 unsigned upper_limit = (std::min)(r, n);
Chris@16 47 if(x != upper_limit)
Chris@16 48 {
Chris@16 49 ++x;
Chris@16 50 result = hypergeometric_pdf<T>(x, r, n, N, pol);
Chris@16 51 T diff = result;
Chris@16 52 while((x <= upper_limit) && (diff > (invert ? T(1) : result) * tools::epsilon<T>()))
Chris@16 53 {
Chris@16 54 diff = T(n - x) * T(r - x) * diff / (T(x + 1) * T((N + x + 1) - n - r));
Chris@16 55 result += diff;
Chris@16 56 ++x;
Chris@16 57 BOOST_MATH_INSTRUMENT_VARIABLE(x);
Chris@16 58 BOOST_MATH_INSTRUMENT_VARIABLE(diff);
Chris@16 59 BOOST_MATH_INSTRUMENT_VARIABLE(result);
Chris@16 60 }
Chris@16 61 }
Chris@16 62 }
Chris@16 63 if(invert)
Chris@16 64 result = 1 - result;
Chris@16 65 return result;
Chris@16 66 #ifdef BOOST_MSVC
Chris@16 67 # pragma warning(pop)
Chris@16 68 #endif
Chris@16 69 }
Chris@16 70
Chris@16 71 template <class T, class Policy>
Chris@16 72 inline T hypergeometric_cdf(unsigned x, unsigned r, unsigned n, unsigned N, bool invert, const Policy&)
Chris@16 73 {
Chris@16 74 BOOST_FPU_EXCEPTION_GUARD
Chris@16 75 typedef typename tools::promote_args<T>::type result_type;
Chris@16 76 typedef typename policies::evaluation<result_type, Policy>::type value_type;
Chris@16 77 typedef typename policies::normalise<
Chris@16 78 Policy,
Chris@16 79 policies::promote_float<false>,
Chris@16 80 policies::promote_double<false>,
Chris@16 81 policies::discrete_quantile<>,
Chris@16 82 policies::assert_undefined<> >::type forwarding_policy;
Chris@16 83
Chris@16 84 value_type result;
Chris@16 85 result = detail::hypergeometric_cdf_imp<value_type>(x, r, n, N, invert, forwarding_policy());
Chris@16 86 if(result > 1)
Chris@16 87 {
Chris@16 88 result = 1;
Chris@16 89 }
Chris@16 90 if(result < 0)
Chris@16 91 {
Chris@16 92 result = 0;
Chris@16 93 }
Chris@16 94 return policies::checked_narrowing_cast<result_type, forwarding_policy>(result, "boost::math::hypergeometric_cdf<%1%>(%1%,%1%,%1%,%1%)");
Chris@16 95 }
Chris@16 96
Chris@16 97 }}} // namespaces
Chris@16 98
Chris@16 99 #endif
Chris@16 100