ulp.hpp 2.0 KB

1234567891011121314151617181920212223242526272829303132333435363738394041424344454647484950515253545556575859606162636465666768697071
  1. // (C) Copyright John Maddock 2015.
  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. #ifndef BOOST_MATH_SPECIAL_ULP_HPP
  6. #define BOOST_MATH_SPECIAL_ULP_HPP
  7. #ifdef _MSC_VER
  8. #pragma once
  9. #endif
  10. #include <boost/math/special_functions/math_fwd.hpp>
  11. #include <boost/math/policies/error_handling.hpp>
  12. #include <boost/math/special_functions/fpclassify.hpp>
  13. #include <boost/math/special_functions/next.hpp>
  14. namespace boost{ namespace math{ namespace detail{
  15. template <class T, class Policy>
  16. T ulp_imp(const T& val, const Policy& pol)
  17. {
  18. BOOST_MATH_STD_USING
  19. int expon;
  20. static const char* function = "ulp<%1%>(%1%)";
  21. int fpclass = (boost::math::fpclassify)(val);
  22. if(fpclass == (int)FP_NAN)
  23. {
  24. return policies::raise_domain_error<T>(
  25. function,
  26. "Argument must be finite, but got %1%", val, pol);
  27. }
  28. else if((fpclass == (int)FP_INFINITE) || (fabs(val) >= tools::max_value<T>()))
  29. {
  30. return (val < 0 ? -1 : 1) * policies::raise_overflow_error<T>(function, 0, pol);
  31. }
  32. else if(fpclass == FP_ZERO)
  33. return detail::get_smallest_value<T>();
  34. //
  35. // This code is almost the same as that for float_next, except for negative integers,
  36. // where we preserve the relation ulp(x) == ulp(-x) as does Java:
  37. //
  38. frexp(fabs(val), &expon);
  39. T diff = ldexp(T(1), expon - tools::digits<T>());
  40. if(diff == 0)
  41. diff = detail::get_smallest_value<T>();
  42. return diff;
  43. }
  44. }
  45. template <class T, class Policy>
  46. inline typename tools::promote_args<T>::type ulp(const T& val, const Policy& pol)
  47. {
  48. typedef typename tools::promote_args<T>::type result_type;
  49. return detail::ulp_imp(static_cast<result_type>(val), pol);
  50. }
  51. template <class T>
  52. inline typename tools::promote_args<T>::type ulp(const T& val)
  53. {
  54. return ulp(val, policies::policy<>());
  55. }
  56. }} // namespaces
  57. #endif // BOOST_MATH_SPECIAL_ULP_HPP