rounded-hw 0.3.0 → 0.4.0
raw patch · 13 files changed
+390/−173 lines, 13 filesdep ~QuickCheckdep ~arraydep ~base
Dependency ranges changed: QuickCheck, array, base, deepseq, doctest, float128, hspec, long-double, primitive, random, tagged, tasty-bench, vector
Files
- ChangeLog.md +6/−0
- benchmark/Benchmark.hs +1/−1
- benchmark/Conversion.hs +1/−1
- benchmark/IGA.hs +1/−1
- cbits/rounded-common.inl +316/−132
- cbits/rounded.c +14/−0
- rounded-hw.cabal +37/−27
- src/Numeric/Rounded/Hardware/Backend/C.hs +0/−5
- src/Numeric/Rounded/Hardware/Backend/Default.hs +2/−0
- src/Numeric/Rounded/Hardware/Backend/FastFFI.hs +5/−5
- src/Numeric/Rounded/Hardware/Backend/ViaRational.hs +1/−1
- src/Numeric/Rounded/Hardware/Interval.hs +3/−0
- src/Numeric/Rounded/Hardware/Interval/NonEmpty.hs +3/−0
ChangeLog.md view
@@ -1,5 +1,11 @@ # Changelog for rounded-hw +## 0.4.0 (2023-11-18)++* Disable `x87-long-double` flag by default.+* Support GHC 9.4 to 9.8 (inclusive).+* Fix `roundedMul` of `ViaRational`.+ ## 0.3.0 (2022-01-08) * Support Clang on AArch64.
benchmark/Benchmark.hs view
@@ -11,13 +11,13 @@ import Data.Ratio import qualified Data.Vector as V import qualified Data.Vector.Unboxed as VU-import Gauge.Main import IGA import Numeric.Rounded.Hardware.Internal import Numeric.Rounded.Hardware.Interval import Numeric.Rounded.Hardware.Interval.Class (makeInterval) import qualified Numeric.Rounded.Hardware.Interval.NonEmpty as NE import qualified Numeric.Rounded.Hardware.Vector.Unboxed as RVU+import Test.Tasty.Bench foreign import ccall unsafe "fma" c_fma_double :: Double -> Double -> Double -> Double
benchmark/Conversion.hs view
@@ -8,13 +8,13 @@ import Data.Int import Data.Ratio import Data.Word-import Gauge.Benchmark import Numeric.Floating.IEEE import qualified Numeric.Floating.IEEE.Internal as IEEE.Internal import Numeric.Rounded.Hardware import qualified Numeric.Rounded.Hardware.Backend.C as C import Numeric.Rounded.Hardware.Class import Numeric.Rounded.Hardware.Interval+import Test.Tasty.Bench word64ToDouble :: RoundingMode -> Word64 -> Double word64ToDouble ToNearest x
benchmark/IGA.hs view
@@ -12,9 +12,9 @@ import qualified Data.Vector.Mutable as VM import qualified Data.Vector.Unboxed as VU import qualified Data.Vector.Unboxed.Mutable as VUM-import Gauge.Benchmark import Numeric.Rounded.Hardware.Interval import qualified Numeric.Rounded.Hardware.Interval.NonEmpty as NE+import Test.Tasty.Bench thawST :: (Ix i, IArray a e) => a i e -> ST s (STArray s i e) thawST = thaw
cbits/rounded-common.inl view
@@ -4,11 +4,14 @@ // double // -static inline double rounded_add_impl_double(native_rounding_mode mode, double a, double b)+static inline double rounded_add_impl_double(native_rounding_mode mode, VOLATILE double a, VOLATILE double b) { fp_reg oldreg = get_fp_reg(); set_rounding(oldreg, mode);- volatile double c = a + b;+ MAY_MODIFY(a);+ MAY_MODIFY(b);+ VOLATILE double c = a + b;+ FORCE_EVAL(c); restore_fp_reg(oldreg); return c; }@@ -21,11 +24,14 @@ extern double rounded_hw_add_double_zero(double a, double b) { return rounded_add_impl_double(ROUND_TOWARDZERO, a, b); } -static inline double rounded_sub_impl_double(native_rounding_mode mode, double a, double b)+static inline double rounded_sub_impl_double(native_rounding_mode mode, VOLATILE double a, VOLATILE double b) { fp_reg oldreg = get_fp_reg(); set_rounding(oldreg, mode);- volatile double c = a - b;+ MAY_MODIFY(a);+ MAY_MODIFY(b);+ VOLATILE double c = a - b;+ FORCE_EVAL(c); restore_fp_reg(oldreg); return c; }@@ -38,11 +44,14 @@ extern double rounded_hw_sub_double_zero(double a, double b) { return rounded_sub_impl_double(ROUND_TOWARDZERO, a, b); } -static inline double rounded_mul_impl_double(native_rounding_mode mode, double a, double b)+static inline double rounded_mul_impl_double(native_rounding_mode mode, VOLATILE double a, VOLATILE double b) { fp_reg oldreg = get_fp_reg(); set_rounding(oldreg, mode);- volatile double c = a * b;+ MAY_MODIFY(a);+ MAY_MODIFY(b);+ VOLATILE double c = a * b;+ FORCE_EVAL(c); restore_fp_reg(oldreg); return c; }@@ -55,11 +64,14 @@ extern double rounded_hw_mul_double_zero(double a, double b) { return rounded_mul_impl_double(ROUND_TOWARDZERO, a, b); } -static inline double rounded_div_impl_double(native_rounding_mode mode, double a, double b)+static inline double rounded_div_impl_double(native_rounding_mode mode, VOLATILE double a, VOLATILE double b) { fp_reg oldreg = get_fp_reg(); set_rounding(oldreg, mode);- volatile double c = a / b;+ MAY_MODIFY(a);+ MAY_MODIFY(b);+ VOLATILE double c = a / b;+ FORCE_EVAL(c); restore_fp_reg(oldreg); return c; }@@ -72,11 +84,19 @@ extern double rounded_hw_div_double_zero(double a, double b) { return rounded_div_impl_double(ROUND_TOWARDZERO, a, b); } -static inline double rounded_sqrt_impl_double(native_rounding_mode mode, double a)+static inline double rounded_sqrt_impl_double(native_rounding_mode mode, VOLATILE double a) { fp_reg oldreg = get_fp_reg(); set_rounding(oldreg, mode);- volatile double c = sqrt(a);+#if defined(__GNUC__) && defined(USE_SSE2)+ // mingw-w64's sqrt may use x87 instead of SSE2+ double c;+ __asm__ volatile("sqrtsd %1, %0" : "=x"(c) : "x"(a));+#else+ MAY_MODIFY(a);+ VOLATILE double c = sqrt(a);+ FORCE_EVAL(c);+#endif restore_fp_reg(oldreg); return c; }@@ -89,11 +109,15 @@ extern double rounded_hw_sqrt_double_zero(double a) { return rounded_sqrt_impl_double(ROUND_TOWARDZERO, a); } -static inline double rounded_fma_impl_double(native_rounding_mode mode, double a, double b, double c)+static inline double rounded_fma_impl_double(native_rounding_mode mode, VOLATILE double a, VOLATILE double b, VOLATILE double c) { fp_reg oldreg = get_fp_reg(); set_rounding(oldreg, mode);- volatile double result = fma(a, b, c);+ MAY_MODIFY(a);+ MAY_MODIFY(b);+ MAY_MODIFY(c);+ VOLATILE double result = fma(a, b, c);+ FORCE_EVAL(result); restore_fp_reg(oldreg); return result; }@@ -106,15 +130,19 @@ extern double rounded_hw_fma_double_zero(double a, double b, double c) { return rounded_fma_impl_double(ROUND_TOWARDZERO, a, b, c); } -static inline double rounded_fma_if_fast_impl_double(native_rounding_mode mode, double a, double b, double c)+static inline double rounded_fma_if_fast_impl_double(native_rounding_mode mode, VOLATILE double a, VOLATILE double b, VOLATILE double c) { fp_reg oldreg = get_fp_reg(); set_rounding(oldreg, mode);+ MAY_MODIFY(a);+ MAY_MODIFY(b);+ MAY_MODIFY(c); #ifdef FP_FAST_FMA- volatile double result = fma(a, b, c);+ VOLATILE double result = fma(a, b, c); #else- volatile double result = a * b + c;+ VOLATILE double result = a * b + c; #endif+ FORCE_EVAL(result); restore_fp_reg(oldreg); return result; }@@ -135,7 +163,8 @@ { fp_reg oldreg = get_fp_reg(); set_rounding(oldreg, mode);- volatile double result = (double)x;+ VOLATILE double result = (double)x;+ FORCE_EVAL(result); restore_fp_reg(oldreg); return result; }@@ -152,7 +181,8 @@ { fp_reg oldreg = get_fp_reg(); set_rounding(oldreg, mode);- volatile double result = (double)x;+ VOLATILE double result = (double)x;+ FORCE_EVAL(result); restore_fp_reg(oldreg); return result; }@@ -188,138 +218,182 @@ return fast_fmin_double(fast_fmin_double(x, y), fast_fmin_double(z, w)); } -extern double rounded_hw_interval_mul_double_up(double lo1, double hi1, double lo2, double hi2)+extern double rounded_hw_interval_mul_double_up(VOLATILE double lo1, VOLATILE double hi1, VOLATILE double lo2, VOLATILE double hi2) { fp_reg oldreg = get_fp_reg(); set_rounding(oldreg, ROUND_UPWARD);- double x = (volatile double)(lo1 * lo2);- double y = (volatile double)(lo1 * hi2);- double z = (volatile double)(hi1 * lo2);- double w = (volatile double)(hi1 * hi2);+ MAY_MODIFY(lo1);+ MAY_MODIFY(hi1);+ MAY_MODIFY(lo2);+ MAY_MODIFY(hi2);+ VOLATILE double x = lo1 * lo2;+ VOLATILE double y = lo1 * hi2;+ VOLATILE double z = hi1 * lo2;+ VOLATILE double w = hi1 * hi2; if (isnan(x)) x = 0.0; /* 0 * inf -> 0 */ if (isnan(y)) y = 0.0; /* 0 * inf -> 0 */ if (isnan(z)) z = 0.0; /* 0 * inf -> 0 */ if (isnan(w)) w = 0.0; /* 0 * inf -> 0 */- double hi = fast_fmax4_double(x, y, z, w);+ VOLATILE double hi = fast_fmax4_double(x, y, z, w);+ FORCE_EVAL(hi); restore_fp_reg(oldreg); return hi; } -extern double rounded_hw_interval_mul_double_down(double lo1, double hi1, double lo2, double hi2)+extern double rounded_hw_interval_mul_double_down(VOLATILE double lo1, VOLATILE double hi1, VOLATILE double lo2, VOLATILE double hi2) { fp_reg oldreg = get_fp_reg(); set_rounding(oldreg, ROUND_DOWNWARD);- double x = (volatile double)(lo1 * lo2);- double y = (volatile double)(lo1 * hi2);- double z = (volatile double)(hi1 * lo2);- double w = (volatile double)(hi1 * hi2);+ MAY_MODIFY(lo1);+ MAY_MODIFY(hi1);+ MAY_MODIFY(lo2);+ MAY_MODIFY(hi2);+ VOLATILE double x = lo1 * lo2;+ VOLATILE double y = lo1 * hi2;+ VOLATILE double z = hi1 * lo2;+ VOLATILE double w = hi1 * hi2; if (isnan(x)) x = 0.0; /* 0 * inf -> 0 */ if (isnan(y)) y = 0.0; /* 0 * inf -> 0 */ if (isnan(z)) z = 0.0; /* 0 * inf -> 0 */ if (isnan(w)) w = 0.0; /* 0 * inf -> 0 */- double lo = fast_fmin4_double(x, y, z, w);+ VOLATILE double lo = fast_fmin4_double(x, y, z, w);+ FORCE_EVAL(lo); restore_fp_reg(oldreg); return lo; } -extern double rounded_hw_interval_mul_add_double_up(double lo1, double hi1, double lo2, double hi2, double hi3)+extern double rounded_hw_interval_mul_add_double_up(VOLATILE double lo1, VOLATILE double hi1, VOLATILE double lo2, VOLATILE double hi2, VOLATILE double hi3) { fp_reg oldreg = get_fp_reg(); set_rounding(oldreg, ROUND_UPWARD);- double x = (volatile double)(lo1 * lo2);- double y = (volatile double)(lo1 * hi2);- double z = (volatile double)(hi1 * lo2);- double w = (volatile double)(hi1 * hi2);+ MAY_MODIFY(lo1);+ MAY_MODIFY(hi1);+ MAY_MODIFY(lo2);+ MAY_MODIFY(hi2);+ MAY_MODIFY(hi3);+ double x = lo1 * lo2;+ double y = lo1 * hi2;+ double z = hi1 * lo2;+ double w = hi1 * hi2; if (isnan(x)) x = 0.0; /* 0 * inf -> 0 */ if (isnan(y)) y = 0.0; /* 0 * inf -> 0 */ if (isnan(z)) z = 0.0; /* 0 * inf -> 0 */ if (isnan(w)) w = 0.0; /* 0 * inf -> 0 */- volatile double hi = fast_fmax4_double(x, y, z, w) + hi3;+ VOLATILE double hi = fast_fmax4_double(x, y, z, w) + hi3;+ FORCE_EVAL(hi); restore_fp_reg(oldreg); return hi; } -extern double rounded_hw_interval_mul_add_double_down(double lo1, double hi1, double lo2, double hi2, double lo3)+extern double rounded_hw_interval_mul_add_double_down(VOLATILE double lo1, VOLATILE double hi1, VOLATILE double lo2, VOLATILE double hi2, VOLATILE double lo3) { fp_reg oldreg = get_fp_reg(); set_rounding(oldreg, ROUND_DOWNWARD);- double x = (volatile double)(lo1 * lo2);- double y = (volatile double)(lo1 * hi2);- double z = (volatile double)(hi1 * lo2);- double w = (volatile double)(hi1 * hi2);+ MAY_MODIFY(lo1);+ MAY_MODIFY(hi1);+ MAY_MODIFY(lo2);+ MAY_MODIFY(hi2);+ MAY_MODIFY(lo3);+ double x = lo1 * lo2;+ double y = lo1 * hi2;+ double z = hi1 * lo2;+ double w = hi1 * hi2; if (isnan(x)) x = 0.0; /* 0 * inf -> 0 */ if (isnan(y)) y = 0.0; /* 0 * inf -> 0 */ if (isnan(z)) z = 0.0; /* 0 * inf -> 0 */ if (isnan(w)) w = 0.0; /* 0 * inf -> 0 */- volatile double lo = fast_fmin4_double(x, y, z, w) + lo3;+ VOLATILE double lo = fast_fmin4_double(x, y, z, w) + lo3;+ FORCE_EVAL(lo); restore_fp_reg(oldreg); return lo; } -extern double rounded_hw_interval_div_double_up(double lo1, double hi1, double lo2, double hi2)+extern double rounded_hw_interval_div_double_up(VOLATILE double lo1, VOLATILE double hi1, VOLATILE double lo2, VOLATILE double hi2) { fp_reg oldreg = get_fp_reg(); set_rounding(oldreg, ROUND_UPWARD);- double x = (volatile double)(lo1 / lo2);- double y = (volatile double)(lo1 / hi2);- double z = (volatile double)(hi1 / lo2);- double w = (volatile double)(hi1 / hi2);+ MAY_MODIFY(lo1);+ MAY_MODIFY(hi1);+ MAY_MODIFY(lo2);+ MAY_MODIFY(hi2);+ double x = lo1 / lo2;+ double y = lo1 / hi2;+ double z = hi1 / lo2;+ double w = hi1 / hi2; if (isnan(x)) x = 0.0; /* 0 / 0, +-inf / +-inf -> 0 */ if (isnan(y)) y = 0.0; /* 0 / 0, +-inf / +-inf -> 0 */ if (isnan(z)) z = 0.0; /* 0 / 0, +-inf / +-inf -> 0 */ if (isnan(w)) w = 0.0; /* 0 / 0, +-inf / +-inf -> 0 */- double hi = fast_fmax4_double(x, y, z, w);+ VOLATILE double hi = fast_fmax4_double(x, y, z, w);+ FORCE_EVAL(hi); restore_fp_reg(oldreg); return hi; } -extern double rounded_hw_interval_div_double_down(double lo1, double hi1, double lo2, double hi2)+extern double rounded_hw_interval_div_double_down(VOLATILE double lo1, VOLATILE double hi1, VOLATILE double lo2, VOLATILE double hi2) { fp_reg oldreg = get_fp_reg(); set_rounding(oldreg, ROUND_DOWNWARD);- double x = (volatile double)(lo1 / lo2);- double y = (volatile double)(lo1 / hi2);- double z = (volatile double)(hi1 / lo2);- double w = (volatile double)(hi1 / hi2);+ MAY_MODIFY(lo1);+ MAY_MODIFY(hi1);+ MAY_MODIFY(lo2);+ MAY_MODIFY(hi2);+ double x = lo1 / lo2;+ double y = lo1 / hi2;+ double z = hi1 / lo2;+ double w = hi1 / hi2; if (isnan(x)) x = 0.0; /* 0 / 0, +-inf / +-inf -> 0 */ if (isnan(y)) y = 0.0; /* 0 / 0, +-inf / +-inf -> 0 */ if (isnan(z)) z = 0.0; /* 0 / 0, +-inf / +-inf -> 0 */ if (isnan(w)) w = 0.0; /* 0 / 0, +-inf / +-inf -> 0 */- double lo = fast_fmin4_double(x, y, z, w);+ VOLATILE double lo = fast_fmin4_double(x, y, z, w);+ FORCE_EVAL(lo); restore_fp_reg(oldreg); return lo; } -extern double rounded_hw_interval_div_add_double_up(double lo1, double hi1, double lo2, double hi2, double hi3)+extern double rounded_hw_interval_div_add_double_up(VOLATILE double lo1, VOLATILE double hi1, VOLATILE double lo2, VOLATILE double hi2, VOLATILE double hi3) { fp_reg oldreg = get_fp_reg(); set_rounding(oldreg, ROUND_UPWARD);- double x = (volatile double)(lo1 / lo2);- double y = (volatile double)(lo1 / hi2);- double z = (volatile double)(hi1 / lo2);- double w = (volatile double)(hi1 / hi2);+ MAY_MODIFY(lo1);+ MAY_MODIFY(hi1);+ MAY_MODIFY(lo2);+ MAY_MODIFY(hi2);+ MAY_MODIFY(hi3);+ double x = lo1 / lo2;+ double y = lo1 / hi2;+ double z = hi1 / lo2;+ double w = hi1 / hi2; if (isnan(x)) x = 0.0; /* 0 / 0, +-inf / +-inf -> 0 */ if (isnan(y)) y = 0.0; /* 0 / 0, +-inf / +-inf -> 0 */ if (isnan(z)) z = 0.0; /* 0 / 0, +-inf / +-inf -> 0 */ if (isnan(w)) w = 0.0; /* 0 / 0, +-inf / +-inf -> 0 */- volatile double hi = fast_fmax4_double(x, y, z, w) + hi3;+ VOLATILE double hi = fast_fmax4_double(x, y, z, w) + hi3;+ FORCE_EVAL(hi); restore_fp_reg(oldreg); return hi; } -extern double rounded_hw_interval_div_add_double_down(double lo1, double hi1, double lo2, double hi2, double lo3)+extern double rounded_hw_interval_div_add_double_down(VOLATILE double lo1, VOLATILE double hi1, VOLATILE double lo2, VOLATILE double hi2, VOLATILE double lo3) { fp_reg oldreg = get_fp_reg(); set_rounding(oldreg, ROUND_DOWNWARD);- double x = (volatile double)(lo1 / lo2);- double y = (volatile double)(lo1 / hi2);- double z = (volatile double)(hi1 / lo2);- double w = (volatile double)(hi1 / hi2);+ MAY_MODIFY(lo1);+ MAY_MODIFY(hi1);+ MAY_MODIFY(lo2);+ MAY_MODIFY(hi2);+ MAY_MODIFY(lo3);+ double x = lo1 / lo2;+ double y = lo1 / hi2;+ double z = hi1 / lo2;+ double w = hi1 / hi2; if (isnan(x)) x = 0.0; /* 0 / 0, +-inf / +-inf -> 0 */ if (isnan(y)) y = 0.0; /* 0 / 0, +-inf / +-inf -> 0 */ if (isnan(z)) z = 0.0; /* 0 / 0, +-inf / +-inf -> 0 */ if (isnan(w)) w = 0.0; /* 0 / 0, +-inf / +-inf -> 0 */- volatile double lo = fast_fmin4_double(x, y, z, w) + lo3;+ VOLATILE double lo = fast_fmin4_double(x, y, z, w) + lo3;+ FORCE_EVAL(lo); restore_fp_reg(oldreg); return lo; }@@ -333,10 +407,12 @@ native_rounding_mode nmode = hs_rounding_mode_to_native(mode); fp_reg oldreg = get_fp_reg(); set_rounding(oldreg, nmode);- volatile double s = 0.0;+ VOLATILE double s = 0.0;+ MAY_MODIFY(s); for (HsInt i = 0; i < length; ++i) { s += a[offset + i]; }+ FORCE_EVAL(s); restore_fp_reg(oldreg); return s; }@@ -401,9 +477,25 @@ native_rounding_mode nmode = hs_rounding_mode_to_native(mode); fp_reg oldreg = get_fp_reg(); set_rounding(oldreg, nmode);+#if defined(__GNUC__) && defined(USE_SSE2)+ // mingw-w64's sqrt may use x87 instead of SSE2+ HsInt i = 0;+ for (; i + 2 <= length; i += 2) {+ __m128d aa = _mm_loadu_pd(&a[offsetA + i]);+ __m128d bb = _mm_sqrt_pd(aa);+ _mm_storeu_pd(&result[offsetR + i], bb);+ }+ for (; i < length; ++i) {+ double aa = a[offsetA + i];+ double bb;+ __asm__ volatile("sqrtsd %1, %0" : "=x"(bb) : "x"(aa));+ result[offsetR + i] = bb;+ }+#else for (HsInt i = 0; i < length; ++i) { result[offsetR + i] = sqrt(a[offsetA + i]); }+#endif restore_fp_reg(oldreg); } @@ -411,11 +503,14 @@ // float // -static inline float rounded_add_impl_float(native_rounding_mode mode, float a, float b)+static inline float rounded_add_impl_float(native_rounding_mode mode, VOLATILE float a, VOLATILE float b) { fp_reg oldreg = get_fp_reg(); set_rounding(oldreg, mode);- volatile float c = a + b;+ MAY_MODIFY(a);+ MAY_MODIFY(b);+ VOLATILE float c = a + b;+ FORCE_EVAL(c); restore_fp_reg(oldreg); return c; }@@ -428,11 +523,14 @@ extern float rounded_hw_add_float_zero(float a, float b) { return rounded_add_impl_float(ROUND_TOWARDZERO, a, b); } -static inline float rounded_sub_impl_float(native_rounding_mode mode, float a, float b)+static inline float rounded_sub_impl_float(native_rounding_mode mode, VOLATILE float a, VOLATILE float b) { fp_reg oldreg = get_fp_reg(); set_rounding(oldreg, mode);- volatile float c = a - b;+ MAY_MODIFY(a);+ MAY_MODIFY(b);+ VOLATILE float c = a - b;+ FORCE_EVAL(c); restore_fp_reg(oldreg); return c; }@@ -445,11 +543,14 @@ extern float rounded_hw_sub_float_zero(float a, float b) { return rounded_sub_impl_float(ROUND_TOWARDZERO, a, b); } -static inline float rounded_mul_impl_float(native_rounding_mode mode, float a, float b)+static inline float rounded_mul_impl_float(native_rounding_mode mode, VOLATILE float a, VOLATILE float b) { fp_reg oldreg = get_fp_reg(); set_rounding(oldreg, mode);- volatile float c = a * b;+ MAY_MODIFY(a);+ MAY_MODIFY(b);+ VOLATILE float c = a * b;+ FORCE_EVAL(c); restore_fp_reg(oldreg); return c; }@@ -462,11 +563,14 @@ extern float rounded_hw_mul_float_zero(float a, float b) { return rounded_mul_impl_float(ROUND_TOWARDZERO, a, b); } -static inline float rounded_div_impl_float(native_rounding_mode mode, float a, float b)+static inline float rounded_div_impl_float(native_rounding_mode mode, VOLATILE float a, VOLATILE float b) { fp_reg oldreg = get_fp_reg(); set_rounding(oldreg, mode);- volatile float c = a / b;+ MAY_MODIFY(a);+ MAY_MODIFY(b);+ VOLATILE float c = a / b;+ FORCE_EVAL(c); restore_fp_reg(oldreg); return c; }@@ -479,11 +583,19 @@ extern float rounded_hw_div_float_zero(float a, float b) { return rounded_div_impl_float(ROUND_TOWARDZERO, a, b); } -static inline float rounded_sqrt_impl_float(native_rounding_mode mode, float a)+static inline float rounded_sqrt_impl_float(native_rounding_mode mode, VOLATILE float a) { fp_reg oldreg = get_fp_reg(); set_rounding(oldreg, mode);- volatile float c = sqrtf(a);+#if defined(__GNUC__) && defined(USE_SSE2)+ // mingw-w64's sqrtf may use x87 instead of SSE2+ float c;+ __asm__ volatile("sqrtss %1, %0" : "=x"(c) : "x"(a));+#else+ MAY_MODIFY(a);+ VOLATILE float c = sqrtf(a);+ FORCE_EVAL(c);+#endif restore_fp_reg(oldreg); return c; }@@ -496,11 +608,15 @@ extern float rounded_hw_sqrt_float_zero(float a) { return rounded_sqrt_impl_float(ROUND_TOWARDZERO, a); } -static inline float rounded_fma_impl_float(native_rounding_mode mode, float a, float b, float c)+static inline float rounded_fma_impl_float(native_rounding_mode mode, VOLATILE float a, VOLATILE float b, VOLATILE float c) { fp_reg oldreg = get_fp_reg(); set_rounding(oldreg, mode);- volatile float result = fmaf(a, b, c);+ MAY_MODIFY(a);+ MAY_MODIFY(b);+ MAY_MODIFY(c);+ VOLATILE float result = fmaf(a, b, c);+ FORCE_EVAL(result); restore_fp_reg(oldreg); return result; }@@ -513,15 +629,19 @@ extern float rounded_hw_fma_float_zero(float a, float b, float c) { return rounded_fma_impl_float(ROUND_TOWARDZERO, a, b, c); } -static inline float rounded_fma_if_fast_impl_float(native_rounding_mode mode, float a, float b, float c)+static inline float rounded_fma_if_fast_impl_float(native_rounding_mode mode, VOLATILE float a, VOLATILE float b, VOLATILE float c) { fp_reg oldreg = get_fp_reg(); set_rounding(oldreg, mode);+ MAY_MODIFY(a);+ MAY_MODIFY(b);+ MAY_MODIFY(c); #ifdef FP_FAST_FMAF- volatile float result = fmaf(a, b, c);+ VOLATILE float result = fmaf(a, b, c); #else- volatile float result = a * b + c;+ VOLATILE float result = a * b + c; #endif+ FORCE_EVAL(result); restore_fp_reg(oldreg); return result; }@@ -542,7 +662,8 @@ { fp_reg oldreg = get_fp_reg(); set_rounding(oldreg, mode);- volatile float result = (float)x;+ VOLATILE float result = (float)x;+ FORCE_EVAL(result); restore_fp_reg(oldreg); return result; }@@ -559,7 +680,8 @@ { fp_reg oldreg = get_fp_reg(); set_rounding(oldreg, mode);- volatile float result = (float)x;+ VOLATILE float result = (float)x;+ FORCE_EVAL(result); restore_fp_reg(oldreg); return result; }@@ -595,138 +717,182 @@ return fast_fmin_float(fast_fmin_float(x, y), fast_fmin_float(z, w)); } -extern float rounded_hw_interval_mul_float_up(float lo1, float hi1, float lo2, float hi2)+extern float rounded_hw_interval_mul_float_up(VOLATILE float lo1, VOLATILE float hi1, VOLATILE float lo2, VOLATILE float hi2) { fp_reg oldreg = get_fp_reg(); set_rounding(oldreg, ROUND_UPWARD);- float x = (volatile float)(lo1 * lo2);- float y = (volatile float)(lo1 * hi2);- float z = (volatile float)(hi1 * lo2);- float w = (volatile float)(hi1 * hi2);+ MAY_MODIFY(lo1);+ MAY_MODIFY(hi1);+ MAY_MODIFY(lo2);+ MAY_MODIFY(hi2);+ VOLATILE float x = lo1 * lo2;+ VOLATILE float y = lo1 * hi2;+ VOLATILE float z = hi1 * lo2;+ VOLATILE float w = hi1 * hi2; if (isnan(x)) x = 0.0f; /* 0 * inf -> 0 */ if (isnan(y)) y = 0.0f; /* 0 * inf -> 0 */ if (isnan(z)) z = 0.0f; /* 0 * inf -> 0 */ if (isnan(w)) w = 0.0f; /* 0 * inf -> 0 */- float hi = fast_fmax4_float(x, y, z, w);+ VOLATILE float hi = fast_fmax4_float(x, y, z, w);+ FORCE_EVAL(hi); restore_fp_reg(oldreg); return hi; } -extern float rounded_hw_interval_mul_float_down(float lo1, float hi1, float lo2, float hi2)+extern float rounded_hw_interval_mul_float_down(VOLATILE float lo1, VOLATILE float hi1, VOLATILE float lo2, VOLATILE float hi2) { fp_reg oldreg = get_fp_reg(); set_rounding(oldreg, ROUND_DOWNWARD);- float x = (volatile float)(lo1 * lo2);- float y = (volatile float)(lo1 * hi2);- float z = (volatile float)(hi1 * lo2);- float w = (volatile float)(hi1 * hi2);+ MAY_MODIFY(lo1);+ MAY_MODIFY(hi1);+ MAY_MODIFY(lo2);+ MAY_MODIFY(hi2);+ VOLATILE float x = lo1 * lo2;+ VOLATILE float y = lo1 * hi2;+ VOLATILE float z = hi1 * lo2;+ VOLATILE float w = hi1 * hi2; if (isnan(x)) x = 0.0f; /* 0 * inf -> 0 */ if (isnan(y)) y = 0.0f; /* 0 * inf -> 0 */ if (isnan(z)) z = 0.0f; /* 0 * inf -> 0 */ if (isnan(w)) w = 0.0f; /* 0 * inf -> 0 */- float lo = fast_fmin4_float(x, y, z, w);+ VOLATILE float lo = fast_fmin4_float(x, y, z, w);+ FORCE_EVAL(lo); restore_fp_reg(oldreg); return lo; } -extern float rounded_hw_interval_mul_add_float_up(float lo1, float hi1, float lo2, float hi2, float hi3)+extern float rounded_hw_interval_mul_add_float_up(VOLATILE float lo1, VOLATILE float hi1, VOLATILE float lo2, VOLATILE float hi2, VOLATILE float hi3) { fp_reg oldreg = get_fp_reg(); set_rounding(oldreg, ROUND_UPWARD);- float x = (volatile float)(lo1 * lo2);- float y = (volatile float)(lo1 * hi2);- float z = (volatile float)(hi1 * lo2);- float w = (volatile float)(hi1 * hi2);+ MAY_MODIFY(lo1);+ MAY_MODIFY(hi1);+ MAY_MODIFY(lo2);+ MAY_MODIFY(hi2);+ MAY_MODIFY(hi3);+ float x = lo1 * lo2;+ float y = lo1 * hi2;+ float z = hi1 * lo2;+ float w = hi1 * hi2; if (isnan(x)) x = 0.0f; /* 0 * inf -> 0 */ if (isnan(y)) y = 0.0f; /* 0 * inf -> 0 */ if (isnan(z)) z = 0.0f; /* 0 * inf -> 0 */ if (isnan(w)) w = 0.0f; /* 0 * inf -> 0 */- volatile float hi = fast_fmax4_float(x, y, z, w) + hi3;+ VOLATILE float hi = fast_fmax4_float(x, y, z, w) + hi3;+ FORCE_EVAL(hi); restore_fp_reg(oldreg); return hi; } -extern float rounded_hw_interval_mul_add_float_down(float lo1, float hi1, float lo2, float hi2, float lo3)+extern float rounded_hw_interval_mul_add_float_down(VOLATILE float lo1, VOLATILE float hi1, VOLATILE float lo2, VOLATILE float hi2, VOLATILE float lo3) { fp_reg oldreg = get_fp_reg(); set_rounding(oldreg, ROUND_DOWNWARD);- float x = (volatile float)(lo1 * lo2);- float y = (volatile float)(lo1 * hi2);- float z = (volatile float)(hi1 * lo2);- float w = (volatile float)(hi1 * hi2);+ MAY_MODIFY(lo1);+ MAY_MODIFY(hi1);+ MAY_MODIFY(lo2);+ MAY_MODIFY(hi2);+ MAY_MODIFY(lo3);+ float x = lo1 * lo2;+ float y = lo1 * hi2;+ float z = hi1 * lo2;+ float w = hi1 * hi2; if (isnan(x)) x = 0.0f; /* 0 * inf -> 0 */ if (isnan(y)) y = 0.0f; /* 0 * inf -> 0 */ if (isnan(z)) z = 0.0f; /* 0 * inf -> 0 */ if (isnan(w)) w = 0.0f; /* 0 * inf -> 0 */- volatile float lo = fast_fmin4_float(x, y, z, w) + lo3;+ VOLATILE float lo = fast_fmin4_float(x, y, z, w) + lo3;+ FORCE_EVAL(lo); restore_fp_reg(oldreg); return lo; } -extern float rounded_hw_interval_div_float_up(float lo1, float hi1, float lo2, float hi2)+extern float rounded_hw_interval_div_float_up(VOLATILE float lo1, VOLATILE float hi1, VOLATILE float lo2, VOLATILE float hi2) { fp_reg oldreg = get_fp_reg(); set_rounding(oldreg, ROUND_UPWARD);- float x = (volatile float)(lo1 / lo2);- float y = (volatile float)(lo1 / hi2);- float z = (volatile float)(hi1 / lo2);- float w = (volatile float)(hi1 / hi2);+ MAY_MODIFY(lo1);+ MAY_MODIFY(hi1);+ MAY_MODIFY(lo2);+ MAY_MODIFY(hi2);+ float x = lo1 / lo2;+ float y = lo1 / hi2;+ float z = hi1 / lo2;+ float w = hi1 / hi2; if (isnan(x)) x = 0.0f; /* 0 / 0, +-inf / +-inf -> 0 */ if (isnan(y)) y = 0.0f; /* 0 / 0, +-inf / +-inf -> 0 */ if (isnan(z)) z = 0.0f; /* 0 / 0, +-inf / +-inf -> 0 */ if (isnan(w)) w = 0.0f; /* 0 / 0, +-inf / +-inf -> 0 */- float hi = fast_fmax4_float(x, y, z, w);+ VOLATILE float hi = fast_fmax4_float(x, y, z, w);+ FORCE_EVAL(hi); restore_fp_reg(oldreg); return hi; } -extern float rounded_hw_interval_div_float_down(float lo1, float hi1, float lo2, float hi2)+extern float rounded_hw_interval_div_float_down(VOLATILE float lo1, VOLATILE float hi1, VOLATILE float lo2, VOLATILE float hi2) { fp_reg oldreg = get_fp_reg(); set_rounding(oldreg, ROUND_DOWNWARD);- float x = (volatile float)(lo1 / lo2);- float y = (volatile float)(lo1 / hi2);- float z = (volatile float)(hi1 / lo2);- float w = (volatile float)(hi1 / hi2);+ MAY_MODIFY(lo1);+ MAY_MODIFY(hi1);+ MAY_MODIFY(lo2);+ MAY_MODIFY(hi2);+ float x = lo1 / lo2;+ float y = lo1 / hi2;+ float z = hi1 / lo2;+ float w = hi1 / hi2; if (isnan(x)) x = 0.0f; /* 0 / 0, +-inf / +-inf -> 0 */ if (isnan(y)) y = 0.0f; /* 0 / 0, +-inf / +-inf -> 0 */ if (isnan(z)) z = 0.0f; /* 0 / 0, +-inf / +-inf -> 0 */ if (isnan(w)) w = 0.0f; /* 0 / 0, +-inf / +-inf -> 0 */- float lo = fast_fmin4_float(x, y, z, w);+ VOLATILE float lo = fast_fmin4_float(x, y, z, w);+ FORCE_EVAL(lo); restore_fp_reg(oldreg); return lo; } -extern float rounded_hw_interval_div_add_float_up(float lo1, float hi1, float lo2, float hi2, float hi3)+extern float rounded_hw_interval_div_add_float_up(VOLATILE float lo1, VOLATILE float hi1, VOLATILE float lo2, VOLATILE float hi2, VOLATILE float hi3) { fp_reg oldreg = get_fp_reg(); set_rounding(oldreg, ROUND_UPWARD);- float x = (volatile float)(lo1 / lo2);- float y = (volatile float)(lo1 / hi2);- float z = (volatile float)(hi1 / lo2);- float w = (volatile float)(hi1 / hi2);+ MAY_MODIFY(lo1);+ MAY_MODIFY(hi1);+ MAY_MODIFY(lo2);+ MAY_MODIFY(hi2);+ MAY_MODIFY(hi3);+ float x = lo1 / lo2;+ float y = lo1 / hi2;+ float z = hi1 / lo2;+ float w = hi1 / hi2; if (isnan(x)) x = 0.0f; /* 0 / 0, +-inf / +-inf -> 0 */ if (isnan(y)) y = 0.0f; /* 0 / 0, +-inf / +-inf -> 0 */ if (isnan(z)) z = 0.0f; /* 0 / 0, +-inf / +-inf -> 0 */ if (isnan(w)) w = 0.0f; /* 0 / 0, +-inf / +-inf -> 0 */- volatile float hi = fast_fmax4_float(x, y, z, w) + hi3;+ VOLATILE float hi = fast_fmax4_float(x, y, z, w) + hi3;+ FORCE_EVAL(hi); restore_fp_reg(oldreg); return hi; } -extern float rounded_hw_interval_div_add_float_down(float lo1, float hi1, float lo2, float hi2, float lo3)+extern float rounded_hw_interval_div_add_float_down(VOLATILE float lo1, VOLATILE float hi1, VOLATILE float lo2, VOLATILE float hi2, VOLATILE float lo3) { fp_reg oldreg = get_fp_reg(); set_rounding(oldreg, ROUND_DOWNWARD);- float x = (volatile float)(lo1 / lo2);- float y = (volatile float)(lo1 / hi2);- float z = (volatile float)(hi1 / lo2);- float w = (volatile float)(hi1 / hi2);+ MAY_MODIFY(lo1);+ MAY_MODIFY(hi1);+ MAY_MODIFY(lo2);+ MAY_MODIFY(hi2);+ MAY_MODIFY(lo3);+ float x = lo1 / lo2;+ float y = lo1 / hi2;+ float z = hi1 / lo2;+ float w = hi1 / hi2; if (isnan(x)) x = 0.0f; /* 0 / 0, +-inf / +-inf -> 0 */ if (isnan(y)) y = 0.0f; /* 0 / 0, +-inf / +-inf -> 0 */ if (isnan(z)) z = 0.0f; /* 0 / 0, +-inf / +-inf -> 0 */ if (isnan(w)) w = 0.0f; /* 0 / 0, +-inf / +-inf -> 0 */- volatile float lo = fast_fmin4_float(x, y, z, w) + lo3;+ VOLATILE float lo = fast_fmin4_float(x, y, z, w) + lo3;+ FORCE_EVAL(lo); restore_fp_reg(oldreg); return lo; }@@ -740,10 +906,12 @@ native_rounding_mode nmode = hs_rounding_mode_to_native(mode); fp_reg oldreg = get_fp_reg(); set_rounding(oldreg, nmode);- volatile float s = 0.0f;+ VOLATILE float s = 0.0f;+ MAY_MODIFY(s); for (HsInt i = 0; i < length; ++i) { s += a[offset + i]; }+ FORCE_EVAL(s); restore_fp_reg(oldreg); return s; }@@ -808,8 +976,24 @@ native_rounding_mode nmode = hs_rounding_mode_to_native(mode); fp_reg oldreg = get_fp_reg(); set_rounding(oldreg, nmode);+#if defined(__GNUC__) && defined(USE_SSE2)+ // mingw-w64's sqrtf may use x87 instead of SSE2+ HsInt i = 0;+ for (; i + 4 <= length; i += 4) {+ __m128 aa = _mm_loadu_ps(&a[offsetA + i]);+ __m128 bb = _mm_sqrt_ps(aa);+ _mm_storeu_ps(&result[offsetR + i], bb);+ }+ for (; i < length; ++i) {+ float aa = a[offsetA + i];+ float bb;+ __asm__ volatile("sqrtss %1, %0" : "=x"(bb) : "x"(aa));+ result[offsetR + i] = bb;+ }+#else for (HsInt i = 0; i < length; ++i) { result[offsetR + i] = sqrtf(a[offsetA + i]); }+#endif restore_fp_reg(oldreg); }
cbits/rounded.c view
@@ -43,6 +43,20 @@ #error "Invalid configuration detected: USE_SSE2 and USE_AVX512 are mutually exclusive" #endif +#if defined(__GNUC__) && defined(__SSE2__)+#define MAY_MODIFY(x) __asm__ volatile("" : "+x"(x))+#define FORCE_EVAL(x) __asm__ volatile("" : : "x"(x))+#define VOLATILE+#elif defined(__GNUC__) && defined(__aarch64__)+#define MAY_MODIFY(x) __asm__ volatile("" : "+w"(x))+#define FORCE_EVAL(x) __asm__ volatile("" : : "w"(x))+#define VOLATILE+#else+#define MAY_MODIFY(x) (void)(x)+#define FORCE_EVAL(x) (void)(x)+#define VOLATILE volatile+#endif+ #if defined(USE_AVX512) #include <x86intrin.h>
rounded-hw.cabal view
@@ -1,7 +1,8 @@-cabal-version: 2.2+cabal-version: 3.0+-- asm-sources is Cabal 3.0 feature name: rounded-hw-version: 0.3.0+version: 0.4.0 synopsis: Directed rounding for built-in floating types description: Please see the README on GitHub at <https://github.com/minoki/haskell-floating-point/tree/master/rounded-hw#readme> category: Numeric, Math@@ -9,18 +10,19 @@ bug-reports: https://github.com/minoki/haskell-floating-point/issues author: ARATA Mizuki maintainer: minorinoki@gmail.com-copyright: 2020 ARATA Mizuki+copyright: 2020-2023 ARATA Mizuki license: BSD-3-Clause license-file: LICENSE build-type: Custom tested-with:- GHC == 8.6.5, GHC == 8.8.4, GHC == 8.10.7, GHC == 9.0.2, GHC == 9.2.1+ GHC == 8.6.5, GHC == 8.8.4, GHC == 8.10.7, GHC == 9.0.2, GHC == 9.2.8, GHC == 9.4.8, GHC == 9.6.3, GHC == 9.8.1 extra-source-files: README.md ChangeLog.md cbits/rounded-common.inl cbits/rounded-avx512.inl cbits/interval-prim-x86_64-sse2.S+ cbits/interval-prim-x86_64-avx512.S source-repository head type: git@@ -30,8 +32,8 @@ -- Custom setup is required to allow assembly sources to #include "ghcconfig.h" custom-setup setup-depends:- Cabal >=1.24- , base >=4.7+ Cabal >=3.0 && <3.11+ , base >=4.12 && <4.20 flag pure-hs description: Disable FFI@@ -56,7 +58,7 @@ flag x87-long-double description: Support x87 "long double" manual: True- default: True+ default: False flag float128 description: Support Float128@@ -65,12 +67,12 @@ common deps build-depends:- array- , base >=4.12 && <5- , deepseq+ array >=0.5.2.0 && <0.6+ , base >=4.12 && <4.20+ , deepseq >=1.4.4.0 && <1.6 , fp-ieee ==0.1.*- , primitive- , vector+ , primitive >=0.6.1.1 && <0.10+ , vector >=0.12.0.1 && <0.14 common options ghc-options: -Wcompat@@ -101,7 +103,7 @@ hs-source-dirs: src build-depends:- tagged+ tagged >=0.8.6 && <0.9 ghc-options: -Wall -- Use FFI when flag(pure-hs) is off if !flag(pure-hs)@@ -124,12 +126,20 @@ exposed-modules: Numeric.Rounded.Hardware.Backend.FastFFI cpp-options: -DUSE_GHC_PRIM- if flag(avx512)- c-sources:- cbits/interval-prim-x86_64-avx512.S+ if impl(ghc >= 9.2)+ if flag(avx512)+ asm-sources:+ cbits/interval-prim-x86_64-avx512.S+ else+ asm-sources:+ cbits/interval-prim-x86_64.S else- c-sources:- cbits/interval-prim-x86_64.S+ if flag(avx512)+ c-sources:+ cbits/interval-prim-x86_64-avx512.S+ else+ c-sources:+ cbits/interval-prim-x86_64.S -- flag(x87-long-double): Support LongDouble on x86 if flag(x87-long-double) && (arch(i386) || arch(x86_64)) other-modules:@@ -138,7 +148,7 @@ c-sources: cbits/rounded-x87longdouble.c build-depends:- long-double+ long-double >=0.1 && <0.2 -- flag(float128): Support Float128 if flag(float128) other-modules:@@ -147,7 +157,7 @@ c-sources: cbits/rounded-float128.c build-depends:- float128+ float128 >=0.1 && <0.2 default-language: Haskell2010 test-suite rounded-hw-doctests@@ -155,8 +165,10 @@ type: exitcode-stdio-1.0 main-is: doctests.hs build-depends:- doctest >=0.8+ doctest ^>=0.22.2 default-language: Haskell2010+ if impl(ghc >= 9.8)+ buildable: False test-suite rounded-hw-test import: deps, options@@ -175,9 +187,9 @@ test ghc-options: -threaded -rtsopts -with-rtsopts=-N build-depends:- QuickCheck- , hspec- , random+ QuickCheck ^>=2.14.3+ , hspec ^>=2.11.7+ , random ^>=1.2.1.1 , rounded-hw if flag(x87-long-double) && (arch(i386) || (arch(x86_64) && !os(windows))) -- Support for 80-bit long double is not good on Win64, so don't test@@ -205,7 +217,5 @@ benchmark build-depends: rounded-hw- , tasty-bench- mixins:- tasty-bench (Test.Tasty.Bench as Gauge, Test.Tasty.Bench as Gauge.Main, Test.Tasty.Bench as Gauge.Benchmark)+ , tasty-bench ^>=0.3.5 default-language: Haskell2010
src/Numeric/Rounded/Hardware/Backend/C.hs view
@@ -90,11 +90,6 @@ (F.roundedFromWord64 r x) {-# INLINE roundedFloatFromWord64 #-} -intervalFloatFromInteger :: Integer -> (Rounded 'TowardNegInf Float, Rounded 'TowardInf Float)-intervalFloatFromInteger x- | -0x1000000 <= x && x <= 0x1000000 {- abs x <= 2^24 -} = (Rounded (fromInteger x), Rounded (fromInteger x))- | otherwise = intervalFromInteger_default x- roundedFloatFromRealFloat :: RealFloat a => RoundingMode -> a -> Float roundedFloatFromRealFloat r x | isNaN x = 0/0 | isInfinite x = if x > 0 then 1/0 else -1/0
src/Numeric/Rounded/Hardware/Backend/Default.hs view
@@ -89,6 +89,7 @@ map_roundedSqrt mode vec = unsafeCoerce (map_roundedSqrt mode (unsafeCoerce vec) :: VS.Vector DoubleImpl) {-# INLINE map_roundedSqrt #-} +#if !MIN_VERSION_base(4, 16, 0) -- orphaned rules {-# RULES "fromIntegral/a->Rounded ToNearest Float"@@ -108,3 +109,4 @@ "fromIntegral/a->Rounded TowardZero Double" fromIntegral = \x -> (Rounded (fromIntegralTowardZero x) :: Rounded 'TowardZero Double) #-}+#endif
src/Numeric/Rounded/Hardware/Backend/FastFFI.hs view
@@ -50,8 +50,8 @@ import Foreign.Storable (Storable) import GHC.Exts import GHC.Generics (Generic)-import GHC.Int (Int64 (I64#))-import GHC.Word (Word64 (W64#))+-- import GHC.Int (Int64 (I64#))+-- import GHC.Word (Word64 (W64#)) import qualified Numeric.Rounded.Hardware.Backend.C as C import Numeric.Rounded.Hardware.Internal.Class import System.IO.Unsafe (unsafePerformIO)@@ -166,6 +166,7 @@ , Double# -- upper, %xmm2 #) +{- foreign import prim "rounded_hw_interval_sqrt" fastIntervalSqrt# :: Double# -- lower 1, %xmm1 -> Double# -- upper 1, %xmm2@@ -173,7 +174,7 @@ , Double# -- upper, %xmm2 #) -#if WORD_SIZE_IN_BITS >= 64+#if WORD_SIZE_IN_BITS >= 64 && !MIN_VERSION_base(4,17,0) type INT64# = Int# type WORD64# = Word# #else@@ -187,7 +188,6 @@ , Double# -- upper, %xmm2 #) -{- foreign import prim "rounded_hw_interval_from_word64" fastIntervalFromWord64# :: WORD64# -- value -> (# Double# -- lower, %xmm1@@ -210,6 +210,7 @@ (# l2, h2 #) -> (D# l2, D# h2) {-# INLINE fastIntervalRecip #-} +{- fastIntervalSqrt :: Double -> Double -> (Double, Double) fastIntervalSqrt (D# l1) (D# h1) = case fastIntervalSqrt# l1 h1 of (# l2, h2 #) -> (D# l2, D# h2)@@ -220,7 +221,6 @@ (# l, h #) -> (D# l, D# h) {-# INLINE fastIntervalFromInt64 #-} -{- fastIntervalFromWord64 :: Word64 -> (Double, Double) fastIntervalFromWord64 (W64# x) = case fastIntervalFromWord64# x of (# l, h #) -> (D# l, D# h)
src/Numeric/Rounded/Hardware/Backend/ViaRational.hs view
@@ -54,7 +54,7 @@ TowardInf -> 0 TowardZero -> 0 roundedMul r (ViaRational x) (ViaRational y)- | isNaN x || isNaN y || isInfinite x || isInfinite y || isNegativeZero x || isNegativeZero y = ViaRational (x * y)+ | isNaN x || isNaN y || isInfinite x || isInfinite y || x == 0 || y == 0 = ViaRational (x * y) | otherwise = roundedFromRational r (toRational x * toRational y) roundedFusedMultiplyAdd r (ViaRational x) (ViaRational y) (ViaRational z) | isFinite x && isFinite y && isFinite z = case toRational x * toRational y + toRational z of
src/Numeric/Rounded/Hardware/Interval.hs view
@@ -1,3 +1,4 @@+{-# LANGUAGE CPP #-} {-# LANGUAGE DataKinds #-} {-# LANGUAGE DeriveGeneric #-} {-# LANGUAGE MultiParamTypeClasses #-}@@ -316,9 +317,11 @@ in pairToInterval (x, y) -- unsafeReplace, unsafeAccum, unsafeAccumArray: Use default +#if !MIN_VERSION_base(4, 16, 0) {-# RULES "fromIntegral/a->Interval Float" fromIntegral = \x -> case intervalFromIntegral x of (l, u) -> I l u :: Interval Float "fromIntegral/a->Interval Double" fromIntegral = \x -> case intervalFromIntegral x of (l, u) -> I l u :: Interval Double #-}+#endif
src/Numeric/Rounded/Hardware/Interval/NonEmpty.hs view
@@ -1,3 +1,4 @@+{-# LANGUAGE CPP #-} {-# LANGUAGE DataKinds #-} {-# LANGUAGE DeriveGeneric #-} {-# LANGUAGE MultiParamTypeClasses #-}@@ -355,9 +356,11 @@ in pairToInterval (x, y) -- unsafeReplace, unsafeAccum, unsafeAccumArray: Use default +#if !MIN_VERSION_base(4, 16, 0) {-# RULES "fromIntegral/a->Interval Float" fromIntegral = \x -> case intervalFromIntegral x of (l, u) -> I l u :: Interval Float "fromIntegral/a->Interval Double" fromIntegral = \x -> case intervalFromIntegral x of (l, u) -> I l u :: Interval Double #-}+#endif