|
| 1 | +// Copyright 2015 The Rust Project Developers. See the COPYRIGHT |
| 2 | +// file at the top-level directory of this distribution and at |
| 3 | +// http://rust-lang.org/COPYRIGHT. |
| 4 | +// |
| 5 | +// Licensed under the Apache License, Version 2.0 <LICENSE-APACHE or |
| 6 | +// http://www.apache.org/licenses/LICENSE-2.0> or the MIT license |
| 7 | +// <LICENSE-MIT or http://opensource.org/licenses/MIT>, at your |
| 8 | +// option. This file may not be copied, modified, or distributed |
| 9 | +// except according to those terms. |
| 10 | + |
| 11 | +//! Custom arbitrary-precision number (bignum) implementation. |
| 12 | +//! |
| 13 | +//! This is designed to avoid the heap allocation at expense of stack memory. |
| 14 | +//! The most used bignum type, `Big32x36`, is limited by 32 × 36 = 1,152 bits |
| 15 | +//! and will take at most 152 bytes of stack memory. This is (barely) enough |
| 16 | +//! for handling all possible finite `f64` values. |
| 17 | +//! |
| 18 | +//! In principle it is possible to have multiple bignum types for different |
| 19 | +//! inputs, but we don't do so to avoid the code bloat. Each bignum is still |
| 20 | +//! tracked for the actual usages, so it normally doesn't matter. |
| 21 | +
|
| 22 | +#![macro_use] |
| 23 | + |
| 24 | +use prelude::*; |
| 25 | +use mem; |
| 26 | +use intrinsics; |
| 27 | + |
| 28 | +/// Arithmetic operations required by bignums. |
| 29 | +pub trait FullOps { |
| 30 | + /// Returns `(carry', v')` such that `carry' * 2^W + v' = self + other + carry`, |
| 31 | + /// where `W` is the number of bits in `Self`. |
| 32 | + fn full_add(self, other: Self, carry: bool) -> (bool /*carry*/, Self); |
| 33 | + |
| 34 | + /// Returns `(carry', v')` such that `carry' * 2^W + v' = self * other + carry`, |
| 35 | + /// where `W` is the number of bits in `Self`. |
| 36 | + fn full_mul(self, other: Self, carry: Self) -> (Self /*carry*/, Self); |
| 37 | + |
| 38 | + /// Returns `(carry', v')` such that `carry' * 2^W + v' = self * other + other2 + carry`, |
| 39 | + /// where `W` is the number of bits in `Self`. |
| 40 | + fn full_mul_add(self, other: Self, other2: Self, carry: Self) -> (Self /*carry*/, Self); |
| 41 | + |
| 42 | + /// Returns `(quo, rem)` such that `borrow * 2^W + self = quo * other + rem` |
| 43 | + /// and `0 <= rem < other`, where `W` is the number of bits in `Self`. |
| 44 | + fn full_div_rem(self, other: Self, borrow: Self) -> (Self /*quotient*/, Self /*remainder*/); |
| 45 | +} |
| 46 | + |
| 47 | +macro_rules! impl_full_ops { |
| 48 | + ($($ty:ty: add($addfn:path), mul/div($bigty:ident);)*) => ( |
| 49 | + $( |
| 50 | + impl FullOps for $ty { |
| 51 | + fn full_add(self, other: $ty, carry: bool) -> (bool, $ty) { |
| 52 | + // this cannot overflow, the output is between 0 and 2*2^nbits - 1 |
| 53 | + // FIXME will LLVM optimize this into ADC or similar??? |
| 54 | + let (v, carry1) = unsafe { $addfn(self, other) }; |
| 55 | + let (v, carry2) = unsafe { $addfn(v, if carry {1} else {0}) }; |
| 56 | + (carry1 || carry2, v) |
| 57 | + } |
| 58 | + |
| 59 | + fn full_mul(self, other: $ty, carry: $ty) -> ($ty, $ty) { |
| 60 | + // this cannot overflow, the output is between 0 and 2^nbits * (2^nbits - 1) |
| 61 | + let nbits = mem::size_of::<$ty>() * 8; |
| 62 | + let v = (self as $bigty) * (other as $bigty) + (carry as $bigty); |
| 63 | + ((v >> nbits) as $ty, v as $ty) |
| 64 | + } |
| 65 | + |
| 66 | + fn full_mul_add(self, other: $ty, other2: $ty, carry: $ty) -> ($ty, $ty) { |
| 67 | + // this cannot overflow, the output is between 0 and 2^(2*nbits) - 1 |
| 68 | + let nbits = mem::size_of::<$ty>() * 8; |
| 69 | + let v = (self as $bigty) * (other as $bigty) + (other2 as $bigty) + |
| 70 | + (carry as $bigty); |
| 71 | + ((v >> nbits) as $ty, v as $ty) |
| 72 | + } |
| 73 | + |
| 74 | + fn full_div_rem(self, other: $ty, borrow: $ty) -> ($ty, $ty) { |
| 75 | + debug_assert!(borrow < other); |
| 76 | + // this cannot overflow, the dividend is between 0 and other * 2^nbits - 1 |
| 77 | + let nbits = mem::size_of::<$ty>() * 8; |
| 78 | + let lhs = ((borrow as $bigty) << nbits) | (self as $bigty); |
| 79 | + let rhs = other as $bigty; |
| 80 | + ((lhs / rhs) as $ty, (lhs % rhs) as $ty) |
| 81 | + } |
| 82 | + } |
| 83 | + )* |
| 84 | + ) |
| 85 | +} |
| 86 | + |
| 87 | +impl_full_ops! { |
| 88 | + u8: add(intrinsics::u8_add_with_overflow), mul/div(u16); |
| 89 | + u16: add(intrinsics::u16_add_with_overflow), mul/div(u32); |
| 90 | + u32: add(intrinsics::u32_add_with_overflow), mul/div(u64); |
| 91 | +// u64: add(intrinsics::u64_add_with_overflow), mul/div(u128); // damn! |
| 92 | +} |
| 93 | + |
| 94 | +macro_rules! define_bignum { |
| 95 | + ($name:ident: type=$ty:ty, n=$n:expr) => ( |
| 96 | + /// Stack-allocated arbitrary-precision (up to certain limit) integer. |
| 97 | + /// |
| 98 | + /// This is backed by an fixed-size array of given type ("digit"). |
| 99 | + /// While the array is not very large (normally some hundred bytes), |
| 100 | + /// copying it recklessly may result in the performance hit. |
| 101 | + /// Thus this is intentionally not `Copy`. |
| 102 | + /// |
| 103 | + /// All operations available to bignums panic in the case of over/underflows. |
| 104 | + /// The caller is responsible to use large enough bignum types. |
| 105 | + pub struct $name { |
| 106 | + /// One plus the offset to the maximum "digit" in the use. |
| 107 | + /// This does not decrease, so be aware of the computation order. |
| 108 | + /// `base[size..]` should be zero. |
| 109 | + size: usize, |
| 110 | + /// Digits. `[a, b, c, ...]` represents `a + b*n + c*n^2 + ...`. |
| 111 | + base: [$ty; $n] |
| 112 | + } |
| 113 | + |
| 114 | + impl $name { |
| 115 | + /// Makes a bignum from one digit. |
| 116 | + pub fn from_small(v: $ty) -> $name { |
| 117 | + let mut base = [0; $n]; |
| 118 | + base[0] = v; |
| 119 | + $name { size: 1, base: base } |
| 120 | + } |
| 121 | + |
| 122 | + /// Makes a bignum from `u64` value. |
| 123 | + pub fn from_u64(mut v: u64) -> $name { |
| 124 | + use mem; |
| 125 | + |
| 126 | + let mut base = [0; $n]; |
| 127 | + let mut sz = 0; |
| 128 | + while v > 0 { |
| 129 | + base[sz] = v as $ty; |
| 130 | + v >>= mem::size_of::<$ty>() * 8; |
| 131 | + sz += 1; |
| 132 | + } |
| 133 | + $name { size: sz, base: base } |
| 134 | + } |
| 135 | + |
| 136 | + /// Returns true if the bignum is zero. |
| 137 | + pub fn is_zero(&self) -> bool { |
| 138 | + self.base[..self.size].iter().all(|&v| v == 0) |
| 139 | + } |
| 140 | + |
| 141 | + /// Adds `other` to itself and returns its own mutable reference. |
| 142 | + pub fn add<'a>(&'a mut self, other: &$name) -> &'a mut $name { |
| 143 | + use cmp; |
| 144 | + use num::flt2dec::bignum::FullOps; |
| 145 | + |
| 146 | + let mut sz = cmp::max(self.size, other.size); |
| 147 | + let mut carry = false; |
| 148 | + for (a, b) in self.base[..sz].iter_mut().zip(other.base[..sz].iter()) { |
| 149 | + let (c, v) = (*a).full_add(*b, carry); |
| 150 | + *a = v; |
| 151 | + carry = c; |
| 152 | + } |
| 153 | + if carry { |
| 154 | + self.base[sz] = 1; |
| 155 | + sz += 1; |
| 156 | + } |
| 157 | + self.size = sz; |
| 158 | + self |
| 159 | + } |
| 160 | + |
| 161 | + /// Subtracts `other` from itself and returns its own mutable reference. |
| 162 | + pub fn sub<'a>(&'a mut self, other: &$name) -> &'a mut $name { |
| 163 | + use cmp; |
| 164 | + use num::flt2dec::bignum::FullOps; |
| 165 | + |
| 166 | + let sz = cmp::max(self.size, other.size); |
| 167 | + let mut noborrow = true; |
| 168 | + for (a, b) in self.base[..sz].iter_mut().zip(other.base[..sz].iter()) { |
| 169 | + let (c, v) = (*a).full_add(!*b, noborrow); |
| 170 | + *a = v; |
| 171 | + noborrow = c; |
| 172 | + } |
| 173 | + assert!(noborrow); |
| 174 | + self.size = sz; |
| 175 | + self |
| 176 | + } |
| 177 | + |
| 178 | + /// Multiplies itself by a digit-sized `other` and returns its own |
| 179 | + /// mutable reference. |
| 180 | + pub fn mul_small<'a>(&'a mut self, other: $ty) -> &'a mut $name { |
| 181 | + use num::flt2dec::bignum::FullOps; |
| 182 | + |
| 183 | + let mut sz = self.size; |
| 184 | + let mut carry = 0; |
| 185 | + for a in self.base[..sz].iter_mut() { |
| 186 | + let (c, v) = (*a).full_mul(other, carry); |
| 187 | + *a = v; |
| 188 | + carry = c; |
| 189 | + } |
| 190 | + if carry > 0 { |
| 191 | + self.base[sz] = carry; |
| 192 | + sz += 1; |
| 193 | + } |
| 194 | + self.size = sz; |
| 195 | + self |
| 196 | + } |
| 197 | + |
| 198 | + /// Multiplies itself by `2^bits` and returns its own mutable reference. |
| 199 | + pub fn mul_pow2<'a>(&'a mut self, bits: usize) -> &'a mut $name { |
| 200 | + use mem; |
| 201 | + |
| 202 | + let digitbits = mem::size_of::<$ty>() * 8; |
| 203 | + let digits = bits / digitbits; |
| 204 | + let bits = bits % digitbits; |
| 205 | + |
| 206 | + assert!(digits < $n); |
| 207 | + debug_assert!(self.base[$n-digits..].iter().all(|&v| v == 0)); |
| 208 | + debug_assert!(bits == 0 || (self.base[$n-digits-1] >> (digitbits - bits)) == 0); |
| 209 | + |
| 210 | + // shift by `digits * digitbits` bits |
| 211 | + for i in (0..self.size).rev() { |
| 212 | + self.base[i+digits] = self.base[i]; |
| 213 | + } |
| 214 | + for i in 0..digits { |
| 215 | + self.base[i] = 0; |
| 216 | + } |
| 217 | + |
| 218 | + // shift by `nbits` bits |
| 219 | + let mut sz = self.size + digits; |
| 220 | + if bits > 0 { |
| 221 | + let last = sz; |
| 222 | + let overflow = self.base[last-1] >> (digitbits - bits); |
| 223 | + if overflow > 0 { |
| 224 | + self.base[last] = overflow; |
| 225 | + sz += 1; |
| 226 | + } |
| 227 | + for i in (digits+1..last).rev() { |
| 228 | + self.base[i] = (self.base[i] << bits) | |
| 229 | + (self.base[i-1] >> (digitbits - bits)); |
| 230 | + } |
| 231 | + self.base[digits] <<= bits; |
| 232 | + // self.base[..digits] is zero, no need to shift |
| 233 | + } |
| 234 | + |
| 235 | + self.size = sz; |
| 236 | + self |
| 237 | + } |
| 238 | + |
| 239 | + /// Multiplies itself by a number described by `other[0] + other[1] * n + |
| 240 | + /// other[2] * n^2 + ...` and returns its own mutable reference. |
| 241 | + pub fn mul_digits<'a>(&'a mut self, other: &[$ty]) -> &'a mut $name { |
| 242 | + // the internal routine. works best when aa.len() <= bb.len(). |
| 243 | + fn mul_inner(ret: &mut [$ty; $n], aa: &[$ty], bb: &[$ty]) -> usize { |
| 244 | + use num::flt2dec::bignum::FullOps; |
| 245 | + |
| 246 | + let mut retsz = 0; |
| 247 | + for (i, &a) in aa.iter().enumerate() { |
| 248 | + if a == 0 { continue; } |
| 249 | + let mut sz = bb.len(); |
| 250 | + let mut carry = 0; |
| 251 | + for (j, &b) in bb.iter().enumerate() { |
| 252 | + let (c, v) = a.full_mul_add(b, ret[i + j], carry); |
| 253 | + ret[i + j] = v; |
| 254 | + carry = c; |
| 255 | + } |
| 256 | + if carry > 0 { |
| 257 | + ret[i + sz] = carry; |
| 258 | + sz += 1; |
| 259 | + } |
| 260 | + if retsz < i + sz { |
| 261 | + retsz = i + sz; |
| 262 | + } |
| 263 | + } |
| 264 | + retsz |
| 265 | + } |
| 266 | + |
| 267 | + let mut ret = [0; $n]; |
| 268 | + let retsz = if self.size < other.len() { |
| 269 | + mul_inner(&mut ret, &self.base[..self.size], other) |
| 270 | + } else { |
| 271 | + mul_inner(&mut ret, other, &self.base[..self.size]) |
| 272 | + }; |
| 273 | + self.base = ret; |
| 274 | + self.size = retsz; |
| 275 | + self |
| 276 | + } |
| 277 | + |
| 278 | + /// Divides itself by a digit-sized `other` and returns its own |
| 279 | + /// mutable reference *and* the remainder. |
| 280 | + pub fn div_rem_small<'a>(&'a mut self, other: $ty) -> (&'a mut $name, $ty) { |
| 281 | + use num::flt2dec::bignum::FullOps; |
| 282 | + |
| 283 | + assert!(other > 0); |
| 284 | + |
| 285 | + let sz = self.size; |
| 286 | + let mut borrow = 0; |
| 287 | + for a in self.base[..sz].iter_mut().rev() { |
| 288 | + let (q, r) = (*a).full_div_rem(other, borrow); |
| 289 | + *a = q; |
| 290 | + borrow = r; |
| 291 | + } |
| 292 | + (self, borrow) |
| 293 | + } |
| 294 | + } |
| 295 | + |
| 296 | + impl ::cmp::PartialEq for $name { |
| 297 | + fn eq(&self, other: &$name) -> bool { self.base[..] == other.base[..] } |
| 298 | + } |
| 299 | + |
| 300 | + impl ::cmp::Eq for $name { |
| 301 | + } |
| 302 | + |
| 303 | + impl ::cmp::PartialOrd for $name { |
| 304 | + fn partial_cmp(&self, other: &$name) -> ::option::Option<::cmp::Ordering> { |
| 305 | + ::option::Option::Some(self.cmp(other)) |
| 306 | + } |
| 307 | + } |
| 308 | + |
| 309 | + impl ::cmp::Ord for $name { |
| 310 | + fn cmp(&self, other: &$name) -> ::cmp::Ordering { |
| 311 | + use cmp::max; |
| 312 | + use iter::order; |
| 313 | + |
| 314 | + let sz = max(self.size, other.size); |
| 315 | + let lhs = self.base[..sz].iter().cloned().rev(); |
| 316 | + let rhs = other.base[..sz].iter().cloned().rev(); |
| 317 | + order::cmp(lhs, rhs) |
| 318 | + } |
| 319 | + } |
| 320 | + |
| 321 | + impl ::clone::Clone for $name { |
| 322 | + fn clone(&self) -> $name { |
| 323 | + $name { size: self.size, base: self.base } |
| 324 | + } |
| 325 | + } |
| 326 | + |
| 327 | + impl ::fmt::Debug for $name { |
| 328 | + fn fmt(&self, f: &mut ::fmt::Formatter) -> ::fmt::Result { |
| 329 | + use mem; |
| 330 | + |
| 331 | + let sz = if self.size < 1 {1} else {self.size}; |
| 332 | + let digitlen = mem::size_of::<$ty>() * 2; |
| 333 | + |
| 334 | + try!(write!(f, "{:#x}", self.base[sz-1])); |
| 335 | + for &v in self.base[..sz-1].iter().rev() { |
| 336 | + try!(write!(f, "_{:01$x}", v, digitlen)); |
| 337 | + } |
| 338 | + ::result::Result::Ok(()) |
| 339 | + } |
| 340 | + } |
| 341 | + ) |
| 342 | +} |
| 343 | + |
| 344 | +/// The digit type for `Big32x36`. |
| 345 | +pub type Digit32 = u32; |
| 346 | + |
| 347 | +define_bignum!(Big32x36: type=Digit32, n=36); |
| 348 | + |
| 349 | +// this one is used for testing only. |
| 350 | +#[doc(hidden)] |
| 351 | +pub mod tests { |
| 352 | + use prelude::*; |
| 353 | + define_bignum!(Big8x3: type=u8, n=3); |
| 354 | +} |
| 355 | + |
0 commit comments