code:
//  Copyright John Maddock 2007.
//  Use, modification and distribution are subject to 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)

#ifndef BOOST_MATH_POLICY_HPP
#define BOOST_MATH_POLICY_HPP

#include <boost/mpl/list.hpp>
#include <boost/mpl/contains.hpp>
#include <boost/mpl/if.hpp>
#include <boost/mpl/find_if.hpp>
#include <boost/mpl/remove_if.hpp>
#include <boost/mpl/vector.hpp>
#include <boost/mpl/push_back.hpp>
#include <boost/mpl/at.hpp>
#include <boost/mpl/size.hpp>
#include <boost/mpl/comparison.hpp>
#include <boost/type_traits/is_same.hpp>
#include <boost/static_assert.hpp>
#include <boost/assert.hpp>
#include <boost/math/tools/config.hpp>
#include <limits>
// Sadly we do need the .h versions of these to be sure of getting
// FLT_MANT_DIG etc.
#include <limits.h>
#include <stdlib.h>
#include <stddef.h>
#include <math.h>

namespace boost{ namespace math{

namespace tools{

template <class T>
BOOST_MATH_CONSTEXPR int digits(BOOST_MATH_EXPLICIT_TEMPLATE_TYPE(T)) BOOST_NOEXCEPT;
template <class T>
BOOST_MATH_CONSTEXPR T epsilon(BOOST_MATH_EXPLICIT_TEMPLATE_TYPE(T)) BOOST_MATH_NOEXCEPT(T);

}

namespace policies{

//
// Define macros for our default policies, if they're not defined already:
//
// Special cases for exceptions disabled first:
//
#ifdef BOOST_NO_EXCEPTIONS
#  ifndef BOOST_MATH_DOMAIN_ERROR_POLICY
#    define BOOST_MATH_DOMAIN_ERROR_POLICY errno_on_error
#  endif
#  ifndef BOOST_MATH_POLE_ERROR_POLICY
#     define BOOST_MATH_POLE_ERROR_POLICY errno_on_error
#  endif
#  ifndef BOOST_MATH_OVERFLOW_ERROR_POLICY
#     define BOOST_MATH_OVERFLOW_ERROR_POLICY errno_on_error
#  endif
#  ifndef BOOST_MATH_EVALUATION_ERROR_POLICY
#     define BOOST_MATH_EVALUATION_ERROR_POLICY errno_on_error
#  endif
#  ifndef BOOST_MATH_ROUNDING_ERROR_POLICY
#     define BOOST_MATH_ROUNDING_ERROR_POLICY errno_on_error
#  endif
#endif
//
// Then the regular cases:
//
#ifndef BOOST_MATH_DOMAIN_ERROR_POLICY
#define BOOST_MATH_DOMAIN_ERROR_POLICY throw_on_error
#endif
#ifndef BOOST_MATH_POLE_ERROR_POLICY
#define BOOST_MATH_POLE_ERROR_POLICY throw_on_error
#endif
#ifndef BOOST_MATH_OVERFLOW_ERROR_POLICY
#define BOOST_MATH_OVERFLOW_ERROR_POLICY throw_on_error
#endif
#ifndef BOOST_MATH_EVALUATION_ERROR_POLICY
#define BOOST_MATH_EVALUATION_ERROR_POLICY throw_on_error
#endif
#ifndef BOOST_MATH_ROUNDING_ERROR_POLICY
#define BOOST_MATH_ROUNDING_ERROR_POLICY throw_on_error
#endif
#ifndef BOOST_MATH_UNDERFLOW_ERROR_POLICY
#define BOOST_MATH_UNDERFLOW_ERROR_POLICY ignore_error
#endif
#ifndef BOOST_MATH_DENORM_ERROR_POLICY
#define BOOST_MATH_DENORM_ERROR_POLICY ignore_error
#endif
#ifndef BOOST_MATH_INDETERMINATE_RESULT_ERROR_POLICY
#define BOOST_MATH_INDETERMINATE_RESULT_ERROR_POLICY ignore_error
#endif
#ifndef BOOST_MATH_DIGITS10_POLICY
#define BOOST_MATH_DIGITS10_POLICY 0
#endif
#ifndef BOOST_MATH_PROMOTE_FLOAT_POLICY
#define BOOST_MATH_PROMOTE_FLOAT_POLICY true
#endif
#ifndef BOOST_MATH_PROMOTE_DOUBLE_POLICY
#ifdef BOOST_MATH_NO_LONG_DOUBLE_MATH_FUNCTIONS
#define BOOST_MATH_PROMOTE_DOUBLE_POLICY false
#else
#define BOOST_MATH_PROMOTE_DOUBLE_POLICY true
#endif
#endif
#ifndef BOOST_MATH_DISCRETE_QUANTILE_POLICY
#define BOOST_MATH_DISCRETE_QUANTILE_POLICY integer_round_outwards
#endif
#ifndef BOOST_MATH_ASSERT_UNDEFINED_POLICY
#define BOOST_MATH_ASSERT_UNDEFINED_POLICY true
#endif
#ifndef BOOST_MATH_MAX_SERIES_ITERATION_POLICY
#define BOOST_MATH_MAX_SERIES_ITERATION_POLICY 1000000
#endif
#ifndef BOOST_MATH_MAX_ROOT_ITERATION_POLICY
#define BOOST_MATH_MAX_ROOT_ITERATION_POLICY 200
#endif

#if !defined(__BORLANDC__)
#define BOOST_MATH_META_INT(type, name, Default)\
   template <type N = Default> struct name : public boost::mpl::int_<N>{};\
   namespace detail{\
   template <type N>\
   char test_is_valid_arg(const name<N>*);\
   char test_is_default_arg(const name<Default>*);\
   template <class T> struct is_##name##_imp\
   {\
      template <type N> static char test(const name<N>*);\
      static double test(...);\
      BOOST_STATIC_CONSTANT(bool, value = sizeof(test(static_cast<T*>(0))) == 1);\
   };\
   }\
   template <class T> struct is_##name : public boost::mpl::bool_< ::boost::math::policies::detail::is_##name##_imp<T>::value>{};

#define BOOST_MATH_META_BOOL(name, Default)\
   template <bool N = Default> struct name : public boost::mpl::bool_<N>{};\
   namespace detail{\
   template <bool N>\
   char test_is_valid_arg(const name<N>*);\
   char test_is_default_arg(const name<Default>*);\
   template <class T> struct is_##name##_imp\
   {\
      template <bool N> static char test(const name<N>*);\
      static double test(...);\
      BOOST_STATIC_CONSTANT(bool, value = sizeof(test(static_cast<T*>(0))) == 1);\
   };\
   }\
   template <class T> struct is_##name : public boost::mpl::bool_< ::boost::math::policies::detail::is_##name##_imp<T>::value>{};
#else
#define BOOST_MATH_META_INT(Type, name, Default)\
   template <Type N = Default> struct name : public boost::mpl::int_<N>{};\
   namespace detail{\
   template <Type N>\
   char test_is_valid_arg(const name<N>*);\
   char test_is_default_arg(const name<Default>*);\
   template <class T> struct is_##name##_tester\
   {\
      template <Type N> static char test(const name<N>&);\
      static double test(...);\
   };\
   template <class T> struct is_##name##_imp\
   {\
      static T inst;\
      BOOST_STATIC_CONSTANT(bool, value = sizeof( ::boost::math::policies::detail::is_##name##_tester<T>::test(inst)) == 1);\
   };\
   }\
   template <class T> struct is_##name : public boost::mpl::bool_< ::boost::math::policies::detail::is_##name##_imp<T>::value>\
   {\
      template <class U> struct apply{ typedef is_##name<U> type; };\
   };

#define BOOST_MATH_META_BOOL(name, Default)\
   template <bool N = Default> struct name : public boost::mpl::bool_<N>{};\
   namespace detail{\
   template <bool N>\
   char test_is_valid_arg(const name<N>*);\
   char test_is_default_arg(const name<Default>*);\
   template <class T> struct is_##name##_tester\
   {\
      template <bool N> static char test(const name<N>&);\
      static double test(...);\
   };\
   template <class T> struct is_##name##_imp\
   {\
      static T inst;\
      BOOST_STATIC_CONSTANT(bool, value = sizeof( ::boost::math::policies::detail::is_##name##_tester<T>::test(inst)) == 1);\
   };\
   }\
   template <class T> struct is_##name : public boost::mpl::bool_< ::boost::math::policies::detail::is_##name##_imp<T>::value>\
   {\
      template <class U> struct apply{ typedef is_##name<U> type;  };\
   };
#endif
//
// Begin by defining policy types for error handling:
//
enum error_policy_type
{
   throw_on_error = 0,
   errno_on_error = 1,
   ignore_error = 2,
   user_error = 3
};

BOOST_MATH_META_INT(error_policy_type, domain_error, BOOST_MATH_DOMAIN_ERROR_POLICY)
BOOST_MATH_META_INT(error_policy_type, pole_error, BOOST_MATH_POLE_ERROR_POLICY)
BOOST_MATH_META_INT(error_policy_type, overflow_error, BOOST_MATH_OVERFLOW_ERROR_POLICY)
BOOST_MATH_META_INT(error_policy_type, underflow_error, BOOST_MATH_UNDERFLOW_ERROR_POLICY)
BOOST_MATH_META_INT(error_policy_type, denorm_error, BOOST_MATH_DENORM_ERROR_POLICY)
BOOST_MATH_META_INT(error_policy_type, evaluation_error, BOOST_MATH_EVALUATION_ERROR_POLICY)
BOOST_MATH_META_INT(error_policy_type, rounding_error, BOOST_MATH_ROUNDING_ERROR_POLICY)
BOOST_MATH_META_INT(error_policy_type, indeterminate_result_error, BOOST_MATH_INDETERMINATE_RESULT_ERROR_POLICY)

//
// Policy types for internal promotion:
//
BOOST_MATH_META_BOOL(promote_float, BOOST_MATH_PROMOTE_FLOAT_POLICY)
BOOST_MATH_META_BOOL(promote_double, BOOST_MATH_PROMOTE_DOUBLE_POLICY)
BOOST_MATH_META_BOOL(assert_undefined, BOOST_MATH_ASSERT_UNDEFINED_POLICY)
//
// Policy types for discrete quantiles:
//
enum discrete_quantile_policy_type
{
   real,
   integer_round_outwards,
   integer_round_inwards,
   integer_round_down,
   integer_round_up,
   integer_round_nearest
};

BOOST_MATH_META_INT(discrete_quantile_policy_type, discrete_quantile, BOOST_MATH_DISCRETE_QUANTILE_POLICY)
//
// Precision:
//
BOOST_MATH_META_INT(int, digits10, BOOST_MATH_DIGITS10_POLICY)
BOOST_MATH_META_INT(int, digits2, 0)
//
// Iterations:
//
BOOST_MATH_META_INT(unsigned long, max_series_iterations, BOOST_MATH_MAX_SERIES_ITERATION_POLICY)
BOOST_MATH_META_INT(unsigned long, max_root_iterations, BOOST_MATH_MAX_ROOT_ITERATION_POLICY)
//
// Define the names for each possible policy:
//
#define BOOST_MATH_PARAMETER(name)\
   BOOST_PARAMETER_TEMPLATE_KEYWORD(name##_name)\
   BOOST_PARAMETER_NAME(name##_name)

struct default_policy{};

namespace detail{
//
// Trait to work out bits precision from digits10 and digits2:
//
template <class Digits10, class Digits2>
struct precision
{
   //
   // Now work out the precision:
   //
   typedef typename mpl::if_c<
      (Digits10::value == 0),
      digits2<0>,
      digits2<((Digits10::value + 1) * 1000L) / 301L>
   >::type digits2_type;
public:
#ifdef __BORLANDC__
   typedef typename mpl::if_c<
      (Digits2::value > ::boost::math::policies::detail::precision<Digits10,Digits2>::digits2_type::value),
      Digits2, digits2_type>::type type;
#else
   typedef typename mpl::if_c<
      (Digits2::value > digits2_type::value),
      Digits2, digits2_type>::type type;
#endif
};

template <class A, class B, bool b>
struct select_result
{
   typedef A type;
};
template <class A, class B>
struct select_result<A, B, false>
{
   typedef typename mpl::deref<B>::type type;
};

template <class Seq, class Pred, class DefaultType>
struct find_arg
{
private:
   typedef typename mpl::find_if<Seq, Pred>::type iter;
   typedef typename mpl::end<Seq>::type end_type;
public:
   typedef typename select_result<
      DefaultType, iter,
      ::boost::is_same<iter, end_type>::value>::type type;
};

double test_is_valid_arg(...);
double test_is_default_arg(...);
char test_is_valid_arg(const default_policy*);
char test_is_default_arg(const default_policy*);

template <class T>
struct is_valid_policy_imp
{
   BOOST_STATIC_CONSTANT(bool, value = sizeof(::boost::math::policies::detail::test_is_valid_arg(static_cast<T*>(0))) == 1);
};

template <class T>
struct is_default_policy_imp
{
   BOOST_STATIC_CONSTANT(bool, value = sizeof(::boost::math::policies::detail::test_is_default_arg(static_cast<T*>(0))) == 1);
};

template <class T> struct is_valid_policy
: public mpl::bool_<
   ::boost::math::policies::detail::is_valid_policy_imp<T>::value>
{};

template <class T> struct is_default_policy
: public mpl::bool_<
   ::boost::math::policies::detail::is_default_policy_imp<T>::value>
{
   template <class U>
   struct apply
   {
      typedef is_default_policy<U> type;
   };
};

template <class Seq, class T, int N>
struct append_N
{
   typedef typename mpl::push_back<Seq, T>::type new_seq;
   typedef typename append_N<new_seq, T, N-1>::type type;
};

template <class Seq, class T>
struct append_N<Seq, T, 0>
{
   typedef Seq type;
};

//
// Traits class to work out what template parameters our default
// policy<> class will have when modified for forwarding:
//
template <bool f, bool d>
struct default_args
{
   typedef promote_float<false> arg1;
   typedef promote_double<false> arg2;
};

template <>
struct default_args<false, false>
{
   typedef default_policy arg1;
   typedef default_policy arg2;
};

template <>
struct default_args<true, false>
{
   typedef promote_float<false> arg1;
   typedef default_policy arg2;
};

template <>
struct default_args<false, true>
{
   typedef promote_double<false> arg1;
   typedef default_policy arg2;
};

typedef default_args<BOOST_MATH_PROMOTE_FLOAT_POLICY, BOOST_MATH_PROMOTE_DOUBLE_POLICY>::arg1 forwarding_arg1;
typedef default_args<BOOST_MATH_PROMOTE_FLOAT_POLICY, BOOST_MATH_PROMOTE_DOUBLE_POLICY>::arg2 forwarding_arg2;

} // detail
//
// Now define the policy type with enough arguments to handle all
// the policies:
//
template <class A1 = default_policy,
          class A2 = default_policy,
          class A3 = default_policy,
          class A4 = default_policy,
          class A5 = default_policy,
          class A6 = default_policy,
          class A7 = default_policy,
          class A8 = default_policy,
          class A9 = default_policy,
          class A10 = default_policy,
          class A11 = default_policy,
          class A12 = default_policy,
          class A13 = default_policy>
struct policy
{
private:
   //
   // Validate all our arguments:
   //
   BOOST_STATIC_ASSERT(::boost::math::policies::detail::is_valid_policy<A1>::value);
   BOOST_STATIC_ASSERT(::boost::math::policies::detail::is_valid_policy<A2>::value);
   BOOST_STATIC_ASSERT(::boost::math::policies::detail::is_valid_policy<A3>::value);
   BOOST_STATIC_ASSERT(::boost::math::policies::detail::is_valid_policy<A4>::value);
   BOOST_STATIC_ASSERT(::boost::math::policies::detail::is_valid_policy<A5>::value);
   BOOST_STATIC_ASSERT(::boost::math::policies::detail::is_valid_policy<A6>::value);
   BOOST_STATIC_ASSERT(::boost::math::policies::detail::is_valid_policy<A7>::value);
   BOOST_STATIC_ASSERT(::boost::math::policies::detail::is_valid_policy<A8>::value);
   BOOST_STATIC_ASSERT(::boost::math::policies::detail::is_valid_policy<A9>::value);
   BOOST_STATIC_ASSERT(::boost::math::policies::detail::is_valid_policy<A10>::value);
   BOOST_STATIC_ASSERT(::boost::math::policies::detail::is_valid_policy<A11>::value);
   BOOST_STATIC_ASSERT(::boost::math::policies::detail::is_valid_policy<A12>::value);
   BOOST_STATIC_ASSERT(::boost::math::policies::detail::is_valid_policy<A13>::value);
   //
   // Typelist of the arguments:
   //
   typedef mpl::list<A1,A2,A3,A4,A5,A6,A7,A8,A9,A10,A11,A12,A13> arg_list;

public:
   typedef typename detail::find_arg<arg_list, is_domain_error<mpl::_1>, domain_error<> >::type domain_error_type;
   typedef typename detail::find_arg<arg_list, is_pole_error<mpl::_1>, pole_error<> >::type pole_error_type;
   typedef typename detail::find_arg<arg_list, is_overflow_error<mpl::_1>, overflow_error<> >::type overflow_error_type;
   typedef typename detail::find_arg<arg_list, is_underflow_error<mpl::_1>, underflow_error<> >::type underflow_error_type;
   typedef typename detail::find_arg<arg_list, is_denorm_error<mpl::_1>, denorm_error<> >::type denorm_error_type;
   typedef typename detail::find_arg<arg_list, is_evaluation_error<mpl::_1>, evaluation_error<> >::type evaluation_error_type;
   typedef typename detail::find_arg<arg_list, is_rounding_error<mpl::_1>, rounding_error<> >::type rounding_error_type;
   typedef typename detail::find_arg<arg_list, is_indeterminate_result_error<mpl::_1>, indeterminate_result_error<> >::type indeterminate_result_error_type;
private:
   //
   // Now work out the precision:
   //
   typedef typename detail::find_arg<arg_list, is_digits10<mpl::_1>, digits10<> >::type digits10_type;
   typedef typename detail::find_arg<arg_list, is_digits2<mpl::_1>, digits2<> >::type bits_precision_type;
public:
   typedef typename detail::precision<digits10_type, bits_precision_type>::type precision_type;
   //
   // Internal promotion:
   //
   typedef typename detail::find_arg<arg_list, is_promote_float<mpl::_1>, promote_float<> >::type promote_float_type;
   typedef typename detail::find_arg<arg_list, is_promote_double<mpl::_1>, promote_double<> >::type promote_double_type;
   //
   // Discrete quantiles:
   //
   typedef typename detail::find_arg<arg_list, is_discrete_quantile<mpl::_1>, discrete_quantile<> >::type discrete_quantile_type;
   //
   // Mathematically undefined properties:
   //
   typedef typename detail::find_arg<arg_list, is_assert_undefined<mpl::_1>, assert_undefined<> >::type assert_undefined_type;
   //
   // Max iterations:
   //
   typedef typename detail::find_arg<arg_list, is_max_series_iterations<mpl::_1>, max_series_iterations<> >::type max_series_iterations_type;
   typedef typename detail::find_arg<arg_list, is_max_root_iterations<mpl::_1>, max_root_iterations<> >::type max_root_iterations_type;
};
//
// These full specializations are defined to reduce the amount of
// template instantiations that have to take place when using the default
// policies, they have quite a large impact on compile times:
//
template <>
struct policy<default_policy, default_policy, default_policy, default_policy, default_policy, default_policy, default_policy, default_policy, default_policy, default_policy, default_policy>
{
public:
   typedef domain_error<> domain_error_type;
   typedef pole_error<> pole_error_type;
   typedef overflow_error<> overflow_error_type;
   typedef underflow_error<> underflow_error_type;
   typedef denorm_error<> denorm_error_type;
   typedef evaluation_error<> evaluation_error_type;
   typedef rounding_error<> rounding_error_type;
   typedef indeterminate_result_error<> indeterminate_result_error_type;
#if BOOST_MATH_DIGITS10_POLICY == 0
   typedef digits2<> precision_type;
#else
   typedef detail::precision<digits10<>, digits2<> >::type precision_type;
#endif
   typedef promote_float<> promote_float_type;
   typedef promote_double<> promote_double_type;
   typedef discrete_quantile<> discrete_quantile_type;
   typedef assert_undefined<> assert_undefined_type;
   typedef max_series_iterations<> max_series_iterations_type;
   typedef max_root_iterations<> max_root_iterations_type;
};

template <>
struct policy<detail::forwarding_arg1, detail::forwarding_arg2, default_policy, default_policy, default_policy, default_policy, default_policy, default_policy, default_policy, default_policy, default_policy>
{
public:
   typedef domain_error<> domain_error_type;
   typedef pole_error<> pole_error_type;
   typedef overflow_error<> overflow_error_type;
   typedef underflow_error<> underflow_error_type;
   typedef denorm_error<> denorm_error_type;
   typedef evaluation_error<> evaluation_error_type;
   typedef rounding_error<> rounding_error_type;
   typedef indeterminate_result_error<> indeterminate_result_error_type;
#if BOOST_MATH_DIGITS10_POLICY == 0
   typedef digits2<> precision_type;
#else
   typedef detail::precision<digits10<>, digits2<> >::type precision_type;
#endif
   typedef promote_float<false> promote_float_type;
   typedef promote_double<false> promote_double_type;
   typedef discrete_quantile<> discrete_quantile_type;
   typedef assert_undefined<> assert_undefined_type;
   typedef max_series_iterations<> max_series_iterations_type;
   typedef max_root_iterations<> max_root_iterations_type;
};

template <class Policy,
          class A1 = default_policy,
          class A2 = default_policy,
          class A3 = default_policy,
          class A4 = default_policy,
          class A5 = default_policy,
          class A6 = default_policy,
          class A7 = default_policy,
          class A8 = default_policy,
          class A9 = default_policy,
          class A10 = default_policy,
          class A11 = default_policy,
          class A12 = default_policy,
          class A13 = default_policy>
struct normalise
{
private:
   typedef mpl::list<A1,A2,A3,A4,A5,A6,A7,A8,A9,A10,A11,A12,A13> arg_list;
   typedef typename detail::find_arg<arg_list, is_domain_error<mpl::_1>, typename Policy::domain_error_type >::type domain_error_type;
   typedef typename detail::find_arg<arg_list, is_pole_error<mpl::_1>, typename Policy::pole_error_type >::type pole_error_type;
   typedef typename detail::find_arg<arg_list, is_overflow_error<mpl::_1>, typename Policy::overflow_error_type >::type overflow_error_type;
   typedef typename detail::find_arg<arg_list, is_underflow_error<mpl::_1>, typename Policy::underflow_error_type >::type underflow_error_type;
   typedef typename detail::find_arg<arg_list, is_denorm_error<mpl::_1>, typename Policy::denorm_error_type >::type denorm_error_type;
   typedef typename detail::find_arg<arg_list, is_evaluation_error<mpl::_1>, typename Policy::evaluation_error_type >::type evaluation_error_type;
   typedef typename detail::find_arg<arg_list, is_rounding_error<mpl::_1>, typename Policy::rounding_error_type >::type rounding_error_type;
   typedef typename detail::find_arg<arg_list, is_indeterminate_result_error<mpl::_1>, typename Policy::indeterminate_result_error_type >::type indeterminate_result_error_type;
   //
   // Now work out the precision:
   //
   typedef typename detail::find_arg<arg_list, is_digits10<mpl::_1>, digits10<> >::type digits10_type;
   typedef typename detail::find_arg<arg_list, is_digits2<mpl::_1>, typename Policy::precision_type >::type bits_precision_type;
   typedef typename detail::precision<digits10_type, bits_precision_type>::type precision_type;
   //
   // Internal promotion:
   //
   typedef typename detail::find_arg<arg_list, is_promote_float<mpl::_1>, typename Policy::promote_float_type >::type promote_float_type;
   typedef typename detail::find_arg<arg_list, is_promote_double<mpl::_1>, typename Policy::promote_double_type >::type promote_double_type;
   //
   // Discrete quantiles:
   //
   typedef typename detail::find_arg<arg_list, is_discrete_quantile<mpl::_1>, typename Policy::discrete_quantile_type >::type discrete_quantile_type;
   //
   // Mathematically undefined properties:
   //
   typedef typename detail::find_arg<arg_list, is_assert_undefined<mpl::_1>, typename Policy::assert_undefined_type >::type assert_undefined_type;
   //
   // Max iterations:
   //
   typedef typename detail::find_arg<arg_list, is_max_series_iterations<mpl::_1>, typename Policy::max_series_iterations_type>::type max_series_iterations_type;
   typedef typename detail::find_arg<arg_list, is_max_root_iterations<mpl::_1>, typename Policy::max_root_iterations_type>::type max_root_iterations_type;
   //
   // Define a typelist of the policies:
   //
   typedef mpl::vector<
      domain_error_type,
      pole_error_type,
      overflow_error_type,
      underflow_error_type,
      denorm_error_type,
      evaluation_error_type,
      rounding_error_type,
      indeterminate_result_error_type,
      precision_type,
      promote_float_type,
      promote_double_type,
      discrete_quantile_type,
      assert_undefined_type,
      max_series_iterations_type,
      max_root_iterations_type> result_list;
   //
   // Remove all the policies that are the same as the default:
   //
   typedef typename mpl::remove_if<result_list, detail::is_default_policy<mpl::_> >::type reduced_list;
   //
   // Pad out the list with defaults:
   //
   typedef typename detail::append_N<reduced_list, default_policy, (14 - ::boost::mpl::size<reduced_list>::value)>::type result_type;
public:
   typedef policy<
      typename mpl::at<result_type, mpl::int_<0> >::type,
      typename mpl::at<result_type, mpl::int_<1> >::type,
      typename mpl::at<result_type, mpl::int_<2> >::type,
      typename mpl::at<result_type, mpl::int_<3> >::type,
      typename mpl::at<result_type, mpl::int_<4> >::type,
      typename mpl::at<result_type, mpl::int_<5> >::type,
      typename mpl::at<result_type, mpl::int_<6> >::type,
      typename mpl::at<result_type, mpl::int_<7> >::type,
      typename mpl::at<result_type, mpl::int_<8> >::type,
      typename mpl::at<result_type, mpl::int_<9> >::type,
      typename mpl::at<result_type, mpl::int_<10> >::type,
      typename mpl::at<result_type, mpl::int_<11> >::type,
      typename mpl::at<result_type, mpl::int_<12> >::type > type;
};
//
// Full specialisation to speed up compilation of the common case:
//
template <>
struct normalise<policy<>,
          promote_float<false>,
          promote_double<false>,
          discrete_quantile<>,
          assert_undefined<>,
          default_policy,
          default_policy,
          default_policy,
          default_policy,
          default_policy,
          default_policy,
          default_policy>
{
   typedef policy<detail::forwarding_arg1, detail::forwarding_arg2> type;
};

template <>
struct normalise<policy<detail::forwarding_arg1, detail::forwarding_arg2>,
          promote_float<false>,
          promote_double<false>,
          discrete_quantile<>,
          assert_undefined<>,
          default_policy,
          default_policy,
          default_policy,
          default_policy,
          default_policy,
          default_policy,
          default_policy>
{
   typedef policy<detail::forwarding_arg1, detail::forwarding_arg2> type;
};

inline BOOST_MATH_CONSTEXPR policy<> make_policy() BOOST_NOEXCEPT
{ return policy<>(); }

template <class A1>
inline BOOST_MATH_CONSTEXPR typename normalise<policy<>, A1>::type make_policy(const A1&) BOOST_NOEXCEPT
{
   typedef typename normalise<policy<>, A1>::type result_type;
   return result_type();
}

template <class A1, class A2>
inline BOOST_MATH_CONSTEXPR typename normalise<policy<>, A1, A2>::type make_policy(const A1&, const A2&) BOOST_NOEXCEPT
{
   typedef typename normalise<policy<>, A1, A2>::type result_type;
   return result_type();
}

template <class A1, class A2, class A3>
inline BOOST_MATH_CONSTEXPR typename normalise<policy<>, A1, A2, A3>::type make_policy(const A1&, const A2&, const A3&) BOOST_NOEXCEPT
{
   typedef typename normalise<policy<>, A1, A2, A3>::type result_type;
   return result_type();
}

template <class A1, class A2, class A3, class A4>
inline BOOST_MATH_CONSTEXPR typename normalise<policy<>, A1, A2, A3, A4>::type make_policy(const A1&, const A2&, const A3&, const A4&) BOOST_NOEXCEPT
{
   typedef typename normalise<policy<>, A1, A2, A3, A4>::type result_type;
   return result_type();
}

template <class A1, class A2, class A3, class A4, class A5>
inline BOOST_MATH_CONSTEXPR typename normalise<policy<>, A1, A2, A3, A4, A5>::type make_policy(const A1&, const A2&, const A3&, const A4&, const A5&) BOOST_NOEXCEPT
{
   typedef typename normalise<policy<>, A1, A2, A3, A4, A5>::type result_type;
   return result_type();
}

template <class A1, class A2, class A3, class A4, class A5, class A6>
inline BOOST_MATH_CONSTEXPR typename normalise<policy<>, A1, A2, A3, A4, A5, A6>::type make_policy(const A1&, const A2&, const A3&, const A4&, const A5&, const A6&) BOOST_NOEXCEPT
{
   typedef typename normalise<policy<>, A1, A2, A3, A4, A5, A6>::type result_type;
   return result_type();
}

template <class A1, class A2, class A3, class A4, class A5, class A6, class A7>
inline BOOST_MATH_CONSTEXPR typename normalise<policy<>, A1, A2, A3, A4, A5, A6, A7>::type make_policy(const A1&, const A2&, const A3&, const A4&, const A5&, const A6&, const A7&) BOOST_NOEXCEPT
{
   typedef typename normalise<policy<>, A1, A2, A3, A4, A5, A6, A7>::type result_type;
   return result_type();
}

template <class A1, class A2, class A3, class A4, class A5, class A6, class A7, class A8>
inline BOOST_MATH_CONSTEXPR typename normalise<policy<>, A1, A2, A3, A4, A5, A6, A7, A8>::type make_policy(const A1&, const A2&, const A3&, const A4&, const A5&, const A6&, const A7&, const A8&) BOOST_NOEXCEPT
{
   typedef typename normalise<policy<>, A1, A2, A3, A4, A5, A6, A7, A8>::type result_type;
   return result_type();
}

template <class A1, class A2, class A3, class A4, class A5, class A6, class A7, class A8, class A9>
inline BOOST_MATH_CONSTEXPR typename normalise<policy<>, A1, A2, A3, A4, A5, A6, A7, A8, A9>::type make_policy(const A1&, const A2&, const A3&, const A4&, const A5&, const A6&, const A7&, const A8&, const A9&) BOOST_NOEXCEPT
{
   typedef typename normalise<policy<>, A1, A2, A3, A4, A5, A6, A7, A8, A9>::type result_type;
   return result_type();
}

template <class A1, class A2, class A3, class A4, class A5, class A6, class A7, class A8, class A9, class A10>
inline BOOST_MATH_CONSTEXPR typename normalise<policy<>, A1, A2, A3, A4, A5, A6, A7, A8, A9, A10>::type make_policy(const A1&, const A2&, const A3&, const A4&, const A5&, const A6&, const A7&, const A8&, const A9&, const A10&) BOOST_NOEXCEPT
{
   typedef typename normalise<policy<>, A1, A2, A3, A4, A5, A6, A7, A8, A9, A10>::type result_type;
   return result_type();
}

template <class A1, class A2, class A3, class A4, class A5, class A6, class A7, class A8, class A9, class A10, class A11>
inline BOOST_MATH_CONSTEXPR typename normalise<policy<>, A1, A2, A3, A4, A5, A6, A7, A8, A9, A10, A11>::type make_policy(const A1&, const A2&, const A3&, const A4&, const A5&, const A6&, const A7&, const A8&, const A9&, const A10&, const A11&) BOOST_NOEXCEPT
{
   typedef typename normalise<policy<>, A1, A2, A3, A4, A5, A6, A7, A8, A9, A10, A11>::type result_type;
   return result_type();
}

//
// Traits class to handle internal promotion:
//
template <class Real, class Policy>
struct evaluation
{
   typedef Real type;
};

template <class Policy>
struct evaluation<float, Policy>
{
   typedef typename mpl::if_<typename Policy::promote_float_type, double, float>::type type;
};

template <class Policy>
struct evaluation<double, Policy>
{
   typedef typename mpl::if_<typename Policy::promote_double_type, long double, double>::type type;
};

#ifdef BOOST_NO_LIMITS_COMPILE_TIME_CONSTANTS

template <class Real>
struct basic_digits : public mpl::int_<0>{ };
template <>
struct basic_digits<float> : public mpl::int_<FLT_MANT_DIG>{ };
template <>
struct basic_digits<double> : public mpl::int_<DBL_MANT_DIG>{ };
template <>
struct basic_digits<long double> : public mpl::int_<LDBL_MANT_DIG>{ };

template <class Real, class Policy>
struct precision
{
   BOOST_STATIC_ASSERT( ::std::numeric_limits<Real>::radix == 2);
   typedef typename Policy::precision_type precision_type;
   typedef basic_digits<Real> digits_t;
   typedef typename mpl::if_<
      mpl::equal_to<digits_t, mpl::int_<0> >,
      // Possibly unknown precision:
      precision_type,
      typename mpl::if_<
         mpl::or_<mpl::less_equal<digits_t, precision_type>, mpl::less_equal<precision_type, mpl::int_<0> > >,
         // Default case, full precision for RealType:
         digits2< ::std::numeric_limits<Real>::digits>,
         // User customised precision:
         precision_type
      >::type
   >::type type;
};

template <class Policy>
struct precision<float, Policy>
{
   typedef digits2<FLT_MANT_DIG> type;
};
template <class Policy>
struct precision<double, Policy>
{
   typedef digits2<DBL_MANT_DIG> type;
};
template <class Policy>
struct precision<long double, Policy>
{
   typedef digits2<LDBL_MANT_DIG> type;
};

#else

template <class Real, class Policy>
struct precision
{
   BOOST_STATIC_ASSERT((::std::numeric_limits<Real>::radix == 2) || ((::std::numeric_limits<Real>::is_specialized == 0) || (::std::numeric_limits<Real>::digits == 0)));
#ifndef __BORLANDC__
   typedef typename Policy::precision_type precision_type;
   typedef typename mpl::if_c<
      ((::std::numeric_limits<Real>::is_specialized == 0) || (::std::numeric_limits<Real>::digits == 0)),
      // Possibly unknown precision:
      precision_type,
      typename mpl::if_c<
         ((::std::numeric_limits<Real>::digits <= precision_type::value)
         || (Policy::precision_type::value <= 0)),
         // Default case, full precision for RealType:
         digits2< ::std::numeric_limits<Real>::digits>,
         // User customised precision:
         precision_type
      >::type
   >::type type;
#else
   typedef typename Policy::precision_type precision_type;
   typedef mpl::int_< ::std::numeric_limits<Real>::digits> digits_t;
   typedef mpl::bool_< ::std::numeric_limits<Real>::is_specialized> spec_t;
   typedef typename mpl::if_<
      mpl::or_<mpl::equal_to<spec_t, mpl::false_>, mpl::equal_to<digits_t, mpl::int_<0> > >,
      // Possibly unknown precision:
      precision_type,
      typename mpl::if_<
         mpl::or_<mpl::less_equal<digits_t, precision_type>, mpl::less_equal<precision_type, mpl::int_<0> > >,
         // Default case, full precision for RealType:
         digits2< ::std::numeric_limits<Real>::digits>,
         // User customised precision:
         precision_type
      >::type
   >::type type;
#endif
};

#endif

#ifdef BOOST_MATH_USE_FLOAT128

template <class Policy>
struct precision<BOOST_MATH_FLOAT128_TYPE, Policy>
{
   typedef mpl::int_<113> type;
};

#endif

namespace detail{

template <class T, class Policy>
inline BOOST_MATH_CONSTEXPR int digits_imp(mpl::true_ const&) BOOST_NOEXCEPT
{
#ifndef BOOST_NO_LIMITS_COMPILE_TIME_CONSTANTS
   BOOST_STATIC_ASSERT( ::std::numeric_limits<T>::is_specialized);
#else
   BOOST_ASSERT(::std::numeric_limits<T>::is_specialized);
#endif
   typedef typename boost::math::policies::precision<T, Policy>::type p_t;
   return p_t::value;
}

template <class T, class Policy>
inline BOOST_MATH_CONSTEXPR int digits_imp(mpl::false_ const&) BOOST_NOEXCEPT
{
   return tools::digits<T>();
}

} // namespace detail

template <class T, class Policy>
inline BOOST_MATH_CONSTEXPR int digits(BOOST_MATH_EXPLICIT_TEMPLATE_TYPE(T)) BOOST_NOEXCEPT
{
   typedef mpl::bool_< std::numeric_limits<T>::is_specialized > tag_type;
   return detail::digits_imp<T, Policy>(tag_type());
}
template <class T, class Policy>
inline BOOST_MATH_CONSTEXPR int digits_base10(BOOST_MATH_EXPLICIT_TEMPLATE_TYPE(T)) BOOST_NOEXCEPT
{
   return boost::math::policies::digits<T, Policy>() * 301 / 1000L;
}

template <class Policy>
inline BOOST_MATH_CONSTEXPR unsigned long get_max_series_iterations() BOOST_NOEXCEPT
{
   typedef typename Policy::max_series_iterations_type iter_type;
   return iter_type::value;
}

template <class Policy>
inline BOOST_MATH_CONSTEXPR unsigned long get_max_root_iterations() BOOST_NOEXCEPT
{
   typedef typename Policy::max_root_iterations_type iter_type;
   return iter_type::value;
}

namespace detail{

template <class T, class Digits, class Small, class Default>
struct series_factor_calc
{
   static T get() BOOST_MATH_NOEXCEPT(T)
   {
      return ldexp(T(1.0), 1 - Digits::value);
   }
};

template <class T, class Digits>
struct series_factor_calc<T, Digits, mpl::true_, mpl::true_>
{
   static BOOST_MATH_CONSTEXPR T get() BOOST_MATH_NOEXCEPT(T)
   {
      return boost::math::tools::epsilon<T>();
   }
};
template <class T, class Digits>
struct series_factor_calc<T, Digits, mpl::true_, mpl::false_>
{
   static BOOST_MATH_CONSTEXPR T get() BOOST_MATH_NOEXCEPT(T)
   {
      return 1 / static_cast<T>(static_cast<boost::uintmax_t>(1u) << (Digits::value - 1));
   }
};
template <class T, class Digits>
struct series_factor_calc<T, Digits, mpl::false_, mpl::true_>
{
   static BOOST_MATH_CONSTEXPR T get() BOOST_MATH_NOEXCEPT(T)
   {
      return boost::math::tools::epsilon<T>();
   }
};

template <class T, class Policy>
inline BOOST_MATH_CONSTEXPR T get_epsilon_imp(mpl::true_ const&) BOOST_MATH_NOEXCEPT(T)
{
#ifndef BOOST_NO_LIMITS_COMPILE_TIME_CONSTANTS
   BOOST_STATIC_ASSERT( ::std::numeric_limits<T>::is_specialized);
   BOOST_STATIC_ASSERT( ::std::numeric_limits<T>::radix == 2);
#else
   BOOST_ASSERT(::std::numeric_limits<T>::is_specialized);
   BOOST_ASSERT(::std::numeric_limits<T>::radix == 2);
#endif
   typedef typename boost::math::policies::precision<T, Policy>::type p_t;
   typedef mpl::bool_<p_t::value <= std::numeric_limits<boost::uintmax_t>::digits> is_small_int;
   typedef mpl::bool_<p_t::value >= std::numeric_limits<T>::digits> is_default_value;
   return series_factor_calc<T, p_t, is_small_int, is_default_value>::get();
}

template <class T, class Policy>
inline BOOST_MATH_CONSTEXPR T get_epsilon_imp(mpl::false_ const&) BOOST_MATH_NOEXCEPT(T)
{
   return tools::epsilon<T>();
}

} // namespace detail

template <class T, class Policy>
inline BOOST_MATH_CONSTEXPR T get_epsilon(BOOST_MATH_EXPLICIT_TEMPLATE_TYPE(T)) BOOST_MATH_NOEXCEPT(T)
{
   typedef mpl::bool_< (std::numeric_limits<T>::is_specialized && (std::numeric_limits<T>::radix == 2)) > tag_type;
   return detail::get_epsilon_imp<T, Policy>(tag_type());
}

namespace detail{

template <class A1,
          class A2,
          class A3,
          class A4,
          class A5,
          class A6,
          class A7,
          class A8,
          class A9,
          class A10,
          class A11>
char test_is_policy(const policy<A1,A2,A3,A4,A5,A6,A7,A8,A9,A10,A11>*);
double test_is_policy(...);

template <class P>
struct is_policy_imp
{
   BOOST_STATIC_CONSTANT(bool, value = (sizeof(::boost::math::policies::detail::test_is_policy(static_cast<P*>(0))) == 1));
};

}

template <class P>
struct is_policy : public mpl::bool_< ::boost::math::policies::detail::is_policy_imp<P>::value> {};

//
// Helper traits class for distribution error handling:
//
template <class Policy>
struct constructor_error_check
{
   typedef typename Policy::domain_error_type domain_error_type;
   typedef typename mpl::if_c<
      (domain_error_type::value == throw_on_error) || (domain_error_type::value == user_error) || (domain_error_type::value == errno_on_error),
      mpl::true_,
      mpl::false_>::type type;
};

template <class Policy>
struct method_error_check
{
   typedef typename Policy::domain_error_type domain_error_type;
   typedef typename mpl::if_c<
      (domain_error_type::value == throw_on_error) && (domain_error_type::value != user_error),
      mpl::false_,
      mpl::true_>::type type;
};
//
// Does the Policy ever throw on error?
//
template <class Policy>
struct is_noexcept_error_policy
{
   typedef typename Policy::domain_error_type               t1;
   typedef typename Policy::pole_error_type                 t2;
   typedef typename Policy::overflow_error_type             t3;
   typedef typename Policy::underflow_error_type            t4;
   typedef typename Policy::denorm_error_type               t5;
   typedef typename Policy::evaluation_error_type           t6;
   typedef typename Policy::rounding_error_type             t7;
   typedef typename Policy::indeterminate_result_error_type t8;

   BOOST_STATIC_CONSTANT(bool, value =
      ((t1::value != throw_on_error) && (t1::value != user_error)
      && (t2::value != throw_on_error) && (t2::value != user_error)
      && (t3::value != throw_on_error) && (t3::value != user_error)
      && (t4::value != throw_on_error) && (t4::value != user_error)
      && (t5::value != throw_on_error) && (t5::value != user_error)
      && (t6::value != throw_on_error) && (t6::value != user_error)
      && (t7::value != throw_on_error) && (t7::value != user_error)
      && (t8::value != throw_on_error) && (t8::value != user_error)));
};

}}} // namespaces

#endif // BOOST_MATH_POLICY_HPP


code:
#ifndef BOOST_MPL_MINUS_HPP_INCLUDED
#define BOOST_MPL_MINUS_HPP_INCLUDED

// Copyright Aleksey Gurtovoy 2000-2004
//
// 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)
//
// See http://www.boost.org/libs/mpl for documentation.

// $Id$
// $Date$
// $Revision$

#define AUX778076_OP_NAME minus
#define AUX778076_OP_TOKEN -
#include <boost/mpl/aux_/arithmetic_op.hpp>

#endif // BOOST_MPL_MINUS_HPP_INCLUDED

license: mit


code:
///////////////////////////////////////////////////////////////////////////////
// Copyright Christopher Kormanyos 2014.
// Copyright John Maddock 2014.
// Copyright Paul Bristow 2014.
// 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)
//

// Implement quadruple-precision <cmath> support.

#ifndef _BOOST_CSTDFLOAT_CMATH_2014_02_15_HPP_
  #define _BOOST_CSTDFLOAT_CMATH_2014_02_15_HPP_

  #include <boost/math/cstdfloat/cstdfloat_types.hpp>
  #include <boost/math/cstdfloat/cstdfloat_limits.hpp>

  #if defined(BOOST_CSTDFLOAT_HAS_INTERNAL_FLOAT128_T) && defined(BOOST_MATH_USE_FLOAT128) && !defined(BOOST_CSTDFLOAT_NO_LIBQUADMATH_SUPPORT)

  #include <cmath>
  #include <stdexcept>
  #include <boost/cstdint.hpp>
  #include <boost/static_assert.hpp>
  #include <boost/throw_exception.hpp>

  #if defined(_WIN32) && defined(__GNUC__)
    // Several versions of Mingw and probably cygwin too have broken
    // libquadmath implementations that segfault as soon as you call
    // expq or any function that depends on it.
    #define BOOST_CSTDFLOAT_BROKEN_FLOAT128_MATH_FUNCTIONS
  #endif

  // Here is a helper function used for raising the value of a given
  // floating-point type to the power of n, where n has integral type.
  namespace boost { namespace math { namespace cstdfloat { namespace detail {

  template<class float_type, class integer_type>
  inline float_type pown(const float_type& x, const integer_type p)
  {
    const bool isneg  = (x < 0);
    const bool isnan  = (x != x);
    const bool isinf  = ((!isneg) ? bool(+x > (std::numeric_limits<float_type>::max)())
                                  : bool(-x > (std::numeric_limits<float_type>::max)()));

    if(isnan) { return x; }

    if(isinf) { return std::numeric_limits<float_type>::quiet_NaN(); }

    const bool       x_is_neg = (x < 0);
    const float_type abs_x    = (x_is_neg ? -x : x);

    if(p < static_cast<integer_type>(0))
    {
      if(abs_x < (std::numeric_limits<float_type>::min)())
      {
        return (x_is_neg ? -std::numeric_limits<float_type>::infinity()
                         : +std::numeric_limits<float_type>::infinity());
      }
      else
      {
        return float_type(1) / pown(x, static_cast<integer_type>(-p));
      }
    }

    if(p == static_cast<integer_type>(0))
    {
      return float_type(1);
    }
    else
    {
      if(p == static_cast<integer_type>(1)) { return x; }

      if(abs_x > (std::numeric_limits<float_type>::max)())
      {
        return (x_is_neg ? -std::numeric_limits<float_type>::infinity()
                         : +std::numeric_limits<float_type>::infinity());
      }

      if     (p == static_cast<integer_type>(2)) { return  (x * x); }
      else if(p == static_cast<integer_type>(3)) { return ((x * x) * x); }
      else if(p == static_cast<integer_type>(4)) { const float_type x2 = (x * x); return (x2 * x2); }
      else
      {
        // The variable xn stores the binary powers of x.
        float_type result(((p % integer_type(2)) != integer_type(0)) ? x : float_type(1));
        float_type xn    (x);

        integer_type p2 = p;

        while(integer_type(p2 /= 2) != integer_type(0))
        {
          // Square xn for each binary power.
          xn *= xn;

          const bool has_binary_power = (integer_type(p2 % integer_type(2)) != integer_type(0));

          if(has_binary_power)
          {
            // Multiply the result with each binary power contained in the exponent.
            result *= xn;
          }
        }

        return result;
      }
    }
  }

  } } } } // boost::math::cstdfloat::detail

  // We will now define preprocessor symbols representing quadruple-precision <cmath> functions.
  #if defined(BOOST_INTEL)
    #define BOOST_CSTDFLOAT_FLOAT128_LDEXP  __ldexpq
    #define BOOST_CSTDFLOAT_FLOAT128_FREXP  __frexpq
    #define BOOST_CSTDFLOAT_FLOAT128_FABS   __fabsq
    #define BOOST_CSTDFLOAT_FLOAT128_FLOOR  __floorq
    #define BOOST_CSTDFLOAT_FLOAT128_CEIL   __ceilq
    #if !defined(BOOST_CSTDFLOAT_FLOAT128_SQRT)
    #define BOOST_CSTDFLOAT_FLOAT128_SQRT   __sqrtq
    #endif
    #define BOOST_CSTDFLOAT_FLOAT128_TRUNC  __truncq
    #define BOOST_CSTDFLOAT_FLOAT128_EXP    __expq
    #define BOOST_CSTDFLOAT_FLOAT128_EXPM1  __expm1q
    #define BOOST_CSTDFLOAT_FLOAT128_POW    __powq
    #define BOOST_CSTDFLOAT_FLOAT128_LOG    __logq
    #define BOOST_CSTDFLOAT_FLOAT128_LOG10  __log10q
    #define BOOST_CSTDFLOAT_FLOAT128_SIN    __sinq
    #define BOOST_CSTDFLOAT_FLOAT128_COS    __cosq
    #define BOOST_CSTDFLOAT_FLOAT128_TAN    __tanq
    #define BOOST_CSTDFLOAT_FLOAT128_ASIN   __asinq
    #define BOOST_CSTDFLOAT_FLOAT128_ACOS   __acosq
    #define BOOST_CSTDFLOAT_FLOAT128_ATAN   __atanq
    #define BOOST_CSTDFLOAT_FLOAT128_SINH   __sinhq
    #define BOOST_CSTDFLOAT_FLOAT128_COSH   __coshq
    #define BOOST_CSTDFLOAT_FLOAT128_TANH   __tanhq
    #define BOOST_CSTDFLOAT_FLOAT128_ASINH  __asinhq
    #define BOOST_CSTDFLOAT_FLOAT128_ACOSH  __acoshq
    #define BOOST_CSTDFLOAT_FLOAT128_ATANH  __atanhq
    #define BOOST_CSTDFLOAT_FLOAT128_FMOD   __fmodq
    #define BOOST_CSTDFLOAT_FLOAT128_ATAN2  __atan2q
    #define BOOST_CSTDFLOAT_FLOAT128_LGAMMA __lgammaq
    #define BOOST_CSTDFLOAT_FLOAT128_TGAMMA __tgammaq
  #elif defined(__GNUC__)
    #define BOOST_CSTDFLOAT_FLOAT128_LDEXP  ldexpq
    #define BOOST_CSTDFLOAT_FLOAT128_FREXP  frexpq
    #define BOOST_CSTDFLOAT_FLOAT128_FABS   fabsq
    #define BOOST_CSTDFLOAT_FLOAT128_FLOOR  floorq
    #define BOOST_CSTDFLOAT_FLOAT128_CEIL   ceilq
    #if !defined(BOOST_CSTDFLOAT_FLOAT128_SQRT)
    #define BOOST_CSTDFLOAT_FLOAT128_SQRT   sqrtq
    #endif
    #define BOOST_CSTDFLOAT_FLOAT128_TRUNC  truncq
    #define BOOST_CSTDFLOAT_FLOAT128_POW    powq
    #define BOOST_CSTDFLOAT_FLOAT128_LOG    logq
    #define BOOST_CSTDFLOAT_FLOAT128_LOG10  log10q
    #define BOOST_CSTDFLOAT_FLOAT128_SIN    sinq
    #define BOOST_CSTDFLOAT_FLOAT128_COS    cosq
    #define BOOST_CSTDFLOAT_FLOAT128_TAN    tanq
    #define BOOST_CSTDFLOAT_FLOAT128_ASIN   asinq
    #define BOOST_CSTDFLOAT_FLOAT128_ACOS   acosq
    #define BOOST_CSTDFLOAT_FLOAT128_ATAN   atanq
    #define BOOST_CSTDFLOAT_FLOAT128_FMOD   fmodq
    #define BOOST_CSTDFLOAT_FLOAT128_ATAN2  atan2q
    #define BOOST_CSTDFLOAT_FLOAT128_LGAMMA lgammaq
    #if !defined(BOOST_CSTDFLOAT_BROKEN_FLOAT128_MATH_FUNCTIONS)
    #define BOOST_CSTDFLOAT_FLOAT128_EXP    expq
    #define BOOST_CSTDFLOAT_FLOAT128_EXPM1  expm1q_internal
    #define BOOST_CSTDFLOAT_FLOAT128_SINH   sinhq
    #define BOOST_CSTDFLOAT_FLOAT128_COSH   coshq
    #define BOOST_CSTDFLOAT_FLOAT128_TANH   tanhq
    #define BOOST_CSTDFLOAT_FLOAT128_ASINH  asinhq
    #define BOOST_CSTDFLOAT_FLOAT128_ACOSH  acoshq
    #define BOOST_CSTDFLOAT_FLOAT128_ATANH  atanhq
    #define BOOST_CSTDFLOAT_FLOAT128_TGAMMA tgammaq
    #else // BOOST_CSTDFLOAT_BROKEN_FLOAT128_MATH_FUNCTIONS
    #define BOOST_CSTDFLOAT_FLOAT128_EXP    expq_patch
    #define BOOST_CSTDFLOAT_FLOAT128_SINH   sinhq_patch
    #define BOOST_CSTDFLOAT_FLOAT128_COSH   coshq_patch
    #define BOOST_CSTDFLOAT_FLOAT128_TANH   tanhq_patch
    #define BOOST_CSTDFLOAT_FLOAT128_ASINH  asinhq_patch
    #define BOOST_CSTDFLOAT_FLOAT128_ACOSH  acoshq_patch
    #define BOOST_CSTDFLOAT_FLOAT128_ATANH  atanhq_patch
    #define BOOST_CSTDFLOAT_FLOAT128_TGAMMA tgammaq_patch
    #endif // BOOST_CSTDFLOAT_BROKEN_FLOAT128_MATH_FUNCTIONS
  #endif

  // Implement quadruple-precision <cmath> functions in the namespace
  // boost::math::cstdfloat::detail. Subsequently inject these into the
  // std namespace via *using* directive.

  // Begin with some forward function declarations. Also implement patches
  // for compilers that have broken float128 exponential functions.

  extern "C" boost::math::cstdfloat::detail::float_internal128_t BOOST_CSTDFLOAT_FLOAT128_LDEXP (boost::math::cstdfloat::detail::float_internal128_t, int)  throw();
  extern "C" boost::math::cstdfloat::detail::float_internal128_t BOOST_CSTDFLOAT_FLOAT128_FREXP (boost::math::cstdfloat::detail::float_internal128_t, int*) throw();
  extern "C" boost::math::cstdfloat::detail::float_internal128_t BOOST_CSTDFLOAT_FLOAT128_FABS  (boost::math::cstdfloat::detail::float_internal128_t) throw();
  extern "C" boost::math::cstdfloat::detail::float_internal128_t BOOST_CSTDFLOAT_FLOAT128_FLOOR (boost::math::cstdfloat::detail::float_internal128_t) throw();
  extern "C" boost::math::cstdfloat::detail::float_internal128_t BOOST_CSTDFLOAT_FLOAT128_CEIL  (boost::math::cstdfloat::detail::float_internal128_t) throw();
  extern "C" boost::math::cstdfloat::detail::float_internal128_t BOOST_CSTDFLOAT_FLOAT128_SQRT  (boost::math::cstdfloat::detail::float_internal128_t) throw();
  extern "C" boost::math::cstdfloat::detail::float_internal128_t BOOST_CSTDFLOAT_FLOAT128_TRUNC (boost::math::cstdfloat::detail::float_internal128_t) throw();
  extern "C" boost::math::cstdfloat::detail::float_internal128_t BOOST_CSTDFLOAT_FLOAT128_POW   (boost::math::cstdfloat::detail::float_internal128_t, boost::math::cstdfloat::detail::float_internal128_t) throw();
  extern "C" boost::math::cstdfloat::detail::float_internal128_t BOOST_CSTDFLOAT_FLOAT128_LOG   (boost::math::cstdfloat::detail::float_internal128_t) throw();
  extern "C" boost::math::cstdfloat::detail::float_internal128_t BOOST_CSTDFLOAT_FLOAT128_LOG10 (boost::math::cstdfloat::detail::float_internal128_t) throw();
  extern "C" boost::math::cstdfloat::detail::float_internal128_t BOOST_CSTDFLOAT_FLOAT128_SIN   (boost::math::cstdfloat::detail::float_internal128_t) throw();
  extern "C" boost::math::cstdfloat::detail::float_internal128_t BOOST_CSTDFLOAT_FLOAT128_COS   (boost::math::cstdfloat::detail::float_internal128_t) throw();
  extern "C" boost::math::cstdfloat::detail::float_internal128_t BOOST_CSTDFLOAT_FLOAT128_TAN   (boost::math::cstdfloat::detail::float_internal128_t) throw();
  extern "C" boost::math::cstdfloat::detail::float_internal128_t BOOST_CSTDFLOAT_FLOAT128_ASIN  (boost::math::cstdfloat::detail::float_internal128_t) throw();
  extern "C" boost::math::cstdfloat::detail::float_internal128_t BOOST_CSTDFLOAT_FLOAT128_ACOS  (boost::math::cstdfloat::detail::float_internal128_t) throw();
  extern "C" boost::math::cstdfloat::detail::float_internal128_t BOOST_CSTDFLOAT_FLOAT128_ATAN  (boost::math::cstdfloat::detail::float_internal128_t) throw();
  extern "C" boost::math::cstdfloat::detail::float_internal128_t BOOST_CSTDFLOAT_FLOAT128_FMOD  (boost::math::cstdfloat::detail::float_internal128_t, boost::math::cstdfloat::detail::float_internal128_t) throw();
  extern "C" boost::math::cstdfloat::detail::float_internal128_t BOOST_CSTDFLOAT_FLOAT128_ATAN2 (boost::math::cstdfloat::detail::float_internal128_t, boost::math::cstdfloat::detail::float_internal128_t) throw();
  extern "C" boost::math::cstdfloat::detail::float_internal128_t BOOST_CSTDFLOAT_FLOAT128_LGAMMA(boost::math::cstdfloat::detail::float_internal128_t) throw();

  #if !defined(BOOST_CSTDFLOAT_BROKEN_FLOAT128_MATH_FUNCTIONS)

  extern "C" boost::math::cstdfloat::detail::float_internal128_t BOOST_CSTDFLOAT_FLOAT128_EXP   (boost::math::cstdfloat::detail::float_internal128_t x) throw();
  extern "C" boost::math::cstdfloat::detail::float_internal128_t BOOST_CSTDFLOAT_FLOAT128_SINH  (boost::math::cstdfloat::detail::float_internal128_t x) throw();
  extern "C" boost::math::cstdfloat::detail::float_internal128_t BOOST_CSTDFLOAT_FLOAT128_COSH  (boost::math::cstdfloat::detail::float_internal128_t x) throw();
  extern "C" boost::math::cstdfloat::detail::float_internal128_t BOOST_CSTDFLOAT_FLOAT128_TANH  (boost::math::cstdfloat::detail::float_internal128_t x) throw();
  extern "C" boost::math::cstdfloat::detail::float_internal128_t BOOST_CSTDFLOAT_FLOAT128_ASINH (boost::math::cstdfloat::detail::float_internal128_t x) throw();
  extern "C" boost::math::cstdfloat::detail::float_internal128_t BOOST_CSTDFLOAT_FLOAT128_ACOSH (boost::math::cstdfloat::detail::float_internal128_t x) throw();
  extern "C" boost::math::cstdfloat::detail::float_internal128_t BOOST_CSTDFLOAT_FLOAT128_ATANH (boost::math::cstdfloat::detail::float_internal128_t x) throw();
  extern "C" boost::math::cstdfloat::detail::float_internal128_t BOOST_CSTDFLOAT_FLOAT128_TGAMMA(boost::math::cstdfloat::detail::float_internal128_t x) throw();

  #else // BOOST_CSTDFLOAT_BROKEN_FLOAT128_MATH_FUNCTIONS

  // Forward declaration of the patched exponent function, exp(x).
  inline     boost::math::cstdfloat::detail::float_internal128_t BOOST_CSTDFLOAT_FLOAT128_EXP   (boost::math::cstdfloat::detail::float_internal128_t x);

  inline     boost::math::cstdfloat::detail::float_internal128_t BOOST_CSTDFLOAT_FLOAT128_EXPM1 (boost::math::cstdfloat::detail::float_internal128_t x)
  {
    // Compute exp(x) - 1 for x small.

    // Use an order-36 polynomial approximation of the exponential function
    // in the range of (-ln2 < x < ln2). Scale the argument to this range
    // and subsequently multiply the result by 2^n accordingly.

    // Derive the polynomial coefficients with Mathematica(R) by generating
    // a table of high-precision values of exp(x) in the range (-ln2 < x < ln2)
    // and subsequently applying the built-in *Fit* function.

    // Table[{x, Exp[x] - 1}, {x, -Log[2], Log[2], 1/180}]
    // N[%, 120]
    // Fit[%, {x, x^2, x^3, x^4, x^5, x^6, x^7, x^8, x^9, x^10, x^11, x^12,
    //         x^13, x^14, x^15, x^16, x^17, x^18, x^19, x^20, x^21, x^22,
    //         x^23, x^24, x^25, x^26, x^27, x^28, x^29, x^30, x^31, x^32,
    //         x^33, x^34, x^35, x^36}, x]

    typedef boost::math::cstdfloat::detail::float_internal128_t float_type;

    float_type sum;

    if(x > BOOST_FLOAT128_C(0.693147180559945309417232121458176568075500134360255))
    {
      sum = ::BOOST_CSTDFLOAT_FLOAT128_EXP(x) - float_type(1);
    }
    else
    {
      // Compute the polynomial approximation of exp(alpha).
      sum = ((((((((((((((((((((((((((((((((((((  float_type(BOOST_FLOAT128_C(2.69291698127774166063293705964720493864630783729857438187365E-42))  * x
                                                + float_type(BOOST_FLOAT128_C(9.70937085471487654794114679403710456028986572118859594614033E-41))) * x
                                                + float_type(BOOST_FLOAT128_C(3.38715585158055097155585505318085512156885389014410753080500E-39))) * x
                                                + float_type(BOOST_FLOAT128_C(1.15162718532861050809222658798662695267019717760563645440433E-37))) * x
                                                + float_type(BOOST_FLOAT128_C(3.80039074689434663295873584133017767349635602413675471702393E-36))) * x
                                                + float_type(BOOST_FLOAT128_C(1.21612504934087520075905434734158045947460467096773246215239E-34))) * x
                                                + float_type(BOOST_FLOAT128_C(3.76998762883139753126119821241037824830069851253295480396224E-33))) * x
                                                + float_type(BOOST_FLOAT128_C(1.13099628863830344684998293828608215735777107850991029729440E-31))) * x
                                                + float_type(BOOST_FLOAT128_C(3.27988923706982293204067897468714277771890104022419696770352E-30))) * x
                                                + float_type(BOOST_FLOAT128_C(9.18368986379558482800593745627556950089950023355628325088207E-29))) * x
                                                + float_type(BOOST_FLOAT128_C(2.47959626322479746949155352659617642905315302382639380521497E-27))) * x
                                                + float_type(BOOST_FLOAT128_C(6.44695028438447337900255966737803112935639344283098705091949E-26))) * x
                                                + float_type(BOOST_FLOAT128_C(1.61173757109611834904452725462599961406036904573072897122957E-24))) * x
                                                + float_type(BOOST_FLOAT128_C(3.86817017063068403772269360016918092488847584660382953555804E-23))) * x
                                                + float_type(BOOST_FLOAT128_C(8.89679139245057328674891109315654704307721758924206107351744E-22))) * x
                                                + float_type(BOOST_FLOAT128_C(1.95729410633912612308475595397946731738088422488032228717097E-20))) * x
                                                + float_type(BOOST_FLOAT128_C(4.11031762331216485847799061511674191805055663711439605760231E-19))) * x
                                                + float_type(BOOST_FLOAT128_C(8.22063524662432971695598123977873600603370758794431071426640E-18))) * x
                                                + float_type(BOOST_FLOAT128_C(1.56192069685862264622163643500633782667263448653185159383285E-16))) * x
                                                + float_type(BOOST_FLOAT128_C(2.81145725434552076319894558300988749849555291507956994126835E-15))) * x
                                                + float_type(BOOST_FLOAT128_C(4.77947733238738529743820749111754320727153728139716409114011E-14))) * x
                                                + float_type(BOOST_FLOAT128_C(7.64716373181981647590113198578807092707697416852226691068627E-13))) * x
                                                + float_type(BOOST_FLOAT128_C(1.14707455977297247138516979786821056670509688396295740818677E-11))) * x
                                                + float_type(BOOST_FLOAT128_C(1.60590438368216145993923771701549479323291461578567184216302E-10))) * x
                                                + float_type(BOOST_FLOAT128_C(2.08767569878680989792100903212014323125428376052986408239620E-09))) * x
                                                + float_type(BOOST_FLOAT128_C(2.50521083854417187750521083854417187750523408006206780016659E-08))) * x
                                                + float_type(BOOST_FLOAT128_C(2.75573192239858906525573192239858906525573195144226062684604E-07))) * x
                                                + float_type(BOOST_FLOAT128_C(2.75573192239858906525573192239858906525573191310049321957902E-06))) * x
                                                + float_type(BOOST_FLOAT128_C(0.00002480158730158730158730158730158730158730158730149317774)))     * x
                                                + float_type(BOOST_FLOAT128_C(0.00019841269841269841269841269841269841269841269841293575920)))     * x
                                                + float_type(BOOST_FLOAT128_C(0.00138888888888888888888888888888888888888888888888889071045)))     * x
                                                + float_type(BOOST_FLOAT128_C(0.00833333333333333333333333333333333333333333333333332986595)))     * x
                                                + float_type(BOOST_FLOAT128_C(0.04166666666666666666666666666666666666666666666666666664876)))     * x
                                                + float_type(BOOST_FLOAT128_C(0.16666666666666666666666666666666666666666666666666666669048)))     * x
                                                + float_type(BOOST_FLOAT128_C(0.50000000000000000000000000000000000000000000000000000000006)))     * x
                                                + float_type(BOOST_FLOAT128_C(0.99999999999999999999999999999999999999999999999999999999995)))     * x);
    }

    return sum;
  }
  inline     boost::math::cstdfloat::detail::float_internal128_t BOOST_CSTDFLOAT_FLOAT128_EXP   (boost::math::cstdfloat::detail::float_internal128_t x)
  {
    // Patch the expq() function for a subset of broken GCC compilers
    // like GCC 4.7, 4.8 on MinGW.

    // Use an order-36 polynomial approximation of the exponential function
    // in the range of (-ln2 < x < ln2). Scale the argument to this range
    // and subsequently multiply the result by 2^n accordingly.

    // Derive the polynomial coefficients with Mathematica(R) by generating
    // a table of high-precision values of exp(x) in the range (-ln2 < x < ln2)
    // and subsequently applying the built-in *Fit* function.

    // Table[{x, Exp[x] - 1}, {x, -Log[2], Log[2], 1/180}]
    // N[%, 120]
    // Fit[%, {x, x^2, x^3, x^4, x^5, x^6, x^7, x^8, x^9, x^10, x^11, x^12,
    //         x^13, x^14, x^15, x^16, x^17, x^18, x^19, x^20, x^21, x^22,
    //         x^23, x^24, x^25, x^26, x^27, x^28, x^29, x^30, x^31, x^32,
    //         x^33, x^34, x^35, x^36}, x]

    typedef boost::math::cstdfloat::detail::float_internal128_t float_type;

    // Scale the argument x to the range (-ln2 < x < ln2).
    BOOST_CONSTEXPR_OR_CONST float_type one_over_ln2 = float_type(BOOST_FLOAT128_C(1.44269504088896340735992468100189213742664595415299));
    const float_type x_over_ln2   = x * one_over_ln2;

    boost::int_fast32_t n;

    if(x != x)
    {
      // The argument is NaN.
      return std::numeric_limits<float_type>::quiet_NaN();
    }
    else if(::BOOST_CSTDFLOAT_FLOAT128_FABS(x) > BOOST_FLOAT128_C(+0.693147180559945309417232121458176568075500134360255))
    {
      // The absolute value of the argument exceeds ln2.
      n = static_cast<boost::int_fast32_t>(::BOOST_CSTDFLOAT_FLOAT128_FLOOR(x_over_ln2));
    }
    else if(::BOOST_CSTDFLOAT_FLOAT128_FABS(x) < BOOST_FLOAT128_C(+0.693147180559945309417232121458176568075500134360255))
    {
      // The absolute value of the argument is less than ln2.
      n = static_cast<boost::int_fast32_t>(0);
    }
    else
    {
      // The absolute value of the argument is exactly equal to ln2 (in the sense of floating-point equality).
      return float_type(2);
    }

    // Check if the argument is very near an integer.
    const float_type floor_of_x = ::BOOST_CSTDFLOAT_FLOAT128_FLOOR(x);

    if(::BOOST_CSTDFLOAT_FLOAT128_FABS(x - floor_of_x) < float_type(BOOST_CSTDFLOAT_FLOAT128_EPS))
    {
      // Return e^n for arguments very near an integer.
      return boost::math::cstdfloat::detail::pown(BOOST_FLOAT128_C(2.71828182845904523536028747135266249775724709369996), static_cast<boost::int_fast32_t>(floor_of_x));
    }

    // Compute the scaled argument alpha.
    const float_type alpha = x - (n * BOOST_FLOAT128_C(0.693147180559945309417232121458176568075500134360255));

    // Compute the polynomial approximation of expm1(alpha) and add to it
    // in order to obtain the scaled result.
    const float_type scaled_result = ::BOOST_CSTDFLOAT_FLOAT128_EXPM1(alpha) + float_type(1);

    // Rescale the result and return it.
    return scaled_result * boost::math::cstdfloat::detail::pown(float_type(2), n);
  }
  inline     boost::math::cstdfloat::detail::float_internal128_t BOOST_CSTDFLOAT_FLOAT128_SINH  (boost::math::cstdfloat::detail::float_internal128_t x)
  {
    // Patch the sinhq() function for a subset of broken GCC compilers
    // like GCC 4.7, 4.8 on MinGW.
    typedef boost::math::cstdfloat::detail::float_internal128_t float_type;

    // Here, we use the following:
    // Set: ex  = exp(x)
    // Set: em1 = expm1(x)
    // Then
    // sinh(x) = (ex - 1/ex) / 2         ; for |x| >= 1
    // sinh(x) = (2em1 + em1^2) / (2ex)  ; for |x| < 1

    const float_type ex = ::BOOST_CSTDFLOAT_FLOAT128_EXP(x);

    if(::BOOST_CSTDFLOAT_FLOAT128_FABS(x) < float_type(+1))
    {
      const float_type em1 = ::BOOST_CSTDFLOAT_FLOAT128_EXPM1(x);

      return ((em1 * 2) + (em1 * em1)) / (ex * 2);
    }
    else
    {
      return (ex - (float_type(1) / ex)) / 2;
    }
  }
  inline     boost::math::cstdfloat::detail::float_internal128_t BOOST_CSTDFLOAT_FLOAT128_COSH  (boost::math::cstdfloat::detail::float_internal128_t x)
  {
    // Patch the coshq() function for a subset of broken GCC compilers
    // like GCC 4.7, 4.8 on MinGW.
    typedef boost::math::cstdfloat::detail::float_internal128_t float_type;
    const float_type ex = ::BOOST_CSTDFLOAT_FLOAT128_EXP(x);
    return (ex + (float_type(1) / ex)) / 2;
  }
  inline     boost::math::cstdfloat::detail::float_internal128_t BOOST_CSTDFLOAT_FLOAT128_TANH  (boost::math::cstdfloat::detail::float_internal128_t x)
  {
    // Patch the tanhq() function for a subset of broken GCC compilers
    // like GCC 4.7, 4.8 on MinGW.
    typedef boost::math::cstdfloat::detail::float_internal128_t float_type;
    const float_type ex_plus  = ::BOOST_CSTDFLOAT_FLOAT128_EXP(x);
    const float_type ex_minus = (float_type(1) / ex_plus);
    return (ex_plus - ex_minus) / (ex_plus + ex_minus);
  }
  inline     boost::math::cstdfloat::detail::float_internal128_t BOOST_CSTDFLOAT_FLOAT128_ASINH(boost::math::cstdfloat::detail::float_internal128_t x) throw()
  {
    // Patch the asinh() function since quadmath does not have it.
    typedef boost::math::cstdfloat::detail::float_internal128_t float_type;
    return ::BOOST_CSTDFLOAT_FLOAT128_LOG(x + ::BOOST_CSTDFLOAT_FLOAT128_SQRT((x * x) + float_type(1)));
  }
  inline     boost::math::cstdfloat::detail::float_internal128_t BOOST_CSTDFLOAT_FLOAT128_ACOSH(boost::math::cstdfloat::detail::float_internal128_t x) throw()
  {
    // Patch the acosh() function since quadmath does not have it.
    typedef boost::math::cstdfloat::detail::float_internal128_t float_type;
    const float_type zp(x + float_type(1));
    const float_type zm(x - float_type(1));

    return ::BOOST_CSTDFLOAT_FLOAT128_LOG(x + (zp * ::BOOST_CSTDFLOAT_FLOAT128_SQRT(zm / zp)));
  }
  inline     boost::math::cstdfloat::detail::float_internal128_t BOOST_CSTDFLOAT_FLOAT128_ATANH(boost::math::cstdfloat::detail::float_internal128_t x) throw()
  {
    // Patch the atanh() function since quadmath does not have it.
    typedef boost::math::cstdfloat::detail::float_internal128_t float_type;
    return (  ::BOOST_CSTDFLOAT_FLOAT128_LOG(float_type(1) + x)
            - ::BOOST_CSTDFLOAT_FLOAT128_LOG(float_type(1) - x)) / 2;
  }
  inline     boost::math::cstdfloat::detail::float_internal128_t BOOST_CSTDFLOAT_FLOAT128_TGAMMA(boost::math::cstdfloat::detail::float_internal128_t x) throw()
  {
    // Patch the tgammaq() function for a subset of broken GCC compilers
    // like GCC 4.7, 4.8 on MinGW.
    typedef boost::math::cstdfloat::detail::float_internal128_t float_type;

    if(x > float_type(0))
    {
      return ::BOOST_CSTDFLOAT_FLOAT128_EXP(::BOOST_CSTDFLOAT_FLOAT128_LGAMMA(x));
    }
    else if(x < float_type(0))
    {
      // For x < 0, compute tgamma(-x) and use the reflection formula.
      const float_type positive_x          = -x;
            float_type gamma_value         = ::BOOST_CSTDFLOAT_FLOAT128_TGAMMA(positive_x);
      const float_type floor_of_positive_x = ::BOOST_CSTDFLOAT_FLOAT128_FLOOR (positive_x);

      // Take the reflection checks (slightly adapted) from <boost/math/gamma.hpp>.
      const bool floor_of_z_is_equal_to_z = (positive_x == ::BOOST_CSTDFLOAT_FLOAT128_FLOOR(positive_x));

      BOOST_CONSTEXPR_OR_CONST float_type my_pi = BOOST_FLOAT128_C(3.14159265358979323846264338327950288419716939937511);

      if(floor_of_z_is_equal_to_z)
      {
        const bool is_odd = ((boost::int32_t(floor_of_positive_x) % boost::int32_t(2)) != boost::int32_t(0));

        return (is_odd ? -std::numeric_limits<float_type>::infinity()
                       : +std::numeric_limits<float_type>::infinity());
      }

      const float_type sinpx_value = x * ::BOOST_CSTDFLOAT_FLOAT128_SIN(my_pi * x);

      gamma_value *= sinpx_value;

      const bool result_is_too_large_to_represent = (   (::BOOST_CSTDFLOAT_FLOAT128_FABS(gamma_value) < float_type(1))
                                                     && (((std::numeric_limits<float_type>::max)() * ::BOOST_CSTDFLOAT_FLOAT128_FABS(gamma_value)) < my_pi));

      if(result_is_too_large_to_represent)
      {
        const bool is_odd = ((boost::int32_t(floor_of_positive_x) % boost::int32_t(2)) != boost::int32_t(0));

        return (is_odd ? -std::numeric_limits<float_type>::infinity()
                       : +std::numeric_limits<float_type>::infinity());
      }

      gamma_value = -my_pi / gamma_value;

      if((gamma_value > float_type(0)) || (gamma_value < float_type(0)))
      {
        return gamma_value;
      }
      else
      {
        // The value of gamma is too small to represent. Return 0.0 here.
        return float_type(0);
      }
    }
    else
    {
      // Gamma of zero is complex infinity. Return NaN here.
      return std::numeric_limits<float_type>::quiet_NaN();
    }
  }
  #endif // BOOST_CSTDFLOAT_BROKEN_FLOAT128_MATH_FUNCTIONS

  // Define the quadruple-precision <cmath> functions in the namespace boost::math::cstdfloat::detail.

  namespace boost { namespace math { namespace cstdfloat { namespace detail {
  inline   boost::math::cstdfloat::detail::float_internal128_t ldexp (boost::math::cstdfloat::detail::float_internal128_t x, int n)                                                 { return ::BOOST_CSTDFLOAT_FLOAT128_LDEXP (x, n); }
  inline   boost::math::cstdfloat::detail::float_internal128_t frexp (boost::math::cstdfloat::detail::float_internal128_t x, int* pn)                                               { return ::BOOST_CSTDFLOAT_FLOAT128_FREXP (x, pn); }
  inline   boost::math::cstdfloat::detail::float_internal128_t fabs  (boost::math::cstdfloat::detail::float_internal128_t x)                                                        { return ::BOOST_CSTDFLOAT_FLOAT128_FABS  (x); }
  inline   boost::math::cstdfloat::detail::float_internal128_t abs   (boost::math::cstdfloat::detail::float_internal128_t x)                                                        { return ::BOOST_CSTDFLOAT_FLOAT128_FABS  (x); }
  inline   boost::math::cstdfloat::detail::float_internal128_t floor (boost::math::cstdfloat::detail::float_internal128_t x)                                                        { return ::BOOST_CSTDFLOAT_FLOAT128_FLOOR (x); }
  inline   boost::math::cstdfloat::detail::float_internal128_t ceil  (boost::math::cstdfloat::detail::float_internal128_t x)                                                        { return ::BOOST_CSTDFLOAT_FLOAT128_CEIL  (x); }
  inline   boost::math::cstdfloat::detail::float_internal128_t sqrt  (boost::math::cstdfloat::detail::float_internal128_t x)                                                        { return ::BOOST_CSTDFLOAT_FLOAT128_SQRT  (x); }
  inline   boost::math::cstdfloat::detail::float_internal128_t trunc (boost::math::cstdfloat::detail::float_internal128_t x)                                                        { return ::BOOST_CSTDFLOAT_FLOAT128_TRUNC (x); }
  inline   boost::math::cstdfloat::detail::float_internal128_t exp   (boost::math::cstdfloat::detail::float_internal128_t x)                                                        { return ::BOOST_CSTDFLOAT_FLOAT128_EXP   (x); }
  inline   boost::math::cstdfloat::detail::float_internal128_t pow   (boost::math::cstdfloat::detail::float_internal128_t x, boost::math::cstdfloat::detail::float_internal128_t a) { return ::BOOST_CSTDFLOAT_FLOAT128_POW   (x, a); }
  inline   boost::math::cstdfloat::detail::float_internal128_t pow   (boost::math::cstdfloat::detail::float_internal128_t x, int a)                                                 { return ::BOOST_CSTDFLOAT_FLOAT128_POW   (x, boost::math::cstdfloat::detail::float_internal128_t(a)); }
  inline   boost::math::cstdfloat::detail::float_internal128_t log   (boost::math::cstdfloat::detail::float_internal128_t x)                                                        { return ::BOOST_CSTDFLOAT_FLOAT128_LOG   (x); }
  inline   boost::math::cstdfloat::detail::float_internal128_t log10 (boost::math::cstdfloat::detail::float_internal128_t x)                                                        { return ::BOOST_CSTDFLOAT_FLOAT128_LOG10 (x); }
  inline   boost::math::cstdfloat::detail::float_internal128_t sin   (boost::math::cstdfloat::detail::float_internal128_t x)                                                        { return ::BOOST_CSTDFLOAT_FLOAT128_SIN   (x); }
  inline   boost::math::cstdfloat::detail::float_internal128_t cos   (boost::math::cstdfloat::detail::float_internal128_t x)                                                        { return ::BOOST_CSTDFLOAT_FLOAT128_COS   (x); }
  inline   boost::math::cstdfloat::detail::float_internal128_t tan   (boost::math::cstdfloat::detail::float_internal128_t x)                                                        { return ::BOOST_CSTDFLOAT_FLOAT128_TAN   (x); }
  inline   boost::math::cstdfloat::detail::float_internal128_t asin  (boost::math::cstdfloat::detail::float_internal128_t x)                                                        { return ::BOOST_CSTDFLOAT_FLOAT128_ASIN  (x); }
  inline   boost::math::cstdfloat::detail::float_internal128_t acos  (boost::math::cstdfloat::detail::float_internal128_t x)                                                        { return ::BOOST_CSTDFLOAT_FLOAT128_ACOS  (x); }
  inline   boost::math::cstdfloat::detail::float_internal128_t atan  (boost::math::cstdfloat::detail::float_internal128_t x)                                                        { return ::BOOST_CSTDFLOAT_FLOAT128_ATAN  (x); }
  inline   boost::math::cstdfloat::detail::float_internal128_t sinh  (boost::math::cstdfloat::detail::float_internal128_t x)                                                        { return ::BOOST_CSTDFLOAT_FLOAT128_SINH  (x); }
  inline   boost::math::cstdfloat::detail::float_internal128_t cosh  (boost::math::cstdfloat::detail::float_internal128_t x)                                                        { return ::BOOST_CSTDFLOAT_FLOAT128_COSH  (x); }
  inline   boost::math::cstdfloat::detail::float_internal128_t tanh  (boost::math::cstdfloat::detail::float_internal128_t x)                                                        { return ::BOOST_CSTDFLOAT_FLOAT128_TANH  (x); }
  inline   boost::math::cstdfloat::detail::float_internal128_t asinh (boost::math::cstdfloat::detail::float_internal128_t x)                                                        { return ::BOOST_CSTDFLOAT_FLOAT128_ASINH (x); }
  inline   boost::math::cstdfloat::detail::float_internal128_t acosh (boost::math::cstdfloat::detail::float_internal128_t x)                                                        { return ::BOOST_CSTDFLOAT_FLOAT128_ACOSH (x); }
  inline   boost::math::cstdfloat::detail::float_internal128_t atanh (boost::math::cstdfloat::detail::float_internal128_t x)                                                        { return ::BOOST_CSTDFLOAT_FLOAT128_ATANH (x); }
  inline   boost::math::cstdfloat::detail::float_internal128_t fmod  (boost::math::cstdfloat::detail::float_internal128_t a, boost::math::cstdfloat::detail::float_internal128_t b) { return ::BOOST_CSTDFLOAT_FLOAT128_FMOD  (a, b); }
  inline   boost::math::cstdfloat::detail::float_internal128_t atan2 (boost::math::cstdfloat::detail::float_internal128_t y, boost::math::cstdfloat::detail::float_internal128_t x) { return ::BOOST_CSTDFLOAT_FLOAT128_ATAN2 (y, x); }
  inline   boost::math::cstdfloat::detail::float_internal128_t lgamma(boost::math::cstdfloat::detail::float_internal128_t x)                                                        { return ::BOOST_CSTDFLOAT_FLOAT128_LGAMMA(x); }
  inline   boost::math::cstdfloat::detail::float_internal128_t tgamma(boost::math::cstdfloat::detail::float_internal128_t x)                                                        { return ::BOOST_CSTDFLOAT_FLOAT128_TGAMMA(x); }
  } } } } // boost::math::cstdfloat::detail

  // We will now inject the quadruple-precision <cmath> functions
  // into the std namespace. This is done via *using* directive.
  namespace std
  {
    using boost::math::cstdfloat::detail::ldexp;
    using boost::math::cstdfloat::detail::frexp;
    using boost::math::cstdfloat::detail::fabs;
    using boost::math::cstdfloat::detail::abs;
    using boost::math::cstdfloat::detail::floor;
    using boost::math::cstdfloat::detail::ceil;
    using boost::math::cstdfloat::detail::sqrt;
    using boost::math::cstdfloat::detail::trunc;
    using boost::math::cstdfloat::detail::exp;
    using boost::math::cstdfloat::detail::pow;
    using boost::math::cstdfloat::detail::log;
    using boost::math::cstdfloat::detail::log10;
    using boost::math::cstdfloat::detail::sin;
    using boost::math::cstdfloat::detail::cos;
    using boost::math::cstdfloat::detail::tan;
    using boost::math::cstdfloat::detail::asin;
    using boost::math::cstdfloat::detail::acos;
    using boost::math::cstdfloat::detail::atan;
    using boost::math::cstdfloat::detail::sinh;
    using boost::math::cstdfloat::detail::cosh;
    using boost::math::cstdfloat::detail::tanh;
    using boost::math::cstdfloat::detail::asinh;
    using boost::math::cstdfloat::detail::acosh;
    using boost::math::cstdfloat::detail::atanh;
    using boost::math::cstdfloat::detail::fmod;
    using boost::math::cstdfloat::detail::atan2;
    using boost::math::cstdfloat::detail::lgamma;
    using boost::math::cstdfloat::detail::tgamma;
  } // namespace std

  // We will now remove the preprocessor symbols representing quadruple-precision <cmath>
  // functions from the preprocessor.

  #undef BOOST_CSTDFLOAT_FLOAT128_LDEXP
  #undef BOOST_CSTDFLOAT_FLOAT128_FREXP
  #undef BOOST_CSTDFLOAT_FLOAT128_FABS
  #undef BOOST_CSTDFLOAT_FLOAT128_FLOOR
  #undef BOOST_CSTDFLOAT_FLOAT128_CEIL
  #undef BOOST_CSTDFLOAT_FLOAT128_SQRT
  #undef BOOST_CSTDFLOAT_FLOAT128_TRUNC
  #undef BOOST_CSTDFLOAT_FLOAT128_EXP
  #undef BOOST_CSTDFLOAT_FLOAT128_EXPM1
  #undef BOOST_CSTDFLOAT_FLOAT128_POW
  #undef BOOST_CSTDFLOAT_FLOAT128_LOG
  #undef BOOST_CSTDFLOAT_FLOAT128_LOG10
  #undef BOOST_CSTDFLOAT_FLOAT128_SIN
  #undef BOOST_CSTDFLOAT_FLOAT128_COS
  #undef BOOST_CSTDFLOAT_FLOAT128_TAN
  #undef BOOST_CSTDFLOAT_FLOAT128_ASIN
  #undef BOOST_CSTDFLOAT_FLOAT128_ACOS
  #undef BOOST_CSTDFLOAT_FLOAT128_ATAN
  #undef BOOST_CSTDFLOAT_FLOAT128_SINH
  #undef BOOST_CSTDFLOAT_FLOAT128_COSH
  #undef BOOST_CSTDFLOAT_FLOAT128_TANH
  #undef BOOST_CSTDFLOAT_FLOAT128_ASINH
  #undef BOOST_CSTDFLOAT_FLOAT128_ACOSH
  #undef BOOST_CSTDFLOAT_FLOAT128_ATANH
  #undef BOOST_CSTDFLOAT_FLOAT128_FMOD
  #undef BOOST_CSTDFLOAT_FLOAT128_ATAN2
  #undef BOOST_CSTDFLOAT_FLOAT128_LGAMMA
  #undef BOOST_CSTDFLOAT_FLOAT128_TGAMMA

  #endif // Not BOOST_CSTDFLOAT_NO_LIBQUADMATH_SUPPORT (i.e., the user would like to have libquadmath support)

#endif // _BOOST_CSTDFLOAT_CMATH_2014_02_15_HPP_


code:
/*
Copyright Rene Rivera 2008-2015
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)
*/
#if !defined(BOOST_PREDEF_LIBRARY_STD_H) || defined(BOOST_PREDEF_INTERNAL_GENERATE_TESTS)
#ifndef BOOST_PREDEF_LIBRARY_STD_H
#define BOOST_PREDEF_LIBRARY_STD_H
#endif

#include <boost/predef/library/std/_prefix.h>

#include <boost/predef/library/std/cxx.h>
#include <boost/predef/library/std/dinkumware.h>
#include <boost/predef/library/std/libcomo.h>
#include <boost/predef/library/std/modena.h>
#include <boost/predef/library/std/msl.h>
#include <boost/predef/library/std/roguewave.h>
#include <boost/predef/library/std/sgi.h>
#include <boost/predef/library/std/stdcpp3.h>
#include <boost/predef/library/std/stlport.h>
#include <boost/predef/library/std/vacpp.h>

#endif


code:
// Copyright Aleksey Gurtovoy 2000-2004
//
// 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)
//

// Preprocessed version of "boost/mpl/apply.hpp" header
// -- DO NOT modify by hand!

namespace boost { namespace mpl {

template<
      typename F
    >
struct apply0

    : apply_wrap0<
          typename lambda<F>::type

        >
{
    BOOST_MPL_AUX_LAMBDA_SUPPORT(
          1
        , apply0
        , (F )
        )
};

/// workaround for ETI bug
template<>
struct apply0<int>
{
    typedef int type;
};

template<
      typename F, typename T1
    >
struct apply1

    : apply_wrap1<
          typename lambda<F>::type
        , T1
        >
{
    BOOST_MPL_AUX_LAMBDA_SUPPORT(
          2
        , apply1
        , (F, T1)
        )
};

/// workaround for ETI bug
template<>
struct apply1< int,int >
{
    typedef int type;
};

template<
      typename F, typename T1, typename T2
    >
struct apply2

    : apply_wrap2<
          typename lambda<F>::type
        , T1, T2
        >
{
    BOOST_MPL_AUX_LAMBDA_SUPPORT(
          3
        , apply2
        , (F, T1, T2)
        )
};

/// workaround for ETI bug
template<>
struct apply2< int,int,int >
{
    typedef int type;
};

template<
      typename F, typename T1, typename T2, typename T3
    >
struct apply3

    : apply_wrap3<
          typename lambda<F>::type
        , T1, T2, T3
        >
{
    BOOST_MPL_AUX_LAMBDA_SUPPORT(
          4
        , apply3
        , (F, T1, T2, T3)
        )
};

/// workaround for ETI bug
template<>
struct apply3< int,int,int,int >
{
    typedef int type;
};

template<
      typename F, typename T1, typename T2, typename T3, typename T4
    >
struct apply4

    : apply_wrap4<
          typename lambda<F>::type
        , T1, T2, T3, T4
        >
{
    BOOST_MPL_AUX_LAMBDA_SUPPORT(
          5
        , apply4
        , (F, T1, T2, T3, T4)
        )
};

/// workaround for ETI bug
template<>
struct apply4< int,int,int,int,int >
{
    typedef int type;
};

template<
      typename F, typename T1, typename T2, typename T3, typename T4
    , typename T5
    >
struct apply5

    : apply_wrap5<
          typename lambda<F>::type
        , T1, T2, T3, T4, T5
        >
{
    BOOST_MPL_AUX_LAMBDA_SUPPORT(
          6
        , apply5
        , (F, T1, T2, T3, T4, T5)
        )
};

/// workaround for ETI bug
template<>
struct apply5< int,int,int,int,int,int >
{
    typedef int type;
};

}}


code:
/*=============================================================================
    Copyright (c) 2001-2011 Joel de Guzman

    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)
==============================================================================*/
#if !defined(FUSION_VECTOR50_05052005_0207)
#define FUSION_VECTOR50_05052005_0207

#include <boost/fusion/support/config.hpp>
#include <boost/fusion/container/vector/detail/cpp03/vector50_fwd.hpp>
#include <boost/fusion/support/sequence_base.hpp>
#include <boost/fusion/support/is_sequence.hpp>
#include <boost/fusion/support/detail/access.hpp>
#include <boost/fusion/iterator/next.hpp>
#include <boost/fusion/iterator/deref.hpp>
#include <boost/fusion/sequence/intrinsic/begin.hpp>
#include <boost/fusion/container/vector/detail/at_impl.hpp>
#include <boost/fusion/container/vector/detail/value_at_impl.hpp>
#include <boost/fusion/container/vector/detail/begin_impl.hpp>
#include <boost/fusion/container/vector/detail/end_impl.hpp>

#include <boost/mpl/void.hpp>
#include <boost/mpl/int.hpp>
#include <boost/mpl/at.hpp>
#include <boost/mpl/bool.hpp>
#include <boost/mpl/vector/vector50.hpp>
#include <boost/type_traits/is_convertible.hpp>
#include <boost/utility/enable_if.hpp>

#include <boost/preprocessor/dec.hpp>
#include <boost/preprocessor/iteration/iterate.hpp>
#include <boost/preprocessor/repetition/enum.hpp>
#include <boost/preprocessor/repetition/enum_shifted.hpp>
#include <boost/preprocessor/repetition/enum_params.hpp>
#include <boost/preprocessor/repetition/enum_binary_params.hpp>
#include <boost/preprocessor/repetition/repeat_from_to.hpp>

#if !defined(BOOST_FUSION_DONT_USE_PREPROCESSED_FILES)
#include <boost/fusion/container/vector/detail/cpp03/preprocessed/vector50.hpp>
#else
#if defined(__WAVE__) && defined(BOOST_FUSION_CREATE_PREPROCESSED_FILES)
#pragma wave option(preserve: 2, line: 0, output: "preprocessed/vector50.hpp")
#endif

/*=============================================================================
    Copyright (c) 2001-2011 Joel de Guzman

    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)

    This is an auto-generated file. Do not edit!
==============================================================================*/

#if defined(__WAVE__) && defined(BOOST_FUSION_CREATE_PREPROCESSED_FILES)
#pragma wave option(preserve: 1)
#endif

namespace boost { namespace fusion
{
    struct vector_tag;
    struct fusion_sequence_tag;
    struct random_access_traversal_tag;

#define FUSION_HASH #
// expand vector41 to vector50
#define BOOST_PP_FILENAME_1 <boost/fusion/container/vector/detail/cpp03/vector_n.hpp>
#define BOOST_PP_ITERATION_LIMITS (41, 50)
#include BOOST_PP_ITERATE()
#undef FUSION_HASH
}}

#if defined(__WAVE__) && defined(BOOST_FUSION_CREATE_PREPROCESSED_FILES)
#pragma wave option(output: null)
#endif

#endif // BOOST_FUSION_DONT_USE_PREPROCESSED_FILES

#endif