vendor/compiler_builtins/libm/src/math/remquo.rs

   1 #[cfg_attr(all(test, assert_no_panic), no_panic::no_panic)]
   2 pub fn remquo(mut x: f64, mut y: f64) -> (f64, i32) {
   3     let ux: u64 = x.to_bits();
   4     let mut uy: u64 = y.to_bits();
   5     let mut ex = ((ux >> 52) & 0x7ff) as i32;
   6     let mut ey = ((uy >> 52) & 0x7ff) as i32;
   7     let sx = (ux >> 63) != 0;
   8     let sy = (uy >> 63) != 0;
   9     let mut q: u32;
  10     let mut i: u64;
  11     let mut uxi: u64 = ux;
  12
  13     if (uy << 1) == 0 || y.is_nan() || ex == 0x7ff {
  14         return ((x * y) / (x * y), 0);
  15     }
  16     if (ux << 1) == 0 {
  17         return (x, 0);
  18     }
  19
  20     /* normalize x and y */
  21     if ex == 0 {
  22         i = uxi << 12;
  23         while (i >> 63) == 0 {
  24             ex -= 1;
  25             i <<= 1;
  26         }
  27         uxi <<= -ex + 1;
  28     } else {
  29         uxi &= (!0) >> 12;
  30         uxi |= 1 << 52;
  31     }
  32     if ey == 0 {
  33         i = uy << 12;
  34         while (i >> 63) == 0 {
  35             ey -= 1;
  36             i <<= 1;
  37         }
  38         uy <<= -ey + 1;
  39     } else {
  40         uy &= (!0) >> 12;
  41         uy |= 1 << 52;
  42     }
  43
  44     q = 0;
  45
  46     if ex + 1 != ey {
  47         if ex < ey {
  48             return (x, 0);
  49         }
  50         /* x mod y */
  51         while ex > ey {
  52             i = uxi.wrapping_sub(uy);
  53             if (i >> 63) == 0 {
  54                 uxi = i;
  55                 q += 1;
  56             }
  57             uxi <<= 1;
  58             q <<= 1;
  59             ex -= 1;
  60         }
  61         i = uxi.wrapping_sub(uy);
  62         if (i >> 63) == 0 {
  63             uxi = i;
  64             q += 1;
  65         }
  66         if uxi == 0 {
  67             ex = -60;
  68         } else {
  69             while (uxi >> 52) == 0 {
  70                 uxi <<= 1;
  71                 ex -= 1;
  72             }
  73         }
  74     }
  75
  76     /* scale result and decide between |x| and |x|-|y| */
  77     if ex > 0 {
  78         uxi -= 1 << 52;
  79         uxi |= (ex as u64) << 52;
  80     } else {
  81         uxi >>= -ex + 1;
  82     }
  83     x = f64::from_bits(uxi);
  84     if sy {
  85         y = -y;
  86     }
  87     if ex == ey || (ex + 1 == ey && (2.0 * x > y || (2.0 * x == y && (q % 2) != 0))) {
  88         x -= y;
  89         // TODO: this matches musl behavior, but it is incorrect
  90         q = q.wrapping_add(1);
  91     }
  92     q &= 0x7fffffff;
  93     let quo = if sx ^ sy { -(q as i32) } else { q as i32 };
  94     if sx {
  95         (-x, quo)
  96     } else {
  97         (x, quo)
  98     }
  99 }
 100
 101 #[cfg(test)]
 102 mod tests {
 103     use super::remquo;
 104
 105     #[test]
 106     fn test_q_overflow() {
 107         // 0xc000000000000001, 0x04c0000000000004
 108         let _ = remquo(-2.0000000000000004, 8.406091369059082e-286);
 109     }
 110 }