More work toward generic

Author: tgross35

src/math/fma.rs

@@ -1,8 +1,6 @@
-/// Floating multiply add (f64)
+/// Fused multiply add (f64)
 ///
-/// Computes `(x*y)+z`, rounded as one ternary operation:
-/// Computes the value (as if) to infinite precision and rounds once to the result format,
-/// according to the rounding mode characterized by the value of FLT_ROUNDS.
+/// Computes `(x*y)+z`, rounded as one ternary operation (i.e. calculated with infinite precision).
 #[cfg_attr(all(test, assert_no_panic), no_panic::no_panic)]
 pub fn fma(x: f64, y: f64, z: f64) -> f64 {
     return super::generic::fma(x, y, z);

src/math/generic/fma.rs

@@ -1,17 +1,13 @@
 use core::{f32, f64};
 
 use super::super::support::{DInt, HInt, IntTy};
-use super::super::{Float, Int, MinInt, scalbn};
+use super::super::{CastFrom, Float, Int, MinInt};
 
 const ZEROINFNAN: i32 = 0x7ff - 0x3ff - 52 - 1;
 
 type F = f64;
 
-/// Floating multiply add (f64)
-///
-/// Computes `(x*y)+z`, rounded as one ternary operation:
-/// Computes the value (as if) to infinite precision and rounds once to the result format,
-/// according to the rounding mode characterized by the value of FLT_ROUNDS.
+/// Fused multiply-add.
 #[cfg_attr(all(test, assert_no_panic), no_panic::no_panic)]
 pub fn fma(x: f64, y: f64, z: f64) -> f64 {
     // let x1p63: f64 = f64::from_bits(0x43e0000000000000); // 0x1p63 === 2 ^ 63
@@ -125,7 +121,8 @@ pub fn fma(x: f64, y: f64, z: f64) -> f64 {
     if neg {
         i = -i;
     }
-    let mut r: f64 = i as f64; /* |r| is in [0x1p62,0x1p63] */
+
+    let mut r: f64 = f64::cast_from_lossy(i); /* |r| is in [0x1p62,0x1p63] */
 
     if e < -(F::EXP_BIAS as i32 - 1) - (sbits - 2) {
         /* result is subnormal before rounding */
@@ -143,11 +140,11 @@ pub fn fma(x: f64, y: f64, z: f64) -> f64 {
             /* one bit is lost when scaled, add another top bit to
              * only round once at conversion if it is inexact */
             if (rhi << F::SIG_BITS) != 0 {
-                i = ((rhi >> 1) | (rhi & 1) | (1 << 62)) as i64;
+                i = ((rhi >> 1) | (rhi & 1) | (1 << 62)).signed();
                 if neg {
                     i = -i;
                 }
-                r = i as f64;
+                r = F::cast_from(i);
                 r = 2.0 * r - c; /* remove top bit */
 
                 /* raise underflow portably, such that it
@@ -161,10 +158,11 @@ pub fn fma(x: f64, y: f64, z: f64) -> f64 {
             if neg {
                 i = -i;
             }
-            r = i as f64;
+            r = f64::cast_from(i);
         }
     }
-    scalbn(r, e)
+
+    F::scalbn(r, e)
 }
 
 struct Norm<F: Float> {
@@ -202,6 +200,7 @@ impl<F: Float> Norm<F> {
 trait RaiseUnderflow {
     fn raise_underflow(self) -> Self;
     fn raise_underflow2(self) -> Self;
+    fn scalbn(self, n: i32) -> Self;
 }
 
 impl RaiseUnderflow for f64 {
@@ -217,4 +216,8 @@ impl RaiseUnderflow for f64 {
         let tiny: f64 = f64::MIN_POSITIVE / f32::MIN_POSITIVE as f64 * self;
         (tiny * tiny) * (self - self)
     }
+
+    fn scalbn(self, n: i32) -> Self {
+        super::super::scalbn(self, n)
+    }
 }

src/math/support/int_traits.rs

@@ -400,6 +400,30 @@ macro_rules! cast_into {
     )*};
 }
 
+macro_rules! cast_into_float {
+    ($ty:ty) => {
+        #[cfg(f16_enabled)]
+        cast_into_float!($ty; f16);
+
+        cast_into_float!($ty; f32, f64);
+
+        #[cfg(f128_enabled)]
+        cast_into_float!($ty; f128);
+    };
+    ($ty:ty; $($into:ty),*) => {$(
+        impl CastInto<$into> for $ty {
+            fn cast(self) -> $into {
+                debug_assert_eq!(self as $into as $ty, self, "inexact float cast");
+                self as $into
+            }
+
+            fn cast_lossy(self) -> $into {
+                self as $into
+            }
+        }
+    )*};
+}
+
 cast_into!(usize);
 cast_into!(isize);
 cast_into!(u8);
@@ -412,3 +436,8 @@ cast_into!(u64);
 cast_into!(i64);
 cast_into!(u128);
 cast_into!(i128);
+
+cast_into_float!(i16);
+cast_into_float!(i32);
+cast_into_float!(i64);
+cast_into_float!(i128);