diff --git a/sprout/math/fdim.hpp b/sprout/math/fdim.hpp index cfb145be..a991fb89 100644 --- a/sprout/math/fdim.hpp +++ b/sprout/math/fdim.hpp @@ -1,9 +1,11 @@ #ifndef SPROUT_MATH_FDIM_HPP #define SPROUT_MATH_FDIM_HPP +#include #include #include #include +#include #include #include @@ -16,9 +18,19 @@ namespace sprout { > inline SPROUT_CONSTEXPR FloatType fdim(FloatType x, FloatType y) { - return x > y ? x - y : 0; + return sprout::math::isnan(y) ? y + : sprout::math::isnan(x) ? x + : x == std::numeric_limits::infinity() ? std::numeric_limits::infinity() + : x == -std::numeric_limits::infinity() ? FloatType(0) + : y == std::numeric_limits::infinity() ? FloatType(0) + : y == -std::numeric_limits::infinity() ? std::numeric_limits::infinity() +#if SPROUT_USE_BUILTIN_CMATH_FUNCTION + : std::fdim(x, y) +#else + : x > y ? x - y : FloatType(0) +#endif + ; } - template< typename ArithmeticType1, typename ArithmeticType2, @@ -33,7 +45,7 @@ namespace sprout { } } // namespace detail - using NS_SPROUT_MATH_DETAIL::fdim; + using sprout::math::detail::fdim; } // namespace math using sprout::math::fdim; diff --git a/sprout/math/fma.hpp b/sprout/math/fma.hpp index d9b4f15f..a936925e 100644 --- a/sprout/math/fma.hpp +++ b/sprout/math/fma.hpp @@ -1,9 +1,13 @@ #ifndef SPROUT_MATH_FMA_HPP #define SPROUT_MATH_FMA_HPP +#include #include #include #include +#include +#include +#include #include #include @@ -16,9 +20,20 @@ namespace sprout { > inline SPROUT_CONSTEXPR FloatType fma(FloatType x, FloatType y, FloatType z) { - return x * y + z; + return sprout::math::isnan(x) ? x + : sprout::math::isnan(y) ? y + : sprout::math::isnan(z) ? z + : sprout::math::isinf(x) && y == 0 ? std::numeric_limits::quiet_NaN() + : sprout::math::isinf(y) && x == 0 ? std::numeric_limits::quiet_NaN() + : (sprout::math::isinf(x) || sprout::math::isinf(y)) && sprout::math::isinf(z) && (sprout::math::signbit(x) ^ sprout::math::signbit(y) ^ sprout::math::signbit(z)) + ? std::numeric_limits::quiet_NaN() +#if SPROUT_USE_BUILTIN_CMATH_FUNCTION + : std::fma(x, y, z) +#else + : x * y + z +#endif + ; } - template< typename ArithmeticType1, typename ArithmeticType2, @@ -40,7 +55,7 @@ namespace sprout { } } // namespace detail - using NS_SPROUT_MATH_DETAIL::fma; + using sprout::math::detail::fma; } // namespace math using sprout::math::fma; diff --git a/sprout/math/fmax.hpp b/sprout/math/fmax.hpp index f8ddf45e..854ea112 100644 --- a/sprout/math/fmax.hpp +++ b/sprout/math/fmax.hpp @@ -1,6 +1,7 @@ #ifndef SPROUT_MATH_FMAX_HPP #define SPROUT_MATH_FMAX_HPP +#include #include #include #include @@ -17,7 +18,19 @@ namespace sprout { > inline SPROUT_CONSTEXPR FloatType fmax(FloatType x, FloatType y) { - return x < y && !sprout::math::isnan(y) ? y : x; + return sprout::math::isnan(y) ? x + : sprout::math::isnan(x) ? y + : y == -std::numeric_limits::infinity() ? x + : x == std::numeric_limits::infinity() ? x + : x == -std::numeric_limits::infinity() ? y + : y == std::numeric_limits::infinity() ? y +#if SPROUT_USE_BUILTIN_CMATH_FUNCTION + : x == 0 && y == 0 ? x + : std::fmax(x, y) +#else + : x < y ? y : x +#endif + ; } template< typename ArithmeticType1, @@ -33,7 +46,7 @@ namespace sprout { } } // namespace detail - using NS_SPROUT_MATH_DETAIL::fmax; + using sprout::math::detail::fmax; } // namespace math using sprout::math::fmax; diff --git a/sprout/math/fmin.hpp b/sprout/math/fmin.hpp index 8429d925..8cee7a82 100644 --- a/sprout/math/fmin.hpp +++ b/sprout/math/fmin.hpp @@ -1,6 +1,7 @@ #ifndef SPROUT_MATH_FMIN_HPP #define SPROUT_MATH_FMIN_HPP +#include #include #include #include @@ -17,7 +18,19 @@ namespace sprout { > inline SPROUT_CONSTEXPR FloatType fmin(FloatType x, FloatType y) { - return y < x && !sprout::math::isnan(y) ? y : x; + return sprout::math::isnan(y) ? x + : sprout::math::isnan(x) ? y + : x == -std::numeric_limits::infinity() ? x + : y == std::numeric_limits::infinity() ? x + : y == -std::numeric_limits::infinity() ? y + : x == std::numeric_limits::infinity() ? y +#if SPROUT_USE_BUILTIN_CMATH_FUNCTION + : x == 0 && y == 0 ? x + : std::fmin(x, y) +#else + : y < x ? y : x +#endif + ; } template< typename ArithmeticType1, @@ -33,7 +46,7 @@ namespace sprout { } } // namespace detail - using NS_SPROUT_MATH_DETAIL::fmin; + using sprout::math::detail::fmin; } // namespace math using sprout::math::fmin; diff --git a/sprout/math/fmod.hpp b/sprout/math/fmod.hpp index 20df9e84..9d3de506 100644 --- a/sprout/math/fmod.hpp +++ b/sprout/math/fmod.hpp @@ -77,9 +77,9 @@ namespace sprout { : std::numeric_limits::quiet_NaN() : y : sprout::math::isnan(x) ? x - : x == 0 && y != 0 ? x : x == std::numeric_limits::infinity() || x == -std::numeric_limits::infinity() || y == 0 ? -std::numeric_limits::quiet_NaN() + : x == 0 ? x : y == std::numeric_limits::infinity() || y == -std::numeric_limits::infinity() ? x : static_cast(sprout::math::detail::fmod_impl( static_cast::type>(x), diff --git a/sprout/math/greater_equal.hpp b/sprout/math/greater_equal.hpp index 05f5201a..98cb694d 100644 --- a/sprout/math/greater_equal.hpp +++ b/sprout/math/greater_equal.hpp @@ -17,7 +17,7 @@ namespace sprout { >::type = sprout::enabler > inline SPROUT_CONSTEXPR bool - greater_equal_equal(FloatType1 x, FloatType2 y) { + greater_equal(FloatType1 x, FloatType2 y) { return sprout::math::equal_to(x, y) || x > y; } @@ -62,7 +62,7 @@ namespace sprout { } } // namespace detail // - // greater_equal_equal + // greater_equal // template< typename T, typename U, diff --git a/sprout/math/iceil.hpp b/sprout/math/iceil.hpp index 9724613a..22385f76 100644 --- a/sprout/math/iceil.hpp +++ b/sprout/math/iceil.hpp @@ -12,7 +12,7 @@ #if SPROUT_USE_BUILTIN_CMATH_FUNCTION # include #else -# include +# include #endif namespace sprout { @@ -42,7 +42,7 @@ namespace sprout { template inline SPROUT_CONSTEXPR To iceil_impl(FloatType x, To x0) { - return sprout::math::equal_to(x, x0) ? x0 + return sprout::math::less_equal(x, x0) ? x0 : x0 + 1 ; } diff --git a/sprout/math/ifloor.hpp b/sprout/math/ifloor.hpp index 4060e48b..0b470168 100644 --- a/sprout/math/ifloor.hpp +++ b/sprout/math/ifloor.hpp @@ -12,7 +12,7 @@ #if SPROUT_USE_BUILTIN_CMATH_FUNCTION # include #else -# include +# include #endif namespace sprout { @@ -42,7 +42,7 @@ namespace sprout { template inline SPROUT_CONSTEXPR To ifloor_impl(FloatType x, To x0) { - return sprout::math::equal_to(x, x0) ? x0 + return sprout::math::greater_equal(x, x0) ? x0 : x0 - 1 ; } diff --git a/sprout/math/less_equal.hpp b/sprout/math/less_equal.hpp index 9599885d..7f70b572 100644 --- a/sprout/math/less_equal.hpp +++ b/sprout/math/less_equal.hpp @@ -18,7 +18,7 @@ namespace sprout { > inline SPROUT_CONSTEXPR bool less_equal(FloatType1 x, FloatType2 y) { - return sprout::math::not_equal_to(x, y) || x < y; + return sprout::math::equal_to(x, y) || x < y; } template< diff --git a/sprout/math/minmax.hpp b/sprout/math/minmax.hpp index a86929bd..1d173a42 100644 --- a/sprout/math/minmax.hpp +++ b/sprout/math/minmax.hpp @@ -2,8 +2,8 @@ #define SPROUT_MATH_MINMAX_HPP #include +#include #include #include -#include #endif // #ifndef SPROUT_MATH_MINMAX_HPP diff --git a/sprout/math/quotient.hpp b/sprout/math/quotient.hpp index 39cb3262..1a1b1287 100644 --- a/sprout/math/quotient.hpp +++ b/sprout/math/quotient.hpp @@ -3,25 +3,17 @@ #include #include -#include #include #include -#include #include #include -#include +#include #include #include namespace sprout { namespace math { namespace detail { - template - inline SPROUT_CONSTEXPR T - quotient_impl(T x, T y) { - return sprout::math::iround(x / y); - } - template< typename R = int, typename FloatType, @@ -29,15 +21,13 @@ namespace sprout { > inline SPROUT_CONSTEXPR R quotient(FloatType x, FloatType y) { - return sprout::math::isnan(y) || sprout::math::isnan(x) ? R(0) - : x == 0 && y != 0 ? R(0) - : x == std::numeric_limits::infinity() || x == -std::numeric_limits::infinity() || y == 0 + return y == 0 ? R(0) + : sprout::math::isnan(y) || sprout::math::isnan(x) ? R(0) + : x == std::numeric_limits::infinity() || x == -std::numeric_limits::infinity() ? R(0) + : x == 0 ? R(0) : y == std::numeric_limits::infinity() || y == -std::numeric_limits::infinity() ? R(0) - : static_cast(sprout::math::detail::quotient_impl( - static_cast::type>(x), - static_cast::type>(y) - )) + : sprout::math::rem_quo(x, y).second ; } template< diff --git a/sprout/math/rem_quo.hpp b/sprout/math/rem_quo.hpp index c3737a44..5a983027 100644 --- a/sprout/math/rem_quo.hpp +++ b/sprout/math/rem_quo.hpp @@ -3,13 +3,14 @@ #include #include -#include #include #include #include #include #include -#include +#include +#include +#include #include #include #include @@ -26,9 +27,58 @@ namespace sprout { template inline SPROUT_CONSTEXPR sprout::pair - rem_quo_impl(T x, T y, R quo) { + rem_quo_impl_7(T x, T y, int k, R iquo, T x1, T z) { typedef sprout::pair type; - return type(x - quo * y, quo); + return (z > y) || ((z == y) && ((iquo & 1) != 0)) + ? x < 0 ? type(-(x1 - y), ((iquo + 1) & 0x7) * k) : type(x1 - y, ((iquo + 1) & 0x7) * k) + : x < 0 ? type(-x1, (iquo & 0x7) * k) : type(x1, (iquo & 0x7) * k) + ; + } + template + inline SPROUT_CONSTEXPR sprout::pair + rem_quo_impl_6(T x, T y, int k, R iquo, T x1) { + return x1 >= y + ? sprout::math::detail::rem_quo_impl_7(x, y, k, iquo + 1, x1 - y, (x1 - y) + (x1 - y)) + : sprout::math::detail::rem_quo_impl_7(x, y, k, iquo, x1, x1 + x1) + ; + } + template + inline SPROUT_CONSTEXPR sprout::pair + rem_quo_impl_5(T x, T y, int k, R iquo, T x1, T y1, T z, int iscy, int idiff, int i) { + return i != idiff + ? z >= 0 + ? sprout::math::detail::rem_quo_impl_5(x, y, k, iquo + iquo + 2, z + z, y1, z + z - y1, iscy, idiff, i + 1) + : sprout::math::detail::rem_quo_impl_5(x, y, k, iquo + iquo, x1 + x1, y1, x1 + x1 - y1, iscy, idiff, i + 1) + : sprout::math::detail::rem_quo_impl_6(x, y, k, iquo, sprout::math::scalbn(x1, iscy)) + ; + } + template + inline SPROUT_CONSTEXPR sprout::pair + rem_quo_impl_4(T x, T y, int k, R iquo, T x1, T y1, int iscy, int idiff, int i) { + return sprout::math::detail::rem_quo_impl_5(x, y, k, iquo, x1, y1, x1 - y1, iscy, idiff, i); + } + template + inline SPROUT_CONSTEXPR sprout::pair + rem_quo_impl_3(T x, T y, int k, R iquo, T x1, int iscx, int iscy, int idiff) { + return idiff < 0 ? sprout::math::detail::rem_quo_impl_7(x, y, k, iquo, x1, x1 + x1) + : idiff ? sprout::math::detail::rem_quo_impl_4(x, y, k, iquo, sprout::math::scalbn(x1, -iscx), sprout::math::scalbn(y, -iscy), iscy, idiff, 0) + : sprout::math::detail::rem_quo_impl_6(x, y, k, iquo, x1) + ; + } + template + inline SPROUT_CONSTEXPR sprout::pair + rem_quo_impl_2(T x, T y, int k, R iquo, T x1, int iscx, int iscy) { + return sprout::math::detail::rem_quo_impl_3(x, y, k, iquo, x1, iscx, iscy, iscx - iscy); + } + template + inline SPROUT_CONSTEXPR sprout::pair + rem_quo_impl_1(T x, T y, int k, R iquo, T x1) { + return sprout::math::detail::rem_quo_impl_2(x, y, k, iquo, x1, sprout::math::ilogb(x1), sprout::math::ilogb(y)); + } + template + inline SPROUT_CONSTEXPR sprout::pair + rem_quo_impl(T x, T y, R iquo) { + return sprout::math::detail::rem_quo_impl_1(x, sprout::math::fabs(y), (sprout::math::signbit(x) ^ sprout::math::signbit(y) ? -1 : 1), iquo, sprout::math::fabs(x)); } template< @@ -39,19 +89,20 @@ namespace sprout { inline SPROUT_CONSTEXPR sprout::pair rem_quo(FloatType x, FloatType y) { typedef sprout::pair type; - return sprout::math::isnan(y) + return y == 0 ? type(-std::numeric_limits::quiet_NaN(), R(0)) + : sprout::math::isnan(y) ? sprout::math::isnan(x) - ? sprout::math::signbit(y) || sprout::math::signbit(x) ? type(-std::numeric_limits::quiet_NaN(), R(0)) - : type(std::numeric_limits::quiet_NaN(), R(0)) + ? type(sprout::math::signbit(y) && sprout::math::signbit(x) ? -std::numeric_limits::quiet_NaN() : std::numeric_limits::quiet_NaN(), R(0)) : type(y, R(0)) - : x == 0 && y != 0 ? type(FloatType(0), R(0)) - : x == std::numeric_limits::infinity() || x == -std::numeric_limits::infinity() || y == 0 - ? type(std::numeric_limits::quiet_NaN(), R(0)) + : sprout::math::isnan(x) ? type(x, R(0)) + : x == std::numeric_limits::infinity() || x == -std::numeric_limits::infinity() + ? type(-std::numeric_limits::quiet_NaN(), R(0)) + : x == 0 ? type(x, R(0)) : y == std::numeric_limits::infinity() || y == -std::numeric_limits::infinity() ? type(x, R(0)) : sprout::math::detail::rem_quo_ret(sprout::math::detail::rem_quo_impl( static_cast::type>(x), static_cast::type>(y), - sprout::math::quotient(x, y) + R(0) )) ; } @@ -72,7 +123,11 @@ namespace sprout { return sprout::math::detail::rem_quo(static_cast(x), static_cast(y)); } } // namespace detail - + // + // issue: + // rem_quo(-NaN, -NaN) returns {-NaN, _} . + // # returns {+NaN, _} . ( same as rem_quo(+NaN, +NaN) ) + // using sprout::math::detail::rem_quo; } // namespace math diff --git a/sprout/math/remainder.hpp b/sprout/math/remainder.hpp index eca81b88..9455f3d3 100644 --- a/sprout/math/remainder.hpp +++ b/sprout/math/remainder.hpp @@ -3,48 +3,38 @@ #include #include -#include #include #include -#include #include #include -#include +#include #include #include namespace sprout { namespace math { namespace detail { - template - inline SPROUT_CONSTEXPR T - remainder_impl(T x, T y) { - return x - sprout::math::round(x / y) * y; - } - template< typename FloatType, typename sprout::enabler_if::value>::type = sprout::enabler > inline SPROUT_CONSTEXPR FloatType remainder(FloatType x, FloatType y) { - return sprout::math::isnan(y) + return y == 0 ? -std::numeric_limits::quiet_NaN() + : sprout::math::isnan(y) ? sprout::math::isnan(x) - ? sprout::math::signbit(y) || sprout::math::signbit(x) ? -std::numeric_limits::quiet_NaN() + ? sprout::math::signbit(y) && sprout::math::signbit(x) ? -std::numeric_limits::quiet_NaN() : std::numeric_limits::quiet_NaN() : y : sprout::math::isnan(x) ? x - : x == 0 && y != 0 ? x - : x == std::numeric_limits::infinity() || x == -std::numeric_limits::infinity() || y == 0 - ? std::numeric_limits::quiet_NaN() + : x == std::numeric_limits::infinity() || x == -std::numeric_limits::infinity() + ? -std::numeric_limits::quiet_NaN() + : x == 0 ? x : y == std::numeric_limits::infinity() || y == -std::numeric_limits::infinity() ? x #if SPROUT_USE_BUILTIN_CMATH_FUNCTION : std::remainder(x, y) #else - : static_cast(sprout::math::detail::remainder_impl( - static_cast::type>(x), - static_cast::type>(y) - )) + : sprout::math::rem_quo(x, y).first #endif ; } @@ -61,8 +51,12 @@ namespace sprout { return sprout::math::detail::remainder(static_cast(x), static_cast(y)); } } // namespace detail - - using NS_SPROUT_MATH_DETAIL::remainder; + // + // issue: + // remainder(-NaN, -NaN) returns -NaN . + // # returns +NaN . ( same as remainder(+NaN, +NaN) ) + // + using sprout::math::detail::remainder; } // namespace math using sprout::math::remainder;