#![doc(hidden)]
use crate::bigint::{Bigint, Limb, LIMB_BITS};
use crate::extended_float::{extended_to_float, ExtendedFloat};
use crate::num::Float;
use crate::number::Number;
use crate::rounding::{round, round_down, round_nearest_tie_even};
use core::cmp;
#[inline]
pub fn slow<'a, F, Iter1, Iter2>(
num: Number,
fp: ExtendedFloat,
integer: Iter1,
fraction: Iter2,
) -> ExtendedFloat
where
F: Float,
Iter1: Iterator<Item = &'a u8> + Clone,
Iter2: Iterator<Item = &'a u8> + Clone,
{
debug_assert!(fp.mant & (1 << 63) != 0);
let sci_exp = scientific_exponent(&num);
let (bigmant, digits) = parse_mantissa(integer, fraction, F::MAX_DIGITS);
let exponent = sci_exp + 1 - digits as i32;
if exponent >= 0 {
positive_digit_comp::<F>(bigmant, exponent)
} else {
negative_digit_comp::<F>(bigmant, fp, exponent)
}
}
pub fn positive_digit_comp<F: Float>(mut bigmant: Bigint, exponent: i32) -> ExtendedFloat {
bigmant.pow(10, exponent as u32).unwrap();
let (mant, is_truncated) = bigmant.hi64();
let exp = bigmant.bit_length() as i32 - 64 + F::EXPONENT_BIAS;
let mut fp = ExtendedFloat {
mant,
exp,
};
round::<F, _>(&mut fp, |f, s| {
round_nearest_tie_even(f, s, |is_odd, is_halfway, is_above| {
is_above || (is_halfway && is_truncated) || (is_odd && is_halfway)
});
});
fp
}
#[allow(clippy::comparison_chain)]
pub fn negative_digit_comp<F: Float>(
bigmant: Bigint,
mut fp: ExtendedFloat,
exponent: i32,
) -> ExtendedFloat {
debug_assert!(fp.mant & (1 << 63) != 0);
let mut real_digits = bigmant;
let real_exp = exponent;
debug_assert!(real_exp < 0);
let mut b = fp;
round::<F, _>(&mut b, round_down);
let b = extended_to_float::<F>(b);
let theor = bh(b);
let mut theor_digits = Bigint::from_u64(theor.mant);
let theor_exp = theor.exp;
let binary_exp = theor_exp - real_exp;
let halfradix_exp = -real_exp;
if halfradix_exp != 0 {
theor_digits.pow(5, halfradix_exp as u32).unwrap();
}
if binary_exp > 0 {
theor_digits.pow(2, binary_exp as u32).unwrap();
} else if binary_exp < 0 {
real_digits.pow(2, (-binary_exp) as u32).unwrap();
}
let ord = real_digits.data.cmp(&theor_digits.data);
round::<F, _>(&mut fp, |f, s| {
round_nearest_tie_even(f, s, |is_odd, _, _| {
match ord {
cmp::Ordering::Greater => true,
cmp::Ordering::Less => false,
cmp::Ordering::Equal if is_odd => true,
cmp::Ordering::Equal => false,
}
});
});
fp
}
macro_rules! add_digit {
($c:ident, $value:ident, $counter:ident, $count:ident) => {{
let digit = $c - b'0';
$value *= 10 as Limb;
$value += digit as Limb;
$counter += 1;
$count += 1;
}};
}
macro_rules! add_temporary {
(@mul $result:ident, $power:expr, $value:expr) => {
$result.data.mul_small($power).unwrap();
$result.data.add_small($value).unwrap();
};
($format:ident, $result:ident, $counter:ident, $value:ident) => {
if $counter != 0 {
let small_power = unsafe { f64::int_pow_fast_path($counter, 10) };
add_temporary!(@mul $result, small_power as Limb, $value);
$counter = 0;
$value = 0;
}
};
(@end $format:ident, $result:ident, $counter:ident, $value:ident) => {
if $counter != 0 {
let small_power = unsafe { f64::int_pow_fast_path($counter, 10) };
add_temporary!(@mul $result, small_power as Limb, $value);
}
};
(@max $format:ident, $result:ident, $counter:ident, $value:ident, $max:ident) => {
add_temporary!(@mul $result, $max, $value);
$counter = 0;
$value = 0;
};
}
macro_rules! round_up_truncated {
($format:ident, $result:ident, $count:ident) => {{
add_temporary!(@mul $result, 10, 1);
$count += 1;
}};
}
macro_rules! round_up_nonzero {
($format:ident, $iter:expr, $result:ident, $count:ident) => {{
for &digit in $iter {
if digit != b'0' {
round_up_truncated!($format, $result, $count);
return ($result, $count);
}
}
}};
}
pub fn parse_mantissa<'a, Iter1, Iter2>(
mut integer: Iter1,
mut fraction: Iter2,
max_digits: usize,
) -> (Bigint, usize)
where
Iter1: Iterator<Item = &'a u8> + Clone,
Iter2: Iterator<Item = &'a u8> + Clone,
{
let mut counter: usize = 0;
let mut count: usize = 0;
let mut value: Limb = 0;
let mut result = Bigint::new();
let step: usize = if LIMB_BITS == 32 {
9
} else {
19
};
let max_native = (10 as Limb).pow(step as u32);
'integer: loop {
while counter < step && count < max_digits {
if let Some(&c) = integer.next() {
add_digit!(c, value, counter, count);
} else {
break 'integer;
}
}
if count == max_digits {
add_temporary!(@end format, result, counter, value);
round_up_nonzero!(format, integer, result, count);
round_up_nonzero!(format, fraction, result, count);
return (result, count);
} else {
add_temporary!(@max format, result, counter, value, max_native);
}
}
if count == 0 {
for &c in &mut fraction {
if c != b'0' {
add_digit!(c, value, counter, count);
break;
}
}
}
'fraction: loop {
while counter < step && count < max_digits {
if let Some(&c) = fraction.next() {
add_digit!(c, value, counter, count);
} else {
break 'fraction;
}
}
if count == max_digits {
add_temporary!(@end format, result, counter, value);
round_up_nonzero!(format, fraction, result, count);
return (result, count);
} else {
add_temporary!(@max format, result, counter, value, max_native);
}
}
add_temporary!(@end format, result, counter, value);
(result, count)
}
#[inline]
pub fn scientific_exponent(num: &Number) -> i32 {
let mut mantissa = num.mantissa;
let mut exponent = num.exponent;
while mantissa >= 10000 {
mantissa /= 10000;
exponent += 4;
}
while mantissa >= 100 {
mantissa /= 100;
exponent += 2;
}
while mantissa >= 10 {
mantissa /= 10;
exponent += 1;
}
exponent as i32
}
#[inline]
pub fn b<F: Float>(float: F) -> ExtendedFloat {
ExtendedFloat {
mant: float.mantissa(),
exp: float.exponent(),
}
}
#[inline]
pub fn bh<F: Float>(float: F) -> ExtendedFloat {
let fp = b(float);
ExtendedFloat {
mant: (fp.mant << 1) + 1,
exp: fp.exp - 1,
}
}