/usr/include/boost/math/special_functions
Edit: /usr/include/boost/math/special_functions/ulp.hpp (3324B)
// (C) Copyright John Maddock 2015.
// 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_SPECIAL_ULP_HPP
#define BOOST_MATH_SPECIAL_ULP_HPP
#ifdef _MSC_VER
#pragma once
#endif
#include
#include
#include
#include
namespace boost{ namespace math{ namespace detail{
template
T ulp_imp(const T& val, const boost::true_type&, const Policy& pol)
{
BOOST_MATH_STD_USING
int expon;
static const char* function = "ulp<%1%>(%1%)";
int fpclass = (boost::math::fpclassify)(val);
if(fpclass == (int)FP_NAN)
{
return policies::raise_domain_error(
function,
"Argument must be finite, but got %1%", val, pol);
}
else if((fpclass == (int)FP_INFINITE) || (fabs(val) >= tools::max_value()))
{
return (val < 0 ? -1 : 1) * policies::raise_overflow_error(function, 0, pol);
}
else if(fpclass == FP_ZERO)
return detail::get_smallest_value();
//
// This code is almost the same as that for float_next, except for negative integers,
// where we preserve the relation ulp(x) == ulp(-x) as does Java:
//
frexp(fabs(val), &expon);
T diff = ldexp(T(1), expon - tools::digits());
if(diff == 0)
diff = detail::get_smallest_value();
return diff;
}
// non-binary version:
template
T ulp_imp(const T& val, const boost::false_type&, const Policy& pol)
{
BOOST_STATIC_ASSERT(std::numeric_limits::is_specialized);
BOOST_STATIC_ASSERT(std::numeric_limits::radix != 2);
BOOST_MATH_STD_USING
int expon;
static const char* function = "ulp<%1%>(%1%)";
int fpclass = (boost::math::fpclassify)(val);
if(fpclass == (int)FP_NAN)
{
return policies::raise_domain_error(
function,
"Argument must be finite, but got %1%", val, pol);
}
else if((fpclass == (int)FP_INFINITE) || (fabs(val) >= tools::max_value()))
{
return (val < 0 ? -1 : 1) * policies::raise_overflow_error(function, 0, pol);
}
else if(fpclass == FP_ZERO)
return detail::get_smallest_value();
//
// This code is almost the same as that for float_next, except for negative integers,
// where we preserve the relation ulp(x) == ulp(-x) as does Java:
//
expon = 1 + ilogb(fabs(val));
T diff = scalbn(T(1), expon - std::numeric_limits::digits);
if(diff == 0)
diff = detail::get_smallest_value();
return diff;
}
}
template
inline typename tools::promote_args::type ulp(const T& val, const Policy& pol)
{
typedef typename tools::promote_args::type result_type;
return detail::ulp_imp(static_cast(val), boost::integral_constant::is_specialized || (std::numeric_limits::radix == 2)>(), pol);
}
template
inline typename tools::promote_args::type ulp(const T& val)
{
return ulp(val, policies::policy<>());
}
}} // namespaces
#endif // BOOST_MATH_SPECIAL_ULP_HPP