use core::marker::PhantomData;
use core::mem;
use core::ops::Neg;
-use smallvec::{SmallVec, smallvec};
+use smallvec::{smallvec, SmallVec};
#[must_use]
pub struct IeeeFloat<S> {
#[derive(Copy, Clone, PartialEq, Eq, Debug)]
enum Loss {
// Example of truncated bits:
- ExactlyZero, // 000000
+ ExactlyZero, // 000000
LessThanHalf, // 0xxxxx x's not all zero
- ExactlyHalf, // 100000
+ ExactlyHalf, // 100000
MoreThanHalf, // 1xxxxx x's not all zero
}
impl<S: Semantics> PartialOrd for IeeeFloat<S> {
fn partial_cmp(&self, rhs: &Self) -> Option<Ordering> {
match (self.category, rhs.category) {
- (Category::NaN, _) |
- (_, Category::NaN) => None,
+ (Category::NaN, _) | (_, Category::NaN) => None,
(Category::Infinity, Category::Infinity) => Some((!self.sign).cmp(&(!rhs.sign))),
(Category::Zero, Category::Zero) => Some(Ordering::Equal),
- (Category::Infinity, _) |
- (Category::Normal, Category::Zero) => Some((!self.sign).cmp(&self.sign)),
+ (Category::Infinity, _) | (Category::Normal, Category::Zero) => {
+ Some((!self.sign).cmp(&self.sign))
+ }
- (_, Category::Infinity) |
- (Category::Zero, Category::Normal) => Some(rhs.sign.cmp(&(!rhs.sign))),
+ (_, Category::Infinity) | (Category::Zero, Category::Normal) => {
+ Some(rhs.sign.cmp(&(!rhs.sign)))
+ }
(Category::Normal, Category::Normal) => {
// Two normal numbers. Do they have the same sign?
impl<S: Semantics> fmt::Debug for IeeeFloat<S> {
fn fmt(&self, f: &mut fmt::Formatter<'_>) -> fmt::Result {
- write!(f, "{}({:?} | {}{:?} * 2^{})",
- self, self.category,
- if self.sign { "-" } else { "+" },
- self.sig,
- self.exp)
+ write!(
+ f,
+ "{}({:?} | {}{:?} * 2^{})",
+ self,
+ self.category,
+ if self.sign { "-" } else { "+" },
+ self.sig,
+ self.exp
+ )
}
}
fn qnan(payload: Option<u128>) -> Self {
IeeeFloat {
- sig: [
- S::QNAN_SIGNIFICAND |
- payload.map_or(0, |payload| {
- // Zero out the excess bits of the significand.
- payload & ((1 << S::QNAN_BIT) - 1)
- }),
- ],
+ sig: [S::QNAN_SIGNIFICAND
+ | payload.map_or(0, |payload| {
+ // Zero out the excess bits of the significand.
+ payload & ((1 << S::QNAN_BIT) - 1)
+ })],
exp: S::MAX_EXP + 1,
category: Category::NaN,
sign: false,
}
// Sign may depend on rounding mode; handled below.
- (_, Category::Zero) |
- (Category::NaN, _) |
- (Category::Infinity, Category::Normal) => Status::OK,
+ (_, Category::Zero) | (Category::NaN, _) | (Category::Infinity, Category::Normal) => {
+ Status::OK
+ }
- (Category::Zero, _) |
- (_, Category::NaN) |
- (_, Category::Infinity) => {
+ (Category::Zero, _) | (_, Category::NaN) | (_, Category::Infinity) => {
self = rhs;
Status::OK
}
// If two numbers add (exactly) to zero, IEEE 754 decrees it is a
// positive zero unless rounding to minus infinity, except that
// adding two like-signed zeroes gives that zero.
- if self.category == Category::Zero &&
- (rhs.category != Category::Zero || self.sign != rhs.sign)
+ if self.category == Category::Zero
+ && (rhs.category != Category::Zero || self.sign != rhs.sign)
{
self.sign = round == Round::TowardNegative;
}
Status::OK.and(self)
}
- (Category::Zero, Category::Infinity) |
- (Category::Infinity, Category::Zero) => Status::INVALID_OP.and(Self::NAN),
+ (Category::Zero, Category::Infinity) | (Category::Infinity, Category::Zero) => {
+ Status::INVALID_OP.and(Self::NAN)
+ }
- (_, Category::Infinity) |
- (Category::Infinity, _) => {
+ (_, Category::Infinity) | (Category::Infinity, _) => {
self.category = Category::Infinity;
Status::OK.and(self)
}
- (Category::Zero, _) |
- (_, Category::Zero) => {
+ (Category::Zero, _) | (_, Category::Zero) => {
self.category = Category::Zero;
Status::OK.and(self)
}
(Category::Normal, Category::Normal) => {
self.exp += rhs.exp;
let mut wide_sig = [0; 2];
- let loss = sig::mul(
- &mut wide_sig,
- &mut self.exp,
- &self.sig,
- &rhs.sig,
- S::PRECISION,
- );
+ let loss =
+ sig::mul(&mut wide_sig, &mut self.exp, &self.sig, &rhs.sig, S::PRECISION);
self.sig = [wide_sig[0]];
let mut status;
self = unpack!(status=, self.normalize(round, loss));
// Extend the addend significand to ext_precision - 1. This guarantees
// that the high bit of the significand is zero (same as wide_sig),
// so the addition will overflow (if it does overflow at all) into the top bit.
- sig::shift_left(
- &mut ext_addend_sig,
- &mut 0,
- ext_precision - 1 - S::PRECISION,
- );
+ sig::shift_left(&mut ext_addend_sig, &mut 0, ext_precision - 1 - S::PRECISION);
loss = sig::add_or_sub(
&mut wide_sig,
&mut self.exp,
// If two numbers add (exactly) to zero, IEEE 754 decrees it is a
// positive zero unless rounding to minus infinity, except that
// adding two like-signed zeroes gives that zero.
- if self.category == Category::Zero && !status.intersects(Status::UNDERFLOW) &&
- self.sign != addend.sign
+ if self.category == Category::Zero
+ && !status.intersects(Status::UNDERFLOW)
+ && self.sign != addend.sign
{
self.sign = round == Round::TowardNegative;
}
Status::OK.and(self)
}
- (Category::Infinity, Category::Infinity) |
- (Category::Zero, Category::Zero) => Status::INVALID_OP.and(Self::NAN),
+ (Category::Infinity, Category::Infinity) | (Category::Zero, Category::Zero) => {
+ Status::INVALID_OP.and(Self::NAN)
+ }
- (Category::Infinity, _) |
- (Category::Zero, _) => Status::OK.and(self),
+ (Category::Infinity, _) | (Category::Zero, _) => Status::OK.and(self),
(Category::Normal, Category::Infinity) => {
self.category = Category::Zero;
fn c_fmod(mut self, rhs: Self) -> StatusAnd<Self> {
match (self.category, rhs.category) {
- (Category::NaN, _) |
- (Category::Zero, Category::Infinity) |
- (Category::Zero, Category::Normal) |
- (Category::Normal, Category::Infinity) => Status::OK.and(self),
+ (Category::NaN, _)
+ | (Category::Zero, Category::Infinity)
+ | (Category::Zero, Category::Normal)
+ | (Category::Normal, Category::Infinity) => Status::OK.and(self),
(_, Category::NaN) => {
self.sign = false;
Status::OK.and(self)
}
- (Category::Infinity, _) |
- (_, Category::Zero) => Status::INVALID_OP.and(Self::NAN),
+ (Category::Infinity, _) | (_, Category::Zero) => Status::INVALID_OP.and(Self::NAN),
(Category::Normal, Category::Normal) => {
- while self.is_finite_non_zero() && rhs.is_finite_non_zero() &&
- self.cmp_abs_normal(rhs) != Ordering::Less
+ while self.is_finite_non_zero()
+ && rhs.is_finite_non_zero()
+ && self.cmp_abs_normal(rhs) != Ordering::Less
{
let mut v = rhs.scalbn(self.ilogb() - rhs.ilogb());
if self.cmp_abs_normal(v) == Ordering::Less {
// 1. exponent != S::MIN_EXP. This implies we are not in the
// smallest binade or are dealing with denormals.
// 2. Our significand excluding the integral bit is all zeros.
- let crossing_binade_boundary = self.exp != S::MIN_EXP &&
- self.sig[0] & sig_mask == 0;
+ let crossing_binade_boundary =
+ self.exp != S::MIN_EXP && self.sig[0] & sig_mask == 0;
// Decrement the significand.
//
// the integral bit to 1, and increment the exponent. If we have a
// denormal always increment since moving denormals and the numbers in the
// smallest normal binade have the same exponent in our representation.
- let crossing_binade_boundary = !self.is_denormal() &&
- self.sig[0] & sig_mask == sig_mask;
+ let crossing_binade_boundary =
+ !self.is_denormal() && self.sig[0] & sig_mask == sig_mask;
if crossing_binade_boundary {
self.sig = [0];
category: Category::Normal,
sign: false,
marker: PhantomData,
- }.normalize(round, Loss::ExactlyZero)
+ }
+ .normalize(round, Loss::ExactlyZero)
}
fn from_str_r(mut s: &str, mut round: Round) -> Result<StatusAnd<Self>, ParseError> {
let mut loss = Loss::ExactlyZero;
if truncated_bits > 0 {
loss = Loss::through_truncation(&self.sig, truncated_bits);
- if loss != Loss::ExactlyZero &&
- self.round_away_from_zero(round, loss, truncated_bits)
+ if loss != Loss::ExactlyZero
+ && self.round_away_from_zero(round, loss, truncated_bits)
{
r = r.wrapping_add(1);
if r == 0 {
assert!(rhs.is_finite_non_zero());
// If exponents are equal, do an unsigned comparison of the significands.
- self.exp.cmp(&rhs.exp).then_with(
- || sig::cmp(&self.sig, &rhs.sig),
- )
+ self.exp.cmp(&rhs.exp).then_with(|| sig::cmp(&self.sig, &rhs.sig))
}
fn bitwise_eq(self, rhs: Self) -> bool {
}
fn is_denormal(self) -> bool {
- self.is_finite_non_zero() && self.exp == S::MIN_EXP &&
- !sig::get_bit(&self.sig, S::PRECISION - 1)
+ self.is_finite_non_zero()
+ && self.exp == S::MIN_EXP
+ && !sig::get_bit(&self.sig, S::PRECISION - 1)
}
fn is_signaling(self) -> bool {
let sig_bits = (S::PRECISION - 1) as ExpInt;
self.exp += sig_bits;
- self = self.normalize(Round::NearestTiesToEven, Loss::ExactlyZero)
- .value;
+ self = self.normalize(Round::NearestTiesToEven, Loss::ExactlyZero).value;
self.exp - sig_bits
}
fn is_x87_double_extended<S: Semantics>() -> bool {
S::QNAN_SIGNIFICAND == X87DoubleExtendedS::QNAN_SIGNIFICAND
}
- let x87_special_nan = is_x87_double_extended::<S>() && !is_x87_double_extended::<T>() &&
- r.category == Category::NaN &&
- (r.sig[0] & S::QNAN_SIGNIFICAND) != S::QNAN_SIGNIFICAND;
+ let x87_special_nan = is_x87_double_extended::<S>()
+ && !is_x87_double_extended::<T>()
+ && r.category == Category::NaN
+ && (r.sig[0] & S::QNAN_SIGNIFICAND) != S::QNAN_SIGNIFICAND;
// If this is a truncation of a denormal number, and the target semantics
// has larger exponent range than the source semantics (this can happen
// OMSB is numbered from 1. We want to place it in the integer
// bit numbered PRECISION if possible, with a compensating change in
// the exponent.
- let mut final_exp = self.exp.saturating_add(
- omsb as ExpInt - S::PRECISION as ExpInt,
- );
+ let mut final_exp = self.exp.saturating_add(omsb as ExpInt - S::PRECISION as ExpInt);
// If the resulting exponent is too high, overflow according to
// the rounding mode.
} else {
dec_exp = dec_exp.saturating_sub((last_sig_digit - dot) as i32);
}
- let significand_digits = last_sig_digit - first_sig_digit + 1 -
- (dot > first_sig_digit && dot < last_sig_digit) as usize;
+ let significand_digits = last_sig_digit - first_sig_digit + 1
+ - (dot > first_sig_digit && dot < last_sig_digit) as usize;
let normalized_exp = dec_exp.saturating_add(significand_digits as i32 - 1);
// Handle the cases where exponents are obviously too large or too
}
// Check for MIN_EXP.
- if normalized_exp.saturating_add(1).saturating_mul(28738) <=
- 8651 * (S::MIN_EXP as i32 - S::PRECISION as i32)
+ if normalized_exp.saturating_add(1).saturating_mul(28738)
+ <= 8651 * (S::MIN_EXP as i32 - S::PRECISION as i32)
{
// Underflow to zero and round.
- let r = if round == Round::TowardPositive {
- IeeeFloat::SMALLEST
- } else {
- IeeeFloat::ZERO
- };
+ let r =
+ if round == Round::TowardPositive { IeeeFloat::SMALLEST } else { IeeeFloat::ZERO };
return Ok((Status::UNDERFLOW | Status::INEXACT).and(r));
}
if power & 1 != 0 {
r_scratch.resize(r.len() + p5.len(), 0);
- let _: Loss = sig::mul(
- &mut r_scratch,
- &mut 0,
- &r,
- &p5,
- (r.len() + p5.len()) * LIMB_BITS,
- );
+ let _: Loss =
+ sig::mul(&mut r_scratch, &mut 0, &r, &p5, (r.len() + p5.len()) * LIMB_BITS);
while r_scratch.last() == Some(&0) {
r_scratch.pop();
}
used_bits = calc_precision.saturating_sub(truncated_bits);
}
// Extra half-ulp lost in reciprocal of exponent.
- half_ulp_err2 = 2 *
- (pow5_status != Status::OK || calc_loss != Loss::ExactlyZero) as Limb;
+ half_ulp_err2 =
+ 2 * (pow5_status != Status::OK || calc_loss != Loss::ExactlyZero) as Limb;
}
// Both sig::mul and sig::div return the
// than the returned value.
//
// See "How to Read Floating Point Numbers Accurately" by William D Clinger.
- assert!(
- half_ulp_err1 < 2 || half_ulp_err2 < 2 || (half_ulp_err1 + half_ulp_err2 < 8)
- );
+ assert!(half_ulp_err1 < 2 || half_ulp_err2 < 2 || (half_ulp_err1 + half_ulp_err2 < 8));
let inexact = (calc_loss != Loss::ExactlyZero) as Limb;
let half_ulp_err = if half_ulp_err1 + half_ulp_err2 == 0 {
/// Implementation details of IeeeFloat significands, such as big integer arithmetic.
/// As a rule of thumb, no functions in this module should dynamically allocate.
mod sig {
+ use super::{limbs_for_bits, ExpInt, Limb, Loss, LIMB_BITS};
use core::cmp::Ordering;
use core::mem;
- use super::{ExpInt, Limb, LIMB_BITS, limbs_for_bits, Loss};
pub(super) fn is_all_zeros(limbs: &[Limb]) -> bool {
limbs.iter().all(|&l| l == 0)
/// One, not zero, based LSB. That is, returns 0 for a zeroed significand.
pub(super) fn olsb(limbs: &[Limb]) -> usize {
- limbs.iter().enumerate().find(|(_, &limb)| limb != 0).map_or(0,
- |(i, limb)| i * LIMB_BITS + limb.trailing_zeros() as usize + 1)
+ limbs
+ .iter()
+ .enumerate()
+ .find(|(_, &limb)| limb != 0)
+ .map_or(0, |(i, limb)| i * LIMB_BITS + limb.trailing_zeros() as usize + 1)
}
/// One, not zero, based MSB. That is, returns 0 for a zeroed significand.
pub(super) fn omsb(limbs: &[Limb]) -> usize {
- limbs.iter().enumerate().rfind(|(_, &limb)| limb != 0).map_or(0,
- |(i, limb)| (i + 1) * LIMB_BITS - limb.leading_zeros() as usize)
+ limbs
+ .iter()
+ .enumerate()
+ .rfind(|(_, &limb)| limb != 0)
+ .map_or(0, |(i, limb)| (i + 1) * LIMB_BITS - limb.leading_zeros() as usize)
}
/// Comparison (unsigned) of two significands.
if precision <= omsb {
extract(dst, src, precision, omsb - precision);
- (
- Loss::through_truncation(src, omsb - precision),
- omsb as ExpInt - 1,
- )
+ (Loss::through_truncation(src, omsb - precision), omsb as ExpInt - 1)
} else {
extract(dst, src, omsb, 0);
(Loss::ExactlyZero, precision as ExpInt - 1)
// caller needs to call IeeeFloat::normalize() if normalized value is
// expected.
let omsb = omsb(dst);
- if omsb <= precision {
- Loss::ExactlyZero
- } else {
- shift_right(dst, exp, omsb - precision)
- }
+ if omsb <= precision { Loss::ExactlyZero } else { shift_right(dst, exp, omsb - precision) }
}
/// `quotient = dividend / divisor`. Returns the lost fraction.
divisor: &mut [Limb],
precision: usize,
) -> Loss {
-
// Normalize the divisor.
let bits = precision - omsb(divisor);
shift_left(divisor, &mut 0, bits);
}
// Helper for figuring out the lost fraction.
- let lost_fraction = |dividend: &[Limb], divisor: &[Limb]| {
- match cmp(dividend, divisor) {
- Ordering::Greater => Loss::MoreThanHalf,
- Ordering::Equal => Loss::ExactlyHalf,
- Ordering::Less => {
- if is_all_zeros(dividend) {
- Loss::ExactlyZero
- } else {
- Loss::LessThanHalf
- }
+ let lost_fraction = |dividend: &[Limb], divisor: &[Limb]| match cmp(dividend, divisor) {
+ Ordering::Greater => Loss::MoreThanHalf,
+ Ordering::Equal => Loss::ExactlyHalf,
+ Ordering::Less => {
+ if is_all_zeros(dividend) {
+ Loss::ExactlyZero
+ } else {
+ Loss::LessThanHalf
}
}
};
return lost_fraction(&[(rem as Limb) << 1], &[divisor as Limb]);
}
- }
+ };
}
try_short_div!(u32, u16, 16);