// Copyright 2020 CoD Technologies Corp.
//
// Licensed under the Apache License, Version 2.0 (the "License");
// you may not use this file except in compliance with the License.
// You may obtain a copy of the License at
//
// http://www.apache.org/licenses/LICENSE-2.0
//
// Unless required by applicable law or agreed to in writing, software
// distributed under the License is distributed on an "AS IS" BASIS,
// WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied.
// See the License for the specific language governing permissions and
// limitations under the License.
//! NumericVar.
use crate::binary::{
NumericDigit, NUMERIC_DIGIT_SIZE, NUMERIC_DSCALE_MAX, NUMERIC_HEADER_NDIGITS, NUMERIC_NAN,
NUMERIC_NEG, NUMERIC_POS, NUMERIC_WEIGHT_MAX, NUMERIC_WEIGHT_MIN, VAR_HEADER_SIZE,
};
use crate::data::{NumericData, NumericDigits};
use crate::num::{NumericBuf, DIVIDE_BY_ZERO_MSG};
use crate::typmod::Typmod;
use lazy_static::lazy_static;
use std::borrow::Cow;
use std::convert::{TryFrom, TryInto};
use std::f64::consts::{LN_10, LOG10_2, LOG10_E};
use std::fmt;
/// Limit on the precision (and hence scale) specifiable in a NUMERIC typmod.
/// Note that the implementation limit on the length of a numeric value is
/// much larger --- beware of what you use this for!
pub const NUMERIC_MAX_PRECISION: i32 = 1000;
// Internal limits on the scales chosen for calculation results
const NUMERIC_MAX_DISPLAY_SCALE: i32 = NUMERIC_MAX_PRECISION;
const NUMERIC_MIN_DISPLAY_SCALE: i32 = 0;
const NUMERIC_MAX_RESULT_SCALE: i32 = NUMERIC_MAX_PRECISION * 2;
/// For inherently inexact calculations such as division and square root,
/// we try to get at least this many significant digits; the idea is to
/// deliver a result no worse than f64 would.
const NUMERIC_MIN_SIG_DIGITS: i32 = 16;
pub const NBASE: i32 = 10000;
const HALF_NBASE: NumericDigit = 5000;
pub const DEC_DIGITS: i32 = 4;
const MUL_GUARD_DIGITS: i32 = 2;
const DIV_GUARD_DIGITS: i32 = 4;
const ROUND_POWERS: [NumericDigit; 4] = [0, 1000, 100, 10];
lazy_static! {
// 0
static ref ZERO: NumericVar<'static> = NumericVar::zero();
// 1
static ref ONE: NumericVar<'static> = NumericVar::borrowed(1, 0, 0, NUMERIC_POS, &[1]);
// 2
static ref TWO: NumericVar<'static> = NumericVar::borrowed(1, 0, 0, NUMERIC_POS, &[2]);
// 10
static ref TEN: NumericVar<'static> = NumericVar::borrowed(1, 0, 0, NUMERIC_POS, &[10]);
// 0.5
static ref ZERO_POINT_FIVE: NumericVar<'static> = NumericVar::borrowed(1, -1, 1, NUMERIC_POS, &[5000]);
// 0.9
static ref ZERO_POINT_NINE: NumericVar<'static> = NumericVar::borrowed(1, -1, 1, NUMERIC_POS, &[9000]);
// 1.1
static ref ONE_POINT_ONE: NumericVar<'static> = NumericVar::borrowed(2, 0, 1, NUMERIC_POS, &[1, 1000]);
}
/// `NumericVar` is the format we use for arithmetic. The `digits`-array part
/// is the same as the `NumericBinary` storage format, but the header is more
/// complex.
///
/// The value represented by a `Numeric` is determined by the `sign`, `weight`,
/// `ndigits`, and `digits[]` array.
///
/// Note: the first digit of a Numeric's value is assumed to be multiplied
/// by NBASE ** weight. Another way to say it is that there are weight+1
/// digits before the decimal point. It is possible to have weight < 0.
///
/// `data.buf` points at the physical start of the digit buffer for the
/// Numeric. `data.offset` points at the first digit in actual use (the one
/// with the specified weight). We normally leave an unused digit or two
/// (preset to zeroes) between buf and digits, so that there is room to store
/// a carry out of the top digit without reallocating space. We just need to
/// decrement digits (and increment weight) to make room for the carry digit.
/// (There is no such extra space in a numeric value stored in the database,
/// only in a Numeric in memory.)
///
/// `dscale`, or display scale, is the nominal precision expressed as number
/// of digits after the decimal point (it must always be >= 0 at present).
/// dscale may be more than the number of physically stored fractional digits,
/// implying that we have suppressed storage of significant trailing zeroes.
/// It should never be less than the number of stored digits, since that would
/// imply hiding digits that are present. NOTE that dscale is always expressed
/// in *decimal* digits, and so it may correspond to a fractional number of
/// base-NBASE digits --- divide by DEC_DIGITS to convert to NBASE digits.
///
/// While we consistently use `weight` to refer to the base-NBASE weight of
/// a numeric value, it is convenient in some scale-related calculations to
/// make use of the base-10 weight (ie, the approximate log10 of the value).
/// To avoid confusion, such a decimal-units weight is called a "dweight".
///
#[derive(Debug, Clone)]
pub(crate) struct NumericVar<'a> {
ndigits: i32,
weight: i32,
dscale: i32,
sign: u16,
digits: Cow<'a, NumericDigits>,
}
impl<'a> NumericVar<'a> {
/// Creates a `NumericVar` with `ndigits` data space.
#[inline]
fn with_ndigits(ndigits: i32) -> Self {
debug_assert!(ndigits >= 0);
NumericVar {
ndigits,
weight: 0,
dscale: 0,
sign: 0,
digits: Cow::Owned(NumericData::with_ndigits(ndigits)),
}
}
/// Creates a `NumericVar` of `NaN`, which has no data space.
#[inline]
pub fn nan() -> Self {
NumericVar::borrowed(0, 0, 0, NUMERIC_NAN, &[])
}
/// Creates a `NumericVar` of zero with given scale.
#[inline]
fn zero_scaled(scale: i32) -> NumericVar<'a> {
debug_assert!(scale >= 0 && scale <= NUMERIC_DSCALE_MAX as i32);
Self::borrowed(0, 0, scale, NUMERIC_POS, &[])
}
/// Creates a `NumericVar` with borrowed data space.
#[inline]
pub fn borrowed(
ndigits: i32,
weight: i32,
dscale: i32,
sign: u16,
digits: &'a [NumericDigit],
) -> Self {
debug_assert_eq!(ndigits as usize, digits.len());
let digits = Cow::Borrowed(unsafe { NumericDigits::from_slice_unchecked(digits) });
NumericVar {
ndigits,
weight,
dscale,
sign,
digits,
}
}
/// Creates a `NumericVar` with given data space.
#[inline]
pub fn owned(ndigits: i32, weight: i32, dscale: i32, sign: u16, digits: NumericData) -> Self {
debug_assert!(digits.offset() + ndigits as u32 <= digits.len());
NumericVar {
ndigits,
weight,
dscale,
sign,
digits: Cow::Owned(digits),
}
}
/// Creates a `NumericVar` of zero.
#[inline]
pub fn zero() -> Self {
NumericVar::borrowed(0, 0, 0, NUMERIC_POS, &[])
}
#[inline]
pub fn ndigits(&self) -> i32 {
self.ndigits
}
#[inline]
pub fn weight(&self) -> i32 {
self.weight
}
#[inline]
pub fn dscale(&self) -> i32 {
self.dscale
}
/// Checks if `self` is `NaN`.
#[inline]
pub const fn is_nan(&self) -> bool {
self.sign == NUMERIC_NAN
}
/// Checks if `self` is positive.
#[inline]
pub const fn is_positive(&self) -> bool {
self.sign == NUMERIC_POS
}
/// Checks if `self` is negative.
#[inline]
pub const fn is_negative(&self) -> bool {
self.sign == NUMERIC_NEG
}
#[inline]
pub fn digits(&self) -> &[NumericDigit] {
&self.digits.as_slice()[0..self.ndigits as usize]
}
#[inline]
pub fn into_numeric_buf(self) -> NumericBuf {
let mut data = self.digits.into_owned();
let header_offset = data.set_header(
self.weight as i16,
self.dscale as u16,
self.sign as u16,
self.ndigits,
);
let (buf, len, _) = data.into_raw_parts();
unsafe {
NumericBuf::from_raw_parts(buf as *const u8, len * NUMERIC_DIGIT_SIZE, header_offset)
}
}
/// Round the value of a variable to no more than rscale decimal digits
/// after the decimal point.
///
/// NOTE: we allow rscale < 0 here, implying rounding before the decimal point.
pub fn round_common(&mut self, rscale: i32) {
debug_assert!(!self.is_nan());
// decimal digits wanted
let di = (self.weight + 1) * DEC_DIGITS + rscale;
// If di = 0, the value loses all digits, but could round up to 1 if its
// first extra digit is >= 5. If di < 0 the result must be 0.
if di < 0 {
self.ndigits = 0;
self.weight = 0;
self.sign = NUMERIC_POS;
} else {
// NBASE digits wanted
let mut ndigits = (di + DEC_DIGITS - 1) / DEC_DIGITS;
// 0, or number of decimal digits to keep in last NBASE digit
let di = di % DEC_DIGITS;
if ndigits < self.ndigits || (ndigits == self.ndigits && di > 0) {
let data = self.digits.to_mut();
if self.ndigits > 0 {
data.reserve_rounding_digit(self.ndigits);
}
// Carry may need one additional digit
debug_assert!(data.offset() > NUMERIC_HEADER_NDIGITS || self.ndigits == 0);
let digits = data.digits_mut(self.ndigits);
self.ndigits = ndigits;
let mut carry: i32 = 0;
if di == 0 {
if digits[ndigits as usize] >= HALF_NBASE {
carry = 1;
}
} else {
// Must round within last NBASE digit
let mut pow10 = ROUND_POWERS[di as usize];
ndigits -= 1;
debug_assert!((ndigits as usize) < digits.len());
let digit = unsafe { digits.get_unchecked_mut(ndigits as usize) };
let extra = *digit % pow10;
*digit -= extra;
if extra >= pow10 / 2 {
pow10 += *digit;
if pow10 >= NBASE as NumericDigit {
pow10 -= NBASE as NumericDigit;
carry = 1;
}
*digit = pow10;
}
}
let offset = data.offset();
// Carry may need one additional digit, so we use buf from start.
let digits = data.as_mut_slice();
digits[offset as usize - 1] = 0;
// Propagate carry if needed
while carry > 0 {
ndigits -= 1;
let i = (offset as i32 + ndigits) as usize;
debug_assert!(i < digits.len());
let digit = unsafe { digits.get_unchecked_mut(i) };
carry += *digit as i32;
if carry >= NBASE as i32 {
*digit = (carry - NBASE as i32) as NumericDigit;
carry = 1;
} else {
*digit = carry as NumericDigit;
carry = 0;
}
}
if ndigits < 0 {
debug_assert_eq!(ndigits, -1);
debug_assert!(data.offset() > 0);
data.offset_sub(1);
self.ndigits += 1;
self.weight += 1;
}
}
}
self.dscale = rscale;
}
/// Truncate (towards zero) the value of a variable at rscale decimal digits
/// after the decimal point.
///
/// NOTE: we allow rscale < 0 here, implying truncation before the decimal point.
pub fn trunc_common(&mut self, rscale: i32) {
debug_assert!(!self.is_nan());
// decimal digits wanted
let di = (self.weight + 1) * DEC_DIGITS + rscale;
// If di <= 0, the value loses all digits.
if di <= 0 {
self.ndigits = 0;
self.weight = 0;
self.sign = NUMERIC_POS;
} else {
// NBASE digits wanted
let mut ndigits = (di + DEC_DIGITS - 1) / DEC_DIGITS;
if ndigits <= self.ndigits {
self.ndigits = ndigits;
// 0, or number of decimal digits to keep in last NBASE digit
let di = di % DEC_DIGITS;
if di > 0 {
let data = self.digits.to_mut();
let digits = data.digits_mut(self.ndigits);
let pow10 = ROUND_POWERS[di as usize];
ndigits -= 1;
let extra = digits[ndigits as usize] % pow10;
digits[ndigits as usize] -= extra;
}
}
}
self.dscale = rscale;
}
/// Return the smallest integer greater than or equal to the argument
/// on variable level
#[inline]
pub fn ceil_common(&self) -> Self {
debug_assert!(!self.is_nan());
let mut result = self.clone();
result.trunc_common(0);
if self.is_positive() && self.cmp_common(&result) != 0 {
result = result.add_common(&ONE);
}
result
}
/// Return the largest integer equal to or less than the argument
/// on variable level
#[inline]
pub fn floor_common(&self) -> Self {
debug_assert!(!self.is_nan());
let mut result = self.clone();
result.trunc_common(0);
if self.is_negative() && self.cmp_common(&result) != 0 {
result = result.sub_common(&ONE);
}
result
}
/// Strips the leading and trailing zeroes, and normalize zero.
pub fn strip(&mut self) {
let data = self.digits.to_mut();
let digits = data.digits(self.ndigits);
let mut ndigits = self.ndigits;
let mut i = 0;
// strip leading zeroes
while ndigits > 0 && unsafe { *digits.get_unchecked(i) } == 0 {
i += 1;
self.weight -= 1;
ndigits -= 1;
}
// strip trailing zeroes
while ndigits > 0 && unsafe { *digits.get_unchecked(i + ndigits as usize - 1) } == 0 {
ndigits -= 1;
}
// if it's zero, normalize the sign and weight
if ndigits == 0 {
self.sign = NUMERIC_POS;
self.weight = 0;
}
data.offset_add(i as u32);
self.ndigits = ndigits;
}
/// Add the absolute values of two variables into result.
pub fn add_abs(&self, other: &Self) -> Self {
debug_assert!(!self.is_nan());
debug_assert!(!other.is_nan());
// copy these values into local vars for speed in inner loop
let var1_ndigits = self.ndigits;
let var2_ndigits = other.ndigits;
let var1_digits = self.digits();
let var2_digits = other.digits();
let res_weight = self.weight.max(other.weight) + 1;
let res_dscale = self.dscale.max(other.dscale);
// Note: here we are figuring rscale in base-NBASE digits
let res_rscale = {
let rscale1 = self.ndigits - self.weight - 1;
let rscale2 = other.ndigits - other.weight - 1;
rscale1.max(rscale2)
};
let res_ndigits = {
let ndigits = res_rscale + res_weight + 1;
if ndigits > 0 {
ndigits
} else {
1
}
};
let mut res = Self::with_ndigits(res_ndigits);
let data = res.digits.to_mut();
let res_digits = data.digits_mut(res_ndigits);
let mut carry: NumericDigit = 0;
let mut i1 = res_rscale + self.weight + 1;
let mut i2 = res_rscale + other.weight + 1;
for i in (0..res_ndigits as usize).rev() {
i1 -= 1;
i2 -= 1;
if i1 >= 0 && i1 < var1_ndigits {
carry += unsafe { *var1_digits.get_unchecked(i1 as usize) };
}
if i2 >= 0 && i2 < var2_ndigits {
carry += unsafe { *var2_digits.get_unchecked(i2 as usize) };
}
let digit = unsafe { res_digits.get_unchecked_mut(i) };
if carry >= NBASE as NumericDigit {
*digit = carry - NBASE as NumericDigit;
carry = 1;
} else {
*digit = carry;
carry = 0;
}
}
debug_assert_eq!(carry, 0); // else we failed to allow for carry out
res.weight = res_weight;
res.dscale = res_dscale;
// Remove leading/trailing zeroes
res.strip();
res
}
/// Subtract the absolute value of `other` from the absolute value of `self`
/// and store in result.
///
/// NOTE: ABS(`self`) MUST BE GREATER OR EQUAL ABS(`other`) !!!
pub fn sub_abs(&self, other: &Self) -> Self {
debug_assert!(!self.is_nan());
debug_assert!(!other.is_nan());
// copy these values into local vars for speed in inner loop
let var1_ndigits = self.ndigits;
let var2_ndigits = other.ndigits;
let var1_digits = self.digits();
let var2_digits = other.digits();
let res_weight = self.weight;
let res_dscale = self.dscale.max(other.dscale);
// Note: here we are figuring rscale in base-NBASE digits
let res_rscale = {
let rscale1 = self.ndigits - self.weight - 1;
let rscale2 = other.ndigits - other.weight - 1;
rscale1.max(rscale2)
};
let res_ndigits = {
let ndigits = res_rscale + res_weight + 1;
if ndigits <= 0 {
1
} else {
ndigits
}
};
let mut res = Self::with_ndigits(res_ndigits);
let data = res.digits.to_mut();
let res_digits = data.digits_mut(res_ndigits);
let mut borrow: NumericDigit = 0;
let mut i1 = res_rscale + self.weight + 1;
let mut i2 = res_rscale + other.weight + 1;
for i in (0..res_ndigits as usize).rev() {
i1 -= 1;
i2 -= 1;
if i1 >= 0 && i1 < var1_ndigits {
borrow += unsafe { *var1_digits.get_unchecked(i1 as usize) };
}
if i2 >= 0 && i2 < var2_ndigits {
borrow -= unsafe { *var2_digits.get_unchecked(i2 as usize) };
}
let digit = unsafe { res_digits.get_unchecked_mut(i) };
if borrow < 0 {
*digit = borrow + NBASE as NumericDigit;
borrow = -1;
} else {
*digit = borrow;
borrow = 0;
}
}
debug_assert_eq!(borrow, 0); // else caller gave us self < other
res.weight = res_weight;
res.dscale = res_dscale;
// Remove leading/trailing zeroes
res.strip();
res
}
/// Compare the absolute values of `self` and `other`
///
/// * -1 for ABS(`self`) < ABS(`other`)
/// * 0 for ABS(`self`) == ABS(`other`)
/// * 1 for ABS(`self`) > ABS(`other`)
pub fn cmp_abs(&self, other: &Self) -> i32 {
debug_assert!(!self.is_nan());
debug_assert!(!other.is_nan());
let var1_ndigits = self.ndigits;
let var1_digits = self.digits();
let mut var1_weight = self.weight;
let var2_ndigits = other.ndigits;
let var2_digits = other.digits();
let mut var2_weight = other.weight;
let mut i1 = 0;
let mut i2 = 0;
// Check any digits before the first common digit
while var1_weight > var2_weight && i1 < var1_ndigits {
if unsafe { *var1_digits.get_unchecked(i1 as usize) } != 0 {
return 1;
}
i1 += 1;
var1_weight -= 1;
}
while var2_weight > var1_weight && i2 < var2_ndigits {
if unsafe { *var2_digits.get_unchecked(i2 as usize) } != 0 {
return -1;
}
i2 += 1;
var2_weight -= 1;
}
// At this point, either var1_weight == var2_weight or we've run out of digits
if var1_weight == var2_weight {
while i1 < var1_ndigits && i2 < var2_ndigits {
let stat = unsafe {
*var1_digits.get_unchecked(i1 as usize)
- *var2_digits.get_unchecked(i2 as usize)
};
if stat != 0 {
return if stat > 0 { 1 } else { -1 };
} else {
i1 += 1;
i2 += 1;
}
}
}
// At this point, we've run out of digits on one side or the other; so any
// remaining nonzero digits imply that side is larger
while i1 < var1_ndigits {
if unsafe { *var1_digits.get_unchecked(i1 as usize) } != 0 {
return 1;
}
i1 += 1;
}
while i2 < var2_ndigits {
if unsafe { *var2_digits.get_unchecked(i2 as usize) } != 0 {
return -1;
}
i2 += 1;
}
0
}
/// Full version of add functionality on variable level (handling signs).
pub fn add_common(&self, other: &Self) -> Self {
debug_assert!(!self.is_nan());
debug_assert!(!other.is_nan());
// Decide on the signs of the two variables what to do
if self.is_positive() {
if other.is_positive() {
// Both are positive
// result = +(ABS(self) + ABS(other))
let mut result = self.add_abs(other);
result.sign = NUMERIC_POS;
result
} else {
let cmp = self.cmp_abs(other);
match cmp {
0 => {
// ABS(self) == ABS(other)
// result = ZERO
Self::zero_scaled(self.dscale.max(other.dscale))
}
1 => {
// ABS(self) > ABS(other)
// result = +(ABS(self) - ABS(other))
let mut result = self.sub_abs(other);
result.sign = NUMERIC_POS;
result
}
-1 => {
// ABS(self) < ABS(other)
// result = -(ABS(other) - ABS(self))
let mut result = other.sub_abs(self);
result.sign = NUMERIC_NEG;
result
}
_ => panic!("invalid comparison result"),
}
}
} else if other.is_positive() {
// self is negative, other is positive
// Must compare absolute values
let cmp = self.cmp_abs(other);
match cmp {
0 => {
// ABS(self) == ABS(other)
// result = ZERO
Self::zero_scaled(self.dscale.max(other.dscale))
}
1 => {
// ABS(self) > ABS(other)
// result = -(ABS(self) - ABS(other))
let mut result = self.sub_abs(other);
result.sign = NUMERIC_NEG;
result
}
-1 => {
// ABS(self) < ABS(other)
// result = +(ABS(other) - ABS(self))
let mut result = other.sub_abs(self);
result.sign = NUMERIC_POS;
result
}
_ => panic!("invalid comparison result"),
}
} else {
// Both are negative
// result = -(ABS(self) + ABS(other))
let mut result = self.add_abs(other);
result.sign = NUMERIC_NEG;
result
}
}
/// Full version of sub functionality on variable level (handling signs).
pub fn sub_common(&self, other: &Self) -> Self {
debug_assert!(!self.is_nan());
debug_assert!(!other.is_nan());
// Decide on the signs of the two variables what to do
if self.is_positive() {
if other.is_negative() {
// self is positive, other is negative
// result = +(ABS(self) + ABS(other))
let mut result = self.add_abs(other);
result.sign = NUMERIC_POS;
result
} else {
// Both are positive
// Must compare absolute values
let cmp = self.cmp_abs(other);
match cmp {
0 => {
// ABS(self) == ABS(other)
// result = ZERO
Self::zero_scaled(self.dscale.max(other.dscale))
}
1 => {
// ABS(self) > ABS(other)
// result = +(ABS(self) - ABS(other))
let mut result = self.sub_abs(other);
result.sign = NUMERIC_POS;
result
}
-1 => {
// ABS(self) < ABS(other)
// result = -(ABS(other) - ABS(self))
let mut result = other.sub_abs(self);
result.sign = NUMERIC_NEG;
result
}
_ => panic!("invalid comparison result"),
}
}
} else if other.is_negative() {
// Both are negative
// Must compare absolute values
let cmp = self.cmp_abs(other);
match cmp {
0 => {
// ABS(self) == ABS(other)
// result = ZERO
Self::zero_scaled(self.dscale.max(other.dscale))
}
1 => {
// ABS(self) > ABS(other)
// result = -(ABS(self) - ABS(other))
let mut result = self.sub_abs(other);
result.sign = NUMERIC_NEG;
result
}
-1 => {
// ABS(self) < ABS(other)
// result = +(ABS(other) - ABS(self))
let mut result = other.sub_abs(self);
result.sign = NUMERIC_POS;
result
}
_ => panic!("invalid comparison result"),
}
} else {
// var1 is negative, var2 is positive
// result = -(ABS(self) + ABS(other))
let mut result = self.add_abs(other);
result.sign = NUMERIC_NEG;
result
}
}
/// Multiplication on variable level.
/// Product of self * other is returned.
/// Result is rounded to no more than rscale fractional digits.
pub fn mul_common(&self, other: &Self, rscale: i32) -> Self {
debug_assert!(!self.is_nan());
debug_assert!(!other.is_nan());
// Arrange for var1 to be the shorter of the two numbers. This improves
// performance because the inner multiplication loop is much simpler than
// the outer loop, so it's better to have a smaller number of iterations
// of the outer loop. This also reduces the number of times that the
// accumulator array needs to be normalized.
let (var1, var2) = if self.ndigits > other.ndigits {
(other, self)
} else {
(self, other)
};
// copy these values into local vars for speed in inner loop
let var1_ndigits = var1.ndigits;
let var1_digits = var1.digits();
let var2_ndigits = var2.ndigits;
let var2_digits = var2.digits();
if var1_ndigits == 0 || var2_ndigits == 0 {
// one or both inputs is zero; so is result
return Self::zero_scaled(rscale);
}
// Determine result sign and (maximum possible) weight
let res_sign = if var1.sign == var2.sign {
NUMERIC_POS
} else {
NUMERIC_NEG
};
let res_weight = var1.weight + var2.weight + 2;
// Determine the number of result digits to compute. If the exact result
// would have more than rscale fractional digits, truncate the computation
// with MUL_GUARD_DIGITS guard digits, i.e., ignore input digits that
// would only contribute to the right of that. (This will give the exact
// rounded-to-rscale answer unless carries out of the ignored positions
// would have propagated through more than MUL_GUARD_DIGITS digits.)
//
// Note: an exact computation could not produce more than var1ndigits +
// var2ndigits digits, but we allocate one extra output digit in case
// rscale-driven rounding produces a carry out of the highest exact digit.
let res_ndigits = {
let ndigits = var1.ndigits + var2.ndigits + 1;
let max_digits =
res_weight + 1 + (rscale + DEC_DIGITS - 1) / DEC_DIGITS + MUL_GUARD_DIGITS;
ndigits.min(max_digits)
};
if res_ndigits < 3 {
// All input digits will be ignored; so result is zero
return Self::zero_scaled(rscale);
}
// We do the arithmetic in an array "dig[]" of signed int32's. Since
// INT32_MAX is noticeably larger than NBASE*NBASE, this gives us headroom
// to avoid normalizing carries immediately.
//
// max_dig tracks the maximum possible value of any dig[] entry; when this
// threatens to exceed INT32_MAX, we take the time to propagate carries.
// Furthermore, we need to ensure that overflow doesn't occur during the
// carry propagation passes either. The carry values could be as much as
// INT32_MAX/NBASE, so really we must normalize when digits threaten to
// exceed INT32_MAX - INT32_MAX/NBASE.
//
// To avoid overflow in max_dig itself, it actually represents the max
// possible value divided by NBASE-1, ie, at the top of the loop it is
// known that no dig[] entry exceeds max_dig * (NBASE-1).
let mut dig = vec![0; res_ndigits as usize * std::mem::size_of::<i32>()];
let mut max_dig = 0i32;
// The least significant digits of var1 should be ignored if they don't
// contribute directly to the first res_ndigits digits of the result that
// we are computing.
//
// Digit i1 of var1 and digit i2 of var2 are multiplied and added to digit
// i1+i2+2 of the accumulator array, so we need only consider digits of
// var1 for which i1 <= res_ndigits - 3.
let bound = (var1_ndigits - 1).min(res_ndigits - 3);
for i1 in (0..=bound).rev() {
let var1_digit = unsafe { *var1_digits.get_unchecked(i1 as usize) } as i32;
if var1_digit == 0 {
continue;
}
// Time to normalize?
max_dig += var1_digit;
if max_dig > ((i32::max_value() - i32::max_value() / NBASE) / (NBASE - 1)) {
// Yes, do it
let mut carry = 0;
for i in (0..res_ndigits).rev() {
let d = unsafe { dig.get_unchecked_mut(i as usize) };
let mut new_dig = *d + carry;
if new_dig >= NBASE {
carry = new_dig / NBASE;
new_dig -= carry * NBASE;
} else {
carry = 0;
}
*d = new_dig;
}
debug_assert_eq!(carry, 0);
// Reset max_dig to indicate new worst-case
max_dig = 1 + var1_digit;
}
// Add the appropriate multiple of var2 into the accumulator.
//
// As above, digits of var2 can be ignored if they don't contribute,
// so we only include digits for which i1+i2+2 <= res_ndigits - 1.
let bound = (var2_ndigits - 1).min(res_ndigits - i1 - 3);
let mut i = i1 + bound + 2;
for i2 in (0..=bound).rev() {
let d = unsafe { dig.get_unchecked_mut(i as usize) };
*d += var1_digit * unsafe { *var2_digits.get_unchecked(i2 as usize) } as i32;
i -= 1;
}
}
// Now we do a final carry propagation pass to normalize the result, which
// we combine with storing the result digits into the output. Note that
// this is still done at full precision w/guard digits.
let mut result = Self::with_ndigits(res_ndigits);
let data = result.digits.to_mut();
let res_digits = data.digits_mut(res_ndigits);
let mut carry = 0;
for i in (0..res_ndigits).rev() {
let mut new_dig = unsafe { dig.get_unchecked(i as usize) } + carry;
if new_dig >= NBASE {
carry = new_dig / NBASE;
new_dig -= carry * NBASE;
} else {
carry = 0;
}
let res_digit = unsafe { res_digits.get_unchecked_mut(i as usize) };
*res_digit = new_dig as NumericDigit;
}
debug_assert_eq!(carry, 0);
// Finally, round the result to the requested precision.
result.weight = res_weight;
result.sign = res_sign;
// Round to target rscale (and set result->dscale)
result.round_common(rscale);
// Strip leading and trailing zeroes
result.strip();
result
}
/// Default scale selection for division
///
/// Returns the appropriate result scale for the division result.
pub fn select_div_scale(&self, other: &Self) -> i32 {
// The result scale of a division isn't specified in any SQL standard. For
// PostgreSQL we select a result scale that will give at least
// NUMERIC_MIN_SIG_DIGITS significant digits, so that numeric gives a
// result no less accurate than f64; but use a scale not less than
// either input's display scale.
// Get the actual (normalized) weight and first digit of each input.
let mut weight1 = 0; // values to use if self is zero
let mut first_digit1 = 0;
let var1_digits = self.digits();
for i in 0..self.ndigits {
first_digit1 = unsafe { *var1_digits.get_unchecked(i as usize) };
if first_digit1 != 0 {
weight1 = self.weight - i;
break;
}
}
let mut weight2 = 0; // values to use if other is zero
let mut first_digit2 = 0;
let var2_digits = other.digits();
for i in 0..other.ndigits {
first_digit2 = unsafe { *var2_digits.get_unchecked(i as usize) };
if first_digit2 != 0 {
weight2 = other.weight - i;
break;
}
}
// Estimate weight of quotient. If the two first digits are equal, we
// can't be sure, but assume that self is less than other.
let qweight = {
let mut w = weight1 - weight2;
if first_digit1 <= first_digit2 {
w -= 1;
}
w
};
// Select result scale
(NUMERIC_MIN_SIG_DIGITS - qweight * DEC_DIGITS)
.max(self.dscale)
.max(other.dscale)
.max(NUMERIC_MIN_DISPLAY_SCALE)
.min(NUMERIC_MAX_DISPLAY_SCALE)
}
/// Division on variable level. Quotient of `self` / `other` is returned.
/// The quotient is figured to exactly rscale fractional digits.
/// If round is true, it is rounded at the rscale'th digit; if false, it
/// is truncated (towards zero) at that digit.
///
/// Returns `None` if `other == 0`.
#[allow(clippy::cognitive_complexity)]
pub fn div_common(&self, other: &Self, rscale: i32, round: bool) -> Option<Self> {
debug_assert!(!self.is_nan());
debug_assert!(!other.is_nan());
// copy these values into local vars for speed in inner loop
let var1_ndigits = self.ndigits;
let var2_ndigits = other.ndigits;
// First of all division by zero check; we must not be handed an
// unnormalized divisor.
if var2_ndigits == 0 || other.digits[0] == 0 {
return None;
}
// Now result zero check
if var1_ndigits == 0 {
return Some(Self::zero_scaled(rscale));
}
// Determine the result sign, weight and number of digits to calculate.
// The weight figured here is correct if the emitted quotient has no
// leading zero digits; otherwise strip() will fix things up.
let res_sign = if self.sign == other.sign {
NUMERIC_POS
} else {
NUMERIC_NEG
};
let res_weight = self.weight - other.weight;
let res_ndigits = {
// The number of accurate result digits we need to produce
let mut ndigits = res_weight + 1 + (rscale + DEC_DIGITS - 1) / DEC_DIGITS;
// ... but always at least 1
ndigits = ndigits.max(1);
// If rounding needed, figure one more digit to ensure correct result
if round {
ndigits += 1;
}
ndigits
};
// The working dividend normally requires res_ndigits + var2_ndigits
// digits, but make it at least var1_ndigits so we can load all of var1
// into it. (There will be an additional digit dividend[0] in the
// dividend space, but for consistency with Knuth's notation we don't
// count that in div_ndigits.)
let div_ndigits = {
let ndigits = res_ndigits + var2_ndigits;
ndigits.max(var1_ndigits)
};
// We need a workspace with room for the working dividend (div_ndigits+1
// digits) plus room for the possibly-normalized divisor (var2_ndigits
// digits). It is convenient also to have a zero at divisor[0] with the
// actual divisor data in divisor[1 .. var2_ndigits]. Transferring the
// digits into the workspace also allows us to realloc the result (which
// might be the same as either input var) before we begin the main loop.
// Note that we use palloc0 to ensure that divisor[0], dividend[0], and
// any additional dividend positions beyond var1_ndigits, start out 0.
let mut workspace = vec![0 as NumericDigit; (div_ndigits + var2_ndigits + 2) as usize];
let (dividend, divisor) = workspace.split_at_mut(div_ndigits as usize + 1);
dividend[1..=var1_ndigits as usize].copy_from_slice(self.digits());
divisor[1..=var2_ndigits as usize].copy_from_slice(other.digits());
// Now we can alloc the result to hold the generated quotient digits.
let mut result = Self::with_ndigits(res_ndigits);
let data = result.digits.to_mut();
let res_digits = data.digits_mut(res_ndigits);
if var2_ndigits == 1 {
// If there's only a single divisor digit, we can use a fast path (cf.
// Knuth section 4.3.1 exercise 16).
let divisor1 = divisor[1] as i32;
let mut carry = 0i32;
for i in 0..res_ndigits as usize {
carry = carry * NBASE + unsafe { *dividend.get_unchecked(i + 1) } as i32;
let res_digit = unsafe { res_digits.get_unchecked_mut(i) };
*res_digit = (carry / divisor1) as NumericDigit;
carry %= divisor1;
}
} else {
// The full multiple-place algorithm is taken from Knuth volume 2,
// Algorithm 4.3.1D.
//
// We need the first divisor digit to be >= NBASE/2. If it isn't,
// make it so by scaling up both the divisor and dividend by the
// factor "d". (The reason for allocating dividend[0] above is to
// leave room for possible carry here.)
if divisor[1] < HALF_NBASE {
let d = NBASE / (divisor[1] + 1) as i32;
let mut carry = 0i32;
for i in (1..=var2_ndigits as usize).rev() {
let div = unsafe { divisor.get_unchecked_mut(i) };
carry += *div as i32 * d;
*div = (carry % NBASE) as NumericDigit;
carry /= NBASE;
}
debug_assert_eq!(carry, 0);
carry = 0;
// at this point only var1_ndigits of dividend can be nonzero
for i in (0..=var1_ndigits as usize).rev() {
let div = unsafe { dividend.get_unchecked_mut(i) };
carry += *div as i32 * d;
*div = (carry % NBASE) as NumericDigit;
carry /= NBASE;
}
debug_assert_eq!(carry, 0);
debug_assert!(divisor[1] >= HALF_NBASE);
}
// First 2 divisor digits are used repeatedly in main loop
let divisor1 = divisor[1];
let divisor2 = divisor[2];
// Begin the main loop. Each iteration of this loop produces the j'th
// quotient digit by dividing dividend[j .. j + var2ndigits] by the
// divisor; this is essentially the same as the common manual
// procedure for long division.
for (j, res_digit) in res_digits.iter_mut().enumerate() {
// Estimate quotient digit from the first two dividend digits
let next2digits = unsafe {
*dividend.get_unchecked(j) as i32 * NBASE
+ *dividend.get_unchecked(j + 1) as i32
};
// If next2digits are 0, then quotient digit must be 0 and there's
// no need to adjust the working dividend. It's worth testing
// here to fall out ASAP when processing trailing zeroes in a
// dividend.
if next2digits == 0 {
*res_digit = 0;
continue;
}
let mut qhat = if unsafe { *dividend.get_unchecked(j) } == divisor1 {
NBASE - 1
} else {
next2digits / divisor1 as i32
};
// Adjust quotient digit if it's too large. Knuth proves that
// after this step, the quotient digit will be either correct or
// just one too large. (Note: it's OK to use dividend[j+2] here
// because we know the divisor length is at least 2.)
while divisor2 as i32 * qhat
> (next2digits - qhat * divisor1 as i32) * NBASE
+ unsafe { *dividend.get_unchecked(j + 2) } as i32
{
qhat -= 1;
}
// As above, need do nothing more when quotient digit is 0
if qhat > 0 {
// Multiply the divisor by qhat, and subtract that from the
// working dividend. "carry" tracks the multiplication,
// "borrow" the subtraction (could we fold these together?)
let mut carry = 0;
let mut borrow = 0;
for i in (0..=var2_ndigits as usize).rev() {
carry += unsafe { *divisor.get_unchecked(i) } as i32 * qhat;
borrow -= carry % NBASE;
carry /= NBASE;
let div = unsafe { dividend.get_unchecked_mut(j + i) };
borrow += *div as i32;
if borrow < 0 {
*div = (borrow + NBASE) as NumericDigit;
borrow = -1;
} else {
*div = borrow as NumericDigit;
borrow = 0;
}
}
debug_assert_eq!(carry, 0);
// If we got a borrow out of the top dividend digit, then
// indeed qhat was one too large. Fix it, and add back the
// divisor to correct the working dividend. (Knuth proves
// that this will occur only about 3/NBASE of the time; hence,
// it's a good idea to test this code with small NBASE to be
// sure this section gets exercised.)
if borrow != 0 {
qhat -= 1;
carry = 0;
for i in (0..=var2_ndigits as usize).rev() {
let div = unsafe { dividend.get_unchecked_mut(j + i) };
carry += *div as i32 + unsafe { *divisor.get_unchecked(i) } as i32;
if carry >= NBASE {
*div = (carry - NBASE) as NumericDigit;
carry = 1;
} else {
*div = carry as NumericDigit;
carry = 0;
}
}
// A carry should occur here to cancel the borrow above
debug_assert_eq!(carry, 1);
}
}
// And we're done with this quotient digit
*res_digit = qhat as NumericDigit;
}
}
// Finally, round or truncate the result to the requested precision.
result.weight = res_weight;
result.sign = res_sign;
// Round or truncate to target rscale (and set result->dscale)
if round {
result.round_common(rscale);
} else {
result.trunc_common(rscale);
}
// Strip leading and trailing zeroes
result.strip();
Some(result)
}
/// This has the same API as `div_common()`, but is implemented using the division
/// algorithm from the "FM" library, rather than Knuth's schoolbook-division
/// approach. This is significantly faster but can produce inaccurate
/// results, because it sometimes has to propagate rounding to the left,
/// and so we can never be entirely sure that we know the requested digits
/// exactly. We compute DIV_GUARD_DIGITS extra digits, but there is
/// no certainty that that's enough. We use this only in the transcendental
/// function calculation routines, where everything is approximate anyway.
///
/// Although we provide a "round" argument for consistency with `div()`,
/// it is unwise to use this function with round=false. In truncation mode
/// it is possible to get a result with no significant digits, for example
/// with rscale=0 we might compute 0.99999... and truncate that to 0 when
/// the correct answer is 1.
///
/// Returns `None` if `other == 0`.
#[allow(clippy::cognitive_complexity)]
pub fn div_fast_common(&self, other: &Self, rscale: i32, round: bool) -> Option<Self> {
debug_assert!(!self.is_nan());
debug_assert!(!other.is_nan());
// copy these values into local vars for speed in inner loop
let var1_ndigits = self.ndigits;
let var1_digits = self.digits();
let var2_ndigits = other.ndigits;
let var2_digits = other.digits();
// First of all division by zero check; we must not be handed an
// unnormalized divisor.
if var2_ndigits == 0 || var2_digits[0] == 0 {
return None;
}
// Now result zero check
if var1_ndigits == 0 {
return Some(Self::zero_scaled(rscale));
}
// Determine the result sign, weight and number of digits to calculate
let res_sign = if self.sign == other.sign {
NUMERIC_POS
} else {
NUMERIC_NEG
};
let res_weight = self.weight - other.weight + 1;
let div_ndigits = {
// The number of accurate result digits we need to produce
let mut ndigits = res_weight + 1 + (rscale + DEC_DIGITS - 1) / DEC_DIGITS;
// Add guard digits for roundoff error
ndigits += DIV_GUARD_DIGITS;
if ndigits < DIV_GUARD_DIGITS {
ndigits = DIV_GUARD_DIGITS;
}
// Must be at least var1_ndigits, too, to simplify data-loading loop
if ndigits < var1_ndigits {
ndigits = var1_ndigits;
}
ndigits
};
// We do the arithmetic in an array "div[]" of signed int32's. Since
// INT32_MAX is noticeably larger than NBASE*NBASE, this gives us headroom
// to avoid normalizing carries immediately.
//
// We start with div[] containing one zero digit followed by the
// dividend's digits (plus appended zeroes to reach the desired precision
// including guard digits). Each step of the main loop computes an
// (approximate) quotient digit and stores it into div[], removing one
// position of dividend space. A final pass of carry propagation takes
// care of any mistaken quotient digits.
let mut div = vec![0i32; div_ndigits as usize + 1];
for i in 0..var1_ndigits as usize {
let d = unsafe { div.get_unchecked_mut(i + 1) };
*d = unsafe { *var1_digits.get_unchecked(i) } as i32;
}
// We estimate each quotient digit using floating-point arithmetic, taking
// the first four digits of the (current) dividend and divisor. This must
// be float to avoid overflow. The quotient digits will generally be off
// by no more than one from the exact answer.
let mut fdivisor = var2_digits[0] as f64;
for i in 1..4 {
fdivisor *= NBASE as f64;
if i < var2_ndigits {
fdivisor += unsafe { *var2_digits.get_unchecked(i as usize) } as f64;
}
}
let fdivisor_inverse = 1.0 / fdivisor;
// max_div tracks the maximum possible absolute value of any div[] entry;
// when this threatens to exceed INT_MAX, we take the time to propagate
// carries. Furthermore, we need to ensure that overflow doesn't occur
// during the carry propagation passes either. The carry values may have
// an absolute value as high as INT_MAX/NBASE + 1, so really we must
// normalize when digits threaten to exceed INT_MAX - INT_MAX/NBASE - 1.
//
// To avoid overflow in max_div itself, it represents the max absolute
// value divided by NBASE-1, ie, at the top of the loop it is known that
// no div[] entry has an absolute value exceeding max_div * (NBASE-1).
//
// Actually, though, that holds good only for div[] entries after div[qi];
// the adjustment done at the bottom of the loop may cause div[qi + 1] to
// exceed the max_div limit, so that div[qi] in the next iteration is
// beyond the limit. This does not cause problems, as explained below.
let mut max_div = 1;
// Outer loop computes next quotient digit, which will go into div[qi]
for qi in 0..div_ndigits as usize {
// Approximate the current dividend value
let mut fdividend = unsafe { *div.get_unchecked(qi) } as f64;
for i in 1..4usize {
fdividend *= NBASE as f64;
if (qi + i) as i32 <= div_ndigits {
fdividend += unsafe { *div.get_unchecked(qi + i) } as f64;
}
}
// Compute the (approximate) quotient digit
let mut fquotient = fdividend * fdivisor_inverse;
let mut qdigit = if fquotient >= 0.0 {
fquotient as i32
} else {
// truncate towards -infinity
fquotient as i32 - 1
};
if qdigit != 0 {
// Do we need to normalize now?
max_div += qdigit.abs();
if max_div > (i32::max_value() - i32::max_value() / NBASE - 1) / (NBASE - 1) {
// Yes, do it
let mut carry = 0;
let mut new_dig;
for i in (qi + 1..=div_ndigits as usize).rev() {
let div_i = unsafe { div.get_unchecked_mut(i) };
new_dig = *div_i + carry;
if new_dig < 0 {
carry = -((-new_dig - 1) / NBASE) - 1;
new_dig -= carry * NBASE;
} else if new_dig >= NBASE {
carry = new_dig / NBASE;
new_dig -= carry * NBASE;
} else {
carry = 0;
}
*div_i = new_dig;
}
let div_qi = unsafe { div.get_unchecked_mut(qi) };
new_dig = *div_qi + carry;
*div_qi = new_dig;
// All the div[] digits except possibly div[qi] are now in the
// range 0..NBASE-1. We do not need to consider div[qi] in
// the max_div value anymore, so we can reset max_div to 1.
max_div = 1;
// Recompute the quotient digit since new info may have
// propagated into the top four dividend digits
fdividend = *div_qi as f64;
for i in 1..4usize {
fdividend *= NBASE as f64;
if (qi + i) as i32 <= div_ndigits {
fdividend += unsafe { *div.get_unchecked(qi + i) } as f64;
}
}
// Compute the (approximate) quotient digit
fquotient = fdividend * fdivisor_inverse;
qdigit = if fquotient >= 0.0 {
fquotient as i32
} else {
// truncate towards -infinity
fquotient as i32 - 1
};
max_div += qdigit.abs();
}
// Subtract off the appropriate multiple of the divisor.
//
// The digits beyond div[qi] cannot overflow, because we know they
// will fall within the max_div limit. As for div[qi] itself, note
// that qdigit is approximately trunc(div[qi] / vardigits[0]),
// which would make the new value simply div[qi] mod vardigits[0].
// The lower-order terms in qdigit can change this result by not
// more than about twice INT_MAX/NBASE, so overflow is impossible.
if qdigit != 0 {
let istop = var2_ndigits.min(div_ndigits - qi as i32 + 1);
for i in 0..istop as usize {
let div_qi_i = unsafe { div.get_unchecked_mut(qi + i) };
*div_qi_i -= qdigit * unsafe { *var2_digits.get_unchecked(i) } as i32;
}
}
}
// The dividend digit we are about to replace might still be nonzero.
// Fold it into the next digit position.
//
// There is no risk of overflow here, although proving that requires
// some care. Much as with the argument for div[qi] not overflowing,
// if we consider the first two terms in the numerator and denominator
// of qdigit, we can see that the final value of div[qi + 1] will be
// approximately a remainder mod (vardigits[0]*NBASE + vardigits[1]).
// Accounting for the lower-order terms is a bit complicated but ends
// up adding not much more than INT_MAX/NBASE to the possible range.
// Thus, div[qi + 1] cannot overflow here, and in its role as div[qi]
// in the next loop iteration, it can't be large enough to cause
// overflow in the carry propagation step (if any), either.
//
// But having said that: div[qi] can be more than INT_MAX/NBASE, as
// noted above, which means that the product div[qi] * NBASE *can*
// overflow. When that happens, adding it to div[qi + 1] will always
// cause a canceling overflow so that the end result is correct. We
// could avoid the intermediate overflow by doing the multiplication
// and addition in int64 arithmetic, but so far there appears no need.
let div_qi = unsafe { *div.get_unchecked(qi) };
let div_qi_1 = unsafe { div.get_unchecked_mut(qi + 1) };
*div_qi_1 += div_qi * NBASE;
let div_qi = unsafe { div.get_unchecked_mut(qi) };
*div_qi = qdigit;
}
// Approximate and store the last quotient digit (div[div_ndigits])
let mut fdividend = div[div_ndigits as usize] as f64;
for _ in 1..4usize {
fdividend *= NBASE as f64;
}
let fquotient = fdividend * fdivisor_inverse;
let qdigit = if fquotient >= 0.0 {
fquotient as i32
} else {
// truncate towards -infinity
fquotient as i32 - 1
};
div[div_ndigits as usize] = qdigit;
// Because the quotient digits might be off by one, some of them might be
// -1 or NBASE at this point. The represented value is correct in a
// mathematical sense, but it doesn't look right. We do a final carry
// propagation pass to normalize the digits, which we combine with storing
// the result digits into the output. Note that this is still done at
// full precision w/guard digits.
let mut result = Self::with_ndigits(div_ndigits + 1);
let data = result.digits.to_mut();
let res_digits = data.digits_mut(result.ndigits);
let mut carry = 0;
for i in (0..=div_ndigits as usize).rev() {
let mut new_dig = unsafe { *div.get_unchecked(i) } + carry;
if new_dig < 0 {
carry = -((-new_dig - 1) / NBASE) - 1;
new_dig -= carry * NBASE;
} else if new_dig >= NBASE {
carry = new_dig / NBASE;
new_dig -= carry * NBASE;
} else {
carry = 0;
}
let res_digit_i = unsafe { res_digits.get_unchecked_mut(i) };
*res_digit_i = new_dig as NumericDigit;
}
debug_assert_eq!(carry, 0);
// Finally, round the result to the requested precision.
result.weight = res_weight;
result.sign = res_sign;
// Round to target rscale (and set result->dscale)
if round {
result.round_common(rscale);
} else {
result.trunc_common(rscale);
}
// Strip leading and trailing zeroes
result.strip();
Some(result)
}
/// Calculate the modulo of two numerics at variable level.
#[inline]
pub fn mod_common(&self, other: &Self) -> Option<Self> {
debug_assert!(!self.is_nan());
debug_assert!(!other.is_nan());
// We do this using the equation
// mod(x,y) = x - trunc(x/y)*y
// div() can be persuaded to give us trunc(x/y) directly.
let mut result = self.div_common(other, 0, false)?;
result = result.mul_common(other, other.dscale);
result = self.sub_common(&result);
Some(result)
}
/// Compare two values on variable level.
/// We assume zeroes have been truncated to no digits.
#[inline]
pub fn cmp_common(&self, other: &Self) -> i32 {
debug_assert!(!self.is_nan());
debug_assert!(!other.is_nan());
if self.ndigits == 0 {
if other.ndigits == 0 {
0
} else if other.is_negative() {
1
} else {
-1
}
} else if other.ndigits == 0 {
if self.is_positive() {
1
} else {
-1
}
} else if self.is_positive() {
if other.is_negative() {
1
} else {
self.cmp_abs(other)
}
} else if other.is_positive() {
-1
} else {
other.cmp_abs(self)
}
}
/// Compute the square root of x using Newton's algorithm.
pub fn sqrt_common(&self, rscale: i32) -> Self {
debug_assert!(self.is_positive());
let local_rscale = rscale + 8;
if self.ndigits == 0 {
return Self::zero_scaled(rscale);
}
// Initialize the result to the first guess
let mut result = Self::with_ndigits(1);
let data = result.digits.to_mut();
data.digits_mut(1)[0] = {
let i = self.digits[0] / 2;
if i == 0 {
1
} else {
i
}
};
result.weight = self.weight / 2;
let mut last_val = result.clone();
loop {
let val = self
.div_fast_common(&result, local_rscale, true)
.expect("should not be zero");
result = result.add_common(&val);
result = result.mul_common(&ZERO_POINT_FIVE, local_rscale);
if result.cmp_common(&last_val) == 0 {
break;
}
last_val = result.clone();
}
// Round to requested precision
result.round_common(rscale);
result
}
/// Raise `self` to the power of exp, where exp is an integer.
///
/// Returns `None` if overflows.
///
/// # Panics
/// Panics if self is zero and exp is less than zero.
pub fn power_int(&self, exp: i32, rscale: i32) -> Option<Self> {
debug_assert!(!self.is_nan());
// Handle some common special cases, as well as corner cases
match exp {
0 => {
// While 0 ^ 0 can be either 1 or indeterminate (error), we treat
// it as 1 because most programming languages do this. SQL:2003
// also requires a return value of 1.
// https://en.wikipedia.org/wiki/Exponentiation#Zero_to_the_zero_power
let mut result = ONE.clone();
result.dscale = rscale; // no need to round
return Some(result);
}
1 => {
let mut result = self.clone();
result.round_common(rscale);
return Some(result);
}
-1 => {
let result = ONE
.div_common(self, rscale, true)
.expect(DIVIDE_BY_ZERO_MSG);
return Some(result);
}
2 => {
let result = self.mul_common(self, rscale);
return Some(result);
}
_ => (),
}
// Handle the special case where the base is zero
if self.ndigits == 0 {
assert!(exp >= 0, DIVIDE_BY_ZERO_MSG);
return Some(Self::zero_scaled(rscale));
}
// The general case repeatedly multiplies base according to the bit
// pattern of exp.
//
// First we need to estimate the weight of the result so that we know how
// many significant digits are needed.
let digits = self.digits();
let mut f = digits[0] as f64;
let mut p = self.weight * DEC_DIGITS;
for (i, &digit) in digits.iter().enumerate().skip(1) {
if (i * DEC_DIGITS as usize) < 16 {
break;
}
f = f * NBASE as f64 + digit as f64;
p -= DEC_DIGITS;
}
// We have base ~= f * 10^p
// so log10(result) = log10(base^exp) ~= exp * (log10(f) + p)
f = exp as f64 * (f.log10() + p as f64);
// Apply crude overflow/underflow tests so we can exit early if the result
// certainly will overflow/underflow.
if f > 3.0 * i16::max_value() as f64 * DEC_DIGITS as f64 {
return None;
}
if f + 1.0 < (-rscale) as f64 || f + 1.0 < (-NUMERIC_MAX_DISPLAY_SCALE) as f64 {
return Some(Self::zero_scaled(rscale));
}
// Approximate number of significant digits in the result. Note that the
// underflow test above means that this is necessarily >= 0.
let mut sig_digits = 1 + rscale + f as i32;
// The multiplications to produce the result may introduce an error of up
// to around log10(abs(exp)) digits, so work with this many extra digits
// of precision (plus a few more for good measure).
sig_digits += (exp.abs() as f64).ln() as i32 + 8;
// Now we can proceed with the multiplications.
let mut neg = exp < 0;
let mut mask = exp.abs();
let mut base_prod = self.clone();
let mut result = if mask & 1 != 0 {
self.clone()
} else {
ONE.clone()
};
loop {
mask >>= 1;
if mask <= 0 {
break;
}
// Do the multiplications using rscales large enough to hold the
// results to the required number of significant digits, but don't
// waste time by exceeding the scales of the numbers themselves.
let local_rscale = (sig_digits - 2 * base_prod.weight * DEC_DIGITS)
.min(2 * base_prod.dscale)
.max(NUMERIC_MIN_DISPLAY_SCALE);
base_prod = base_prod.mul_common(&base_prod, local_rscale);
if mask & 1 != 0 {
let local_rscale = (sig_digits - (base_prod.weight + result.weight) * DEC_DIGITS)
.min(base_prod.dscale + result.dscale)
.max(NUMERIC_MIN_DISPLAY_SCALE);
result = base_prod.mul_common(&result, local_rscale);
}
// When abs(base) > 1, the number of digits to the left of the decimal
// point in base_prod doubles at each iteration, so if exp is large we
// could easily spend large amounts of time and memory space doing the
// multiplications. But once the weight exceeds what will fit in
// int16, the final result is guaranteed to overflow (or underflow, if
// exp < 0), so we can give up before wasting too many cycles.
if base_prod.weight > i16::max_value() as i32 || result.weight > i16::max_value() as i32
{
// overflow, unless neg, in which case result should be 0
if !neg {
return None;
}
result = ZERO.clone();
neg = false;
break;
}
}
// Compensate for input sign, and round to requested rscale
if neg {
result = ONE
.div_fast_common(&result, rscale, true)
.expect(DIVIDE_BY_ZERO_MSG);
} else {
result.round_common(rscale);
}
Some(result)
}
/// Compute the natural log of `self`
pub fn ln_common(&self, rscale: i32) -> Self {
debug_assert!(!self.is_nan());
debug_assert!(self.cmp_common(&ZERO) > 0);
let mut x = self.clone();
let mut fact = TWO.clone();
// Reduce input into range 0.9 < x < 1.1 with repeated sqrt() operations.
//
// The final logarithm will have up to around rscale+6 significant digits.
// Each sqrt() will roughly halve the weight of x, so adjust the local
// rscale as we work so that we keep this many significant digits at each
// step (plus a few more for good measure).
while x.cmp_common(&ZERO_POINT_NINE) <= 0 {
let mut local_rscale = rscale - x.weight * DEC_DIGITS / 2 + 8;
local_rscale = local_rscale.max(NUMERIC_MIN_DISPLAY_SCALE);
x = x.sqrt_common(local_rscale);
fact = fact.mul_common(&TWO, 0);
}
while x.cmp_common(&ONE_POINT_ONE) >= 0 {
let mut local_rscale = rscale - x.weight * DEC_DIGITS / 2 + 8;
local_rscale = local_rscale.max(NUMERIC_MIN_DISPLAY_SCALE);
x = x.sqrt_common(local_rscale);
fact = fact.mul_common(&TWO, 0);
}
// We use the Taylor series for 0.5 * ln((1+z)/(1-z)),
//
// z + z^3/3 + z^5/5 + ...
//
// where z = (x-1)/(x+1) is in the range (approximately) -0.053 .. 0.048
// due to the above range-reduction of x.
//
// The convergence of this is not as fast as one would like, but is
// tolerable given that z is small.
let local_rscale = rscale + 8;
let mut result = x.sub_common(&ONE);
let mut elem = x.add_common(&ONE);
result = result
.div_fast_common(&elem, local_rscale, true)
.expect(DIVIDE_BY_ZERO_MSG);
let mut xx = result.clone();
x = result.mul_common(&result, local_rscale);
let mut ni = ONE.clone();
loop {
ni = ni.add_common(&TWO);
xx = xx.mul_common(&x, local_rscale);
elem = xx
.div_fast_common(&ni, local_rscale, true)
.expect(DIVIDE_BY_ZERO_MSG);
if elem.ndigits == 0 {
break;
}
result = result.add_common(&elem);
if elem.weight < (result.weight - local_rscale * 2 / DEC_DIGITS) {
break;
}
}
// Compensate for argument range reduction, round to requested rscale
result = result.mul_common(&fact, rscale);
result
}
/// Estimate the dweight of the most significant decimal digit of the natural
/// logarithm of a number.
///
/// Essentially, we're approximating `log10(abs(ln(self)))`. This is used to
/// determine the appropriate rscale when computing natural logarithms.
pub fn estimate_ln_dweight(&self) -> i32 {
debug_assert!(!self.is_nan());
let ln_dweight: i32;
if self.cmp_common(&ZERO_POINT_NINE) >= 0 && self.cmp_common(&ONE_POINT_ONE) <= 0 {
// 0.9 <= self <= 1.1
//
// ln(self) has a negative weight (possibly very large). To get a
// reasonably accurate result, estimate it using ln(1+x) ~= x.
let x = self.sub_common(&ONE);
if x.ndigits > 0 {
// Use weight of most significant decimal digit of x
ln_dweight = x.weight * DEC_DIGITS + (x.digits()[0] as f64).log10() as i32;
} else {
// x = 0. Since ln(1) = 0 exactly, we don't need extra digits
ln_dweight = 0;
}
} else {
// Estimate the logarithm using the first couple of digits from the
// input number. This will give an accurate result whenever the input
// is not too close to 1.
if self.ndigits > 0 {
let d = self.digits();
let mut digits = d[0] as i32;
let mut dweight = self.weight * DEC_DIGITS;
if self.ndigits > 1 {
digits = digits * NBASE + d[1] as i32;
dweight -= DEC_DIGITS;
}
// We have self ~= digits * 10^dweight
// so ln(self) ~= ln(digits) + dweight * ln(10)
let ln_var = (digits as f64).ln() + dweight as f64 * LN_10;
ln_dweight = ln_var.abs().log10() as i32;
} else {
ln_dweight = 0;
}
}
ln_dweight
}
/// Compute the logarithm of `self` in a given base.
///
/// Note: this routine chooses dscale of the result.
pub fn log_common(&self, base: &Self) -> Self {
debug_assert!(!self.is_nan());
debug_assert!(!base.is_nan());
debug_assert!(self.cmp_common(&ZERO) > 0);
// Estimated dweights of ln(base), ln(self) and the final result
let ln_base_dweight = base.estimate_ln_dweight();
let ln_num_dweight = self.estimate_ln_dweight();
let result_dweight = ln_num_dweight - ln_base_dweight;
// Select the scale of the result so that it will have at least
// NUMERIC_MIN_SIG_DIGITS significant digits and is not less than either
// input's display scale.
let rscale = (NUMERIC_MIN_SIG_DIGITS - result_dweight)
.max(base.dscale)
.max(self.dscale)
.max(NUMERIC_MIN_DISPLAY_SCALE)
.min(NUMERIC_MAX_DISPLAY_SCALE);
// Set the scales for ln(base) and ln(num) so that they each have more
// significant digits than the final result.
let ln_base_rscale =
(rscale + result_dweight - ln_base_dweight + 8).max(NUMERIC_MIN_DISPLAY_SCALE);
let ln_num_rscale =
(rscale + result_dweight - ln_num_dweight + 8).max(NUMERIC_MIN_DISPLAY_SCALE);
// Form natural logarithms
let ln_base = base.ln_common(ln_base_rscale);
let ln_num = self.ln_common(ln_num_rscale);
// Divide and round to the required scale
ln_num
.div_fast_common(&ln_base, rscale, true)
.expect(DIVIDE_BY_ZERO_MSG)
}
/// Raise e to the power of x, computed to rscale fractional digits.
///
/// Returns `None` if overflows.
pub fn exp_common(&self, rscale: i32) -> Option<Self> {
debug_assert!(!self.is_nan());
let mut x = self.clone();
// Estimate the dweight of the result using floating point arithmetic, so
// that we can choose an appropriate local rscale for the calculation.
let mut val: f64 = TryFrom::try_from(&x).unwrap();
// Guard against overflow
// If you change this limit, see also power_common()'s limit
if val.abs() >= NUMERIC_MAX_RESULT_SCALE as f64 * 3.0 {
// value overflows numeric format
return None;
}
// decimal weight = log10(e^x) = x * log10(e)
let dweight = (val * LOG10_E) as i32;
let mut ndiv2: i32;
// Reduce x to the range -0.01 <= x <= 0.01 (approximately) by dividing by
// 2^n, to improve the convergence rate of the Taylor series.
if val.abs() > 0.01 {
let mut tmp = TWO.clone();
ndiv2 = 1;
val /= 2.0;
while val.abs() > 0.01 {
ndiv2 += 1;
val /= 2.0;
tmp = tmp.add_common(&tmp);
}
let local_rscale = x.dscale + ndiv2;
x = x
.div_fast_common(&tmp, local_rscale, true)
.expect(DIVIDE_BY_ZERO_MSG);
} else {
ndiv2 = 0;
}
// Set the scale for the Taylor series expansion. The final result has
// (dweight + rscale + 1) significant digits. In addition, we have to
// raise the Taylor series result to the power 2^ndiv2, which introduces
// an error of up to around log10(2^ndiv2) digits, so work with this many
// extra digits of precision (plus a few more for good measure).
let mut sig_digits = 1 + dweight + rscale + (ndiv2 as f64 * LOG10_2) as i32;
sig_digits = sig_digits.max(0) + 8;
let local_rscale = sig_digits - 1;
// Use the Taylor series
//
// exp(x) = 1 + x + x^2/2! + x^3/3! + ...
//
// Given the limited range of x, this should converge reasonably quickly.
// We run the series until the terms fall below the local_rscale limit.
let mut result = ONE.add_common(&x);
let mut elem = x.mul_common(&x, local_rscale);
let mut ni = TWO.clone();
elem = elem
.div_fast_common(&ni, local_rscale, true)
.expect(DIVIDE_BY_ZERO_MSG);
while elem.ndigits != 0 {
result = result.add_common(&elem);
elem = elem.mul_common(&x, local_rscale);
ni = ni.add_common(&ONE);
elem = elem
.div_fast_common(&ni, local_rscale, true)
.expect(DIVIDE_BY_ZERO_MSG);
}
// Compensate for the argument range reduction. Since the weight of the
// result doubles with each multiplication, we can reduce the local rscale
// as we proceed.
for _ in 1..=ndiv2 {
let mut local_rscale = sig_digits - result.weight * 2 * DEC_DIGITS;
local_rscale = local_rscale.max(NUMERIC_MIN_DISPLAY_SCALE);
result = result.mul_common(&result, local_rscale);
}
// Round to requested rscale
result.round_common(rscale);
Some(result)
}
/// Raise `self` to the power of `exp`
///
/// Returns `None` if overflows.
///
/// Note: this routine chooses dscale of the result.
///
/// # Panics
/// Panics if self is zero and exp is less than zero.
pub fn power_common(&self, exp: &Self) -> Option<Self> {
debug_assert!(!self.is_nan());
debug_assert!(!exp.is_nan());
// If exp can be represented as an integer, use power_int()
if exp.ndigits == 0 || exp.ndigits <= exp.weight + 1 {
// exact integer, but does it fit in i32?
if let Ok(exp_val) = TryInto::<i32>::try_into(exp) {
// Okay, select rscale
let rscale = NUMERIC_MIN_SIG_DIGITS
.max(self.dscale)
.max(NUMERIC_MIN_DISPLAY_SCALE)
.min(NUMERIC_MAX_DISPLAY_SCALE);
let result = self.power_int(exp_val, rscale);
return result;
}
}
// This avoids log(0) for cases of 0 raised to a non-integer. 0 ^ 0 is
// handled by power_int().
if self.cmp_common(&ZERO) == 0 {
return Some(Self::zero_scaled(NUMERIC_MIN_SIG_DIGITS));
}
// Decide on the scale for the ln() calculation. For this we need an
// estimate of the weight of the result, which we obtain by doing an
// initial low-precision calculation of exp * ln(base).
//
// We want result = e ^ (exp * ln(base))
// so result dweight = log10(result) = exp * ln(base) * log10(e)
//
// We also perform a crude overflow test here so that we can exit early if
// the full-precision result is sure to overflow, and to guard against
// integer overflow when determining the scale for the real calculation.
// exp_common() supports inputs up to NUMERIC_MAX_RESULT_SCALE * 3, so the
// result will overflow if exp * ln(base) >= NUMERIC_MAX_RESULT_SCALE * 3.
// Since the values here are only approximations, we apply a small fuzz
// factor to this overflow test and let exp_common() determine the exact
// overflow threshold so that it is consistent for all inputs.
let ln_dweight = self.estimate_ln_dweight();
let local_rscale = (8 - ln_dweight)
.max(NUMERIC_MIN_DISPLAY_SCALE)
.min(NUMERIC_MAX_DISPLAY_SCALE);
let mut ln_base = self.ln_common(local_rscale);
let mut ln_num = ln_base.mul_common(exp, local_rscale);
let mut val: f64 = TryFrom::try_from(&ln_num).unwrap();
// initial overflow test with fuzz factor
if val.abs() > NUMERIC_MAX_RESULT_SCALE as f64 * 3.01 {
// value overflows numeric format
return None;
}
val *= LOG10_E; // approximate decimal result weight
// choose the result scale
let rscale = (NUMERIC_MIN_SIG_DIGITS - val as i32)
.max(self.dscale)
.max(exp.dscale)
.max(NUMERIC_MIN_DISPLAY_SCALE)
.min(NUMERIC_MAX_DISPLAY_SCALE);
// set the scale for the real exp * ln(base) calculation
let local_rscale = (rscale + val as i32 - ln_dweight + 8).max(NUMERIC_MIN_DISPLAY_SCALE);
// and do the real calculation
ln_base = self.ln_common(local_rscale);
ln_num = ln_base.mul_common(exp, local_rscale);
ln_num.exp_common(rscale)
}
/// Convert `self` to text representation.
/// `self` is displayed to the number of digits indicated by its dscale.
pub fn write<W: fmt::Write>(&self, f: &mut W) -> Result<(), fmt::Error> {
if self.is_nan() {
return f.write_str("NaN");
}
// Output a dash for negative values.
if self.sign == NUMERIC_NEG {
f.write_char('-')?;
}
// Output all digits before the decimal point.
if self.weight < 0 {
f.write_char('0')?;
} else {
let digits = self.digits();
debug_assert_eq!(digits.len(), self.ndigits as usize);
for d in 0..=self.weight {
let dig = if d < self.ndigits {
digits[d as usize]
} else {
0
};
debug_assert!(dig >= 0);
// In the first digit, suppress extra leading decimal zeroes.
if d > 0 {
write!(f, "{:>0width$}", dig, width = DEC_DIGITS as usize)?;
} else {
write!(f, "{}", dig)?;
}
}
}
// If requested, output a decimal point and all the digits that follow it.
if self.dscale > 0 {
f.write_char('.')?;
let digits = self.digits();
debug_assert_eq!(digits.len(), self.ndigits as usize);
let mut d = self.weight + 1;
for scale in (0..self.dscale).step_by(DEC_DIGITS as usize) {
let dig = if d >= 0 && d < self.ndigits {
digits[d as usize]
} else {
0
};
if scale + DEC_DIGITS <= self.dscale {
write!(f, "{:>0width$}", dig, width = DEC_DIGITS as usize)?;
} else {
// truncate the last digit
let width = (self.dscale - scale) as usize;
let dig = (0..DEC_DIGITS as usize - width).fold(dig, |acc, _| acc / 10);
write!(f, "{:>0width$}", dig, width = width)?;
}
d += 1;
}
}
Ok(())
}
/// Convert `self` to a normalised scientific notation text representation.
///
/// This notation has the general form a * 10^b, where a is known as the
/// "significand" and b is known as the "exponent".
///
/// Because we can't do superscript in ASCII (and because we want to copy
/// printf's behaviour) we display the exponent using E notation, with a
/// minimum of two exponent digits.
///
/// For example, the value 1234 could be output as 1.2e+03.
///
/// We assume that the exponent can fit into an int32.
///
/// `rscale` is the number of decimal digits desired after the decimal point in
/// the output, negative values will be treated as meaning zero.
///
/// `lower_exp` indicates use 'e' if true or else use 'E'.
pub fn write_sci<W: fmt::Write>(
&self,
f: &mut W,
rscale: i32,
lower_exp: bool,
) -> Result<(), fmt::Error> {
if self.is_nan() {
return write!(f, "NaN");
}
let rscale = if rscale < 0 { 0 } else { rscale };
// Determine the exponent of this number in normalised form.
//
// This is the exponent required to represent the number with only one
// significant digit before the decimal place.
let exponent = if self.ndigits > 0 {
let mut exp = (self.weight + 1) * DEC_DIGITS;
// Compensate for leading decimal zeroes in the first numeric digit by
// decrementing the exponent.
exp -= DEC_DIGITS - (self.digits[0] as f64).log10() as i32;
exp
} else {
// If has no digits, then it must be zero.
//
// Zero doesn't technically have a meaningful exponent in normalised
// notation, but we just display the exponent as zero for consistency
// of output.
0
};
// The denominator is set to 10 raised to the power of the exponent.
//
// We then divide var by the denominator to get the significand, rounding
// to rscale decimal digits in the process.
let denom_scale = if exponent < 0 { -exponent } else { 0 };
let denominator = TEN
.power_int(exponent, denom_scale)
.expect("attempt to multiply with overflow");
let significand = self
.div_common(&denominator, rscale, true)
.expect(DIVIDE_BY_ZERO_MSG);
if lower_exp {
write!(f, "{}e{:<+03}", significand, exponent)
} else {
write!(f, "{}E{:<+03}", significand, exponent)
}
}
/// Returns the appropriate result scale for scientific notation representation.
pub fn select_sci_scale(&self) -> i32 {
// 1 => (1, 0)
// 10 => (1, 1)
// 11 => (2, 0)
// 100 => (1, 2)
// 101 => (3, 0)
// 1010 => (3, 1)
fn count_zeros(digit: NumericDigit) -> (i32, i32) {
let mut val = digit;
let mut n = 0;
let mut zero = 0;
for _ in 0..DEC_DIGITS {
let d = val % 10;
val /= 10;
if d == 0 && n == 0 {
// all previous d are zeros.
zero += 1;
} else {
n += 1;
}
if val == 0 {
break;
}
}
(n, zero)
}
let digits = self.digits();
// find first non-zero digit from front to end
let (i, digit) = match digits.iter().enumerate().find(|(_, &d)| d != 0) {
Some((i, &digit)) => (i, digit),
None => {
// all digits are 0
return 0;
}
};
// find first non-zero digit from end to front
let (ri, rdigit) = match digits.iter().enumerate().rfind(|(_, &d)| d != 0) {
Some((ri, &rdigit)) => (ri, rdigit),
None => {
// all digits are 0, actually unreachable!
return 0;
}
};
debug_assert!(i <= ri);
if i == ri {
// only one digit
let (n, _) = count_zeros(digit);
return n - 1;
}
let (n, zero) = count_zeros(digit);
let (_, rzero) = count_zeros(rdigit);
let front = n + zero;
let end = DEC_DIGITS - rzero;
front + end + (ri - i - 1) as i32 * DEC_DIGITS - 1
}
/// Make `self` to be a result numeric.
/// We assume that `self` is not overflowed.
#[inline]
pub fn make_result_no_overflow(mut self) -> NumericBuf {
debug_assert!(!self.is_nan());
debug_assert!(
self.weight <= NUMERIC_WEIGHT_MAX as i32
|| self.weight >= NUMERIC_WEIGHT_MIN as i32
|| self.dscale <= NUMERIC_DSCALE_MAX as i32
|| self.dscale >= 0
);
self.strip();
self.into_numeric_buf()
}
/// Make `self` to be a result numeric.
/// Returns `None` if overflows.
#[inline]
pub fn make_result(self) -> Option<NumericBuf> {
debug_assert!(!self.is_nan());
debug_assert!(self.dscale >= 0);
if self.weight > NUMERIC_WEIGHT_MAX as i32
|| self.weight < NUMERIC_WEIGHT_MIN as i32
|| self.dscale > NUMERIC_DSCALE_MAX as i32
{
return None;
}
Some(self.make_result_no_overflow())
}
/// Negate this value.
#[inline]
pub fn negate(&mut self) {
debug_assert!(!self.is_nan());
if self.ndigits > 0 {
if self.is_positive() {
self.sign = NUMERIC_NEG;
} else if self.is_negative() {
self.sign = NUMERIC_POS;
}
}
}
/// Returns a numeric that represents the sign of self.
/// * -1 if `self` is less than 0
/// * 0 if `self` is equal to 0
/// * 1 if `self` is greater than zero
/// * `NaN` if `self` is `NaN`
#[inline]
pub fn signum(&self) -> Self {
debug_assert!(!self.is_nan());
if self.ndigits == 0 {
ZERO.clone()
} else {
let mut result = ONE.clone();
result.sign = self.sign;
result
}
}
/// Increment `self` by one.
#[inline]
pub fn inc(&self) -> Self {
debug_assert!(!self.is_nan());
// Compute the result and return it
self.add_common(&ONE)
}
/// Checked numeric division.
/// Computes `self / other`, returning `None` if `other == 0`.
#[inline]
pub fn checked_div(&self, other: &Self) -> Option<Self> {
debug_assert!(!self.is_nan());
debug_assert!(!other.is_nan());
// Select scale for division result
let rscale = self.select_div_scale(&other);
self.div_common(&other, rscale, true)
}
/// Computes `self / other`, truncating the result to an integer.
///
/// Returns `None` if `other == 0`.
#[inline]
pub fn checked_div_trunc(&self, other: &Self) -> Option<Self> {
debug_assert!(!self.is_nan());
debug_assert!(!other.is_nan());
self.div_common(&other, 0, false)
}
/// Round a value to have `scale` digits after the decimal point.
/// We allow negative `scale`, implying rounding before the decimal
/// point --- Oracle interprets rounding that way.
#[inline]
pub fn round(&mut self, scale: i32) {
debug_assert!(!self.is_nan());
// Limit the scale value to avoid possible overflow in calculations
let rscale = scale
.max(-NUMERIC_MAX_DISPLAY_SCALE)
.min(NUMERIC_MAX_DISPLAY_SCALE);
self.round_common(rscale);
// We don't allow negative output dscale
if rscale < 0 {
self.dscale = 0;
}
}
/// Truncate a value to have `scale` digits after the decimal point.
/// We allow negative `scale`, implying a truncation before the decimal
/// point --- Oracle interprets truncation that way.
#[inline]
pub fn trunc(&mut self, scale: i32) {
debug_assert!(!self.is_nan());
// Limit the scale value to avoid possible overflow in calculations
let rscale = scale
.max(-NUMERIC_MAX_DISPLAY_SCALE)
.min(NUMERIC_MAX_DISPLAY_SCALE);
self.trunc_common(rscale);
// We don't allow negative output dscale
if rscale < 0 {
self.dscale = 0;
}
}
/// Return the smallest integer greater than or equal to the argument.
#[inline]
pub fn ceil(&self) -> Self {
debug_assert!(!self.is_nan());
self.ceil_common()
}
/// Return the largest integer equal to or less than the argument.
#[inline]
pub fn floor(&self) -> Self {
debug_assert!(!self.is_nan());
self.floor_common()
}
/// Compute the absolute value of `self`.
#[inline]
pub fn abs(&mut self) {
debug_assert!(!self.is_nan());
if self.is_negative() {
self.sign = NUMERIC_POS;
}
}
/// Compute the square root of a numeric.
#[inline]
pub fn sqrt(&self) -> Self {
debug_assert!(!self.is_negative());
debug_assert!(!self.is_nan());
// Determine the result scale.
// We choose a scale to give at least NUMERIC_MIN_SIG_DIGITS significant digits;
// but in any case not less than the input's dscale.
// Assume the input was normalized, so arg.weight is accurate
let sweight = (self.weight + 1) * DEC_DIGITS / 2 - 1;
let rscale = (NUMERIC_MIN_SIG_DIGITS - sweight)
.max(self.dscale)
.max(NUMERIC_MIN_DISPLAY_SCALE)
.min(NUMERIC_MAX_DISPLAY_SCALE);
self.sqrt_common(rscale)
}
/// Compute the natural logarithm of `self`.
///
/// # Panics
/// Panics if `self <= 0`.
#[inline]
pub fn ln(&self) -> Self {
debug_assert!(!self.is_nan());
let cmp = self.cmp_common(&ZERO);
assert_ne!(cmp, 0, "cannot take logarithm of zero");
assert!(cmp > 0, "cannot take logarithm of a negative number");
// Estimated dweight of logarithm
let ln_dweight = self.estimate_ln_dweight();
let rscale = (NUMERIC_MIN_SIG_DIGITS - ln_dweight)
.max(self.dscale)
.max(NUMERIC_MIN_DISPLAY_SCALE)
.min(NUMERIC_MAX_DISPLAY_SCALE);
self.ln_common(rscale)
}
/// Compute the logarithm of `self` in a given base.
///
/// # Panics
/// Panics if `self <= 0` or `base <= 0`.
#[inline]
pub fn log(&self, base: &Self) -> Self {
debug_assert!(!self.is_nan());
debug_assert!(!base.is_nan());
let cmp = self.cmp_common(&ZERO);
assert_ne!(cmp, 0, "cannot take logarithm of zero");
assert!(cmp > 0, "cannot take logarithm of a negative number");
let cmp = base.cmp_common(&ZERO);
assert_ne!(cmp, 0, "cannot take logarithm of zero");
assert!(cmp > 0, "cannot take logarithm of a negative number");
// Call log_common() to compute and return the result;
// note it handles scale selection itself.
self.log_common(&base)
}
/// Compute the base 2 logarithm of `self`.
///
/// # Panics
/// Panics if `self <= 0`.
#[inline]
pub fn log2(&self) -> Self {
self.log(&TWO)
}
/// Compute the base 10 logarithm of `self`.
///
/// # Panics
/// Panics if `self <= 0`.
#[inline]
pub fn log10(&self) -> Self {
self.log(&TEN)
}
/// Raise e to the power of `self` (`e^(self)`).
///
/// Returns `None` if overflows.
///
#[inline]
pub fn exp(&self) -> Option<Self> {
debug_assert!(!self.is_nan());
// Determine the result scale. We choose a scale
// to give at least NUMERIC_MIN_SIG_DIGITS significant digits;
// but in any case not less than the input's dscale.
let mut val: f64 = TryFrom::try_from(self).unwrap();
// log10(result) = num * log10(e), so this is approximately the decimal
// weight of the result:
val *= LOG10_E;
// limit to something that won't cause integer overflow
val = val
.max(-NUMERIC_MAX_RESULT_SCALE as f64)
.min(NUMERIC_MAX_RESULT_SCALE as f64);
let rscale = (NUMERIC_MIN_SIG_DIGITS - val as i32)
.max(self.dscale)
.max(NUMERIC_MIN_DISPLAY_SCALE)
.min(NUMERIC_MAX_DISPLAY_SCALE);
// Let exp_common() do the calculation and return the result.
self.exp_common(rscale)
}
/// Raise `self` to the power of `exp`.
///
/// Returns `None` if overflows.
///
/// # Panics
/// if arguments are invalid:
/// - `self` is zero and `exp` is less than zero
/// - `self` is less than zero and `exp` is not a integer.
#[inline]
pub fn pow(&self, exp: &Self) -> Option<Self> {
// Handle NaN cases. We follow the POSIX spec for pow(3), which says that
// NaN ^ 0 = 1, and 1 ^ NaN = 1, while all other cases with NaN inputs
// yield NaN (with no error).
if self.is_nan() {
if !exp.is_nan() && exp.cmp_common(&ZERO) == 0 {
return Some(ONE.clone());
}
return Some(Self::nan());
} else if exp.is_nan() {
if self.cmp_common(&ONE) == 0 {
return Some(ONE.clone());
}
return Some(Self::nan());
}
if self.cmp_common(&ZERO) == 0 && exp.cmp_common(&ZERO) < 0 {
panic!("zero raised to a negative power is undefined")
}
let mut exp_trunc = exp.clone();
exp_trunc.trunc_common(0);
if self.cmp_common(&ZERO) < 0 && exp.cmp_common(&exp_trunc) != 0 {
panic!("a negative number raised to a non-integer power yields a complex result")
}
// Call power_common() to compute and return the result; note it handles
// scale selection itself.
self.power_common(&exp)
}
/// Do bounds checking and rounding according to `typmod`.
///
/// Returns true if overflows.
///
/// Notes that no matter whether overflows, `self` will be rounded.
pub fn apply_typmod(&mut self, typmod: Typmod) -> bool {
// Do nothing if we have a default typmod (-1)
if typmod.value() < VAR_HEADER_SIZE {
return false;
}
let (precision, scale) = typmod.extract();
let max_digits = precision - scale;
// Round to target scale (and set self.dscale)
self.round_common(scale);
// Check for overflow - note we can't do this before rounding, because
// rounding could raise the weight. Also note that the self's weight could
// be inflated by leading zeroes, which will be stripped before storage
// but perhaps might not have been yet. In any case, we must recognize a
// true zero, whose weight doesn't mean anything.
let mut ddigits = (self.weight + 1) * DEC_DIGITS;
if ddigits > max_digits {
// Determine true weight; and check for all-zero result
for &dig in self.digits().iter() {
if dig != 0 {
// Adjust for any high-order decimal zero digits
debug_assert_eq!(DEC_DIGITS, 4);
if dig < 10 {
ddigits -= 3;
} else if dig < 100 {
ddigits -= 2;
} else if dig < 1000 {
ddigits -= 1;
}
if ddigits > max_digits {
return true;
}
break;
}
ddigits -= DEC_DIGITS;
}
}
false
}
/// Returns a normalized string, suppressing insignificant
/// trailing zeroes and then any trailing decimal point.
///
/// The intent of this is to produce strings that are equal
/// if and only if the input numeric values compare equal.
pub fn normalize(&self) -> String {
if self.is_nan() {
return String::from("NaN");
}
let mut s = self.to_string();
// If there's no decimal point, there's certainly nothing to remove.
if s.find('.').is_some() {
// Back up over trailing fractional zeroes. Since there is a decimal
// point, this loop will terminate safely.
let mut new_len = s.len();
let bytes = s.as_bytes();
let count_zeros = bytes.iter().rev().take_while(|i| **i == b'0').count();
new_len -= count_zeros;
// We want to get rid of the decimal point too, if it's now last.
if bytes[new_len - 1] == b'.' {
new_len -= 1;
}
// Delete whatever we backed up over.
s.truncate(new_len);
}
s
}
/// Compute factorial.
///
/// Returns `None` if overflows.
pub fn factorial(num: i64) -> Option<Self> {
if num <= 1 {
return Some(ONE.clone());
}
// Fail immediately if the result would overflow
if num > 32177 {
// value overflows numeric format
return None;
}
let mut result: NumericVar = From::from(num);
for n in (2..num).rev() {
let fact = From::from(n);
result = result.mul_common(&fact, 0);
}
Some(result)
}
}
impl<'a> fmt::Display for NumericVar<'a> {
#[inline]
fn fmt(&self, f: &mut fmt::Formatter) -> Result<(), fmt::Error> {
self.write(f)
}
}