dtoa/
diyfp.rs

1// Licensed under the Apache License, Version 2.0 <LICENSE-APACHE or
2// http://www.apache.org/licenses/LICENSE-2.0> or the MIT license
3// <LICENSE-MIT or http://opensource.org/licenses/MIT>, at your
4// option. This file may not be copied, modified, or distributed
5// except according to those terms.
6//
7// ---
8//
9// The C++ implementation preserved here in comments is licensed as follows:
10//
11// Tencent is pleased to support the open source community by making RapidJSON
12// available.
13//
14// Copyright (C) 2015 THL A29 Limited, a Tencent company, and Milo Yip. All
15// rights reserved.
16//
17// Licensed under the MIT License (the "License"); you may not use this file
18// except in compliance with the License. You may obtain a copy of the License
19// at
20//
21// http://opensource.org/licenses/MIT
22//
23// Unless required by applicable law or agreed to in writing, software
24// distributed under the License is distributed on an "AS IS" BASIS, WITHOUT
25// WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied. See the
26// License for the specific language governing permissions and limitations under
27// the License.
28
29use core::ops::{Mul, Sub};
30#[cfg(feature = "no-panic")]
31use no_panic::no_panic;
32
33#[derive(Copy, Clone, Debug)]
34pub struct DiyFp<F, E> {
35    pub f: F,
36    pub e: E,
37}
38
39impl<F, E> DiyFp<F, E> {
40    #[cfg_attr(feature = "no-panic", no_panic)]
41    pub fn new(f: F, e: E) -> Self {
42        DiyFp { f, e }
43    }
44}
45
46impl<F, E> Sub for DiyFp<F, E>
47where
48    F: Sub<F, Output = F>,
49{
50    type Output = Self;
51
52    #[cfg_attr(feature = "no-panic", no_panic)]
53    fn sub(self, rhs: Self) -> Self {
54        DiyFp {
55            f: self.f - rhs.f,
56            e: self.e,
57        }
58    }
59}
60
61impl Mul for DiyFp<u32, i32> {
62    type Output = Self;
63
64    #[cfg_attr(feature = "no-panic", no_panic)]
65    fn mul(self, rhs: Self) -> Self {
66        let mut tmp = self.f as u64 * rhs.f as u64;
67        tmp += 1u64 << 31; // mult_round
68        DiyFp {
69            f: (tmp >> 32) as u32,
70            e: self.e + rhs.e + 32,
71        }
72    }
73}
74
75impl Mul for DiyFp<u64, isize> {
76    type Output = Self;
77
78    #[cfg_attr(feature = "no-panic", no_panic)]
79    fn mul(self, rhs: Self) -> Self {
80        let m32 = 0xFFFFFFFFu64;
81        let a = self.f >> 32;
82        let b = self.f & m32;
83        let c = rhs.f >> 32;
84        let d = rhs.f & m32;
85        let ac = a * c;
86        let bc = b * c;
87        let ad = a * d;
88        let bd = b * d;
89        let mut tmp = (bd >> 32) + (ad & m32) + (bc & m32);
90        tmp += 1u64 << 31; // mult_round
91        DiyFp {
92            f: ac + (ad >> 32) + (bc >> 32) + (tmp >> 32),
93            e: self.e + rhs.e + 64,
94        }
95    }
96}
97
98macro_rules! diyfp {
99    (
100        floating_type: $fty:ty,
101        significand_type: $sigty:ty,
102        exponent_type: $expty:ty,
103
104        diy_significand_size: $diy_significand_size:expr,
105        significand_size: $significand_size:expr,
106        exponent_bias: $exponent_bias:expr,
107        mask_type: $mask_type:ty,
108        exponent_mask: $exponent_mask:expr,
109        significand_mask: $significand_mask:expr,
110        hidden_bit: $hidden_bit:expr,
111        cached_powers_f: $cached_powers_f:expr,
112        cached_powers_e: $cached_powers_e:expr,
113        min_power: $min_power:expr,
114    ) => {
115        type DiyFp = diyfp::DiyFp<$sigty, $expty>;
116
117        impl DiyFp {
118            // Preconditions:
119            // `d` must have a positive sign and must not be infinity or NaN.
120            /*
121            explicit DiyFp(double d) {
122                union {
123                    double d;
124                    uint64_t u64;
125                } u = { d };
126
127                int biased_e = static_cast<int>((u.u64 & kDpExponentMask) >> kDpSignificandSize);
128                uint64_t significand = (u.u64 & kDpSignificandMask);
129                if (biased_e != 0) {
130                    f = significand + kDpHiddenBit;
131                    e = biased_e - kDpExponentBias;
132                }
133                else {
134                    f = significand;
135                    e = kDpMinExponent + 1;
136                }
137            }
138            */
139            #[cfg_attr(feature = "no-panic", no_panic)]
140            unsafe fn from(d: $fty) -> Self {
141                let u: $mask_type = mem::transmute(d);
142
143                let biased_e = ((u & $exponent_mask) >> $significand_size) as $expty;
144                let significand = u & $significand_mask;
145                if biased_e != 0 {
146                    DiyFp {
147                        f: significand + $hidden_bit,
148                        e: biased_e - $exponent_bias - $significand_size,
149                    }
150                } else {
151                    DiyFp {
152                        f: significand,
153                        e: 1 - $exponent_bias - $significand_size,
154                    }
155                }
156            }
157
158            // Normalizes so that the highest bit of the diy significand is 1.
159            /*
160            DiyFp Normalize() const {
161                DiyFp res = *this;
162                while (!(res.f & (static_cast<uint64_t>(1) << 63))) {
163                    res.f <<= 1;
164                    res.e--;
165                }
166                return res;
167            }
168            */
169            #[cfg_attr(feature = "no-panic", no_panic)]
170            fn normalize(self) -> DiyFp {
171                let mut res = self;
172                while (res.f & (1 << ($diy_significand_size - 1))) == 0 {
173                    res.f <<= 1;
174                    res.e -= 1;
175                }
176                res
177            }
178
179            // Normalizes so that the highest bit of the diy significand is 1.
180            //
181            // Precondition:
182            // `self.f` must be no more than 2 bits longer than the f64 significand.
183            /*
184            DiyFp NormalizeBoundary() const {
185                DiyFp res = *this;
186                while (!(res.f & (kDpHiddenBit << 1))) {
187                    res.f <<= 1;
188                    res.e--;
189                }
190                res.f <<= (kDiySignificandSize - kDpSignificandSize - 2);
191                res.e = res.e - (kDiySignificandSize - kDpSignificandSize - 2);
192                return res;
193            }
194            */
195            #[cfg_attr(feature = "no-panic", no_panic)]
196            fn normalize_boundary(self) -> DiyFp {
197                let mut res = self;
198                while (res.f & $hidden_bit << 1) == 0 {
199                    res.f <<= 1;
200                    res.e -= 1;
201                }
202                res.f <<= $diy_significand_size - $significand_size - 2;
203                res.e -= $diy_significand_size - $significand_size - 2;
204                res
205            }
206
207            // Normalizes `self - e` and `self + e` where `e` is half of the least
208            // significant digit of `self`. The plus is normalized so that the highest
209            // bit of the diy significand is 1. The minus is normalized so that it has
210            // the same exponent as the plus.
211            //
212            // Preconditions:
213            // `self` must have been returned directly from `DiyFp::from_f64`.
214            // `self.f` must not be zero.
215            /*
216            void NormalizedBoundaries(DiyFp* minus, DiyFp* plus) const {
217                DiyFp pl = DiyFp((f << 1) + 1, e - 1).NormalizeBoundary();
218                DiyFp mi = (f == kDpHiddenBit) ? DiyFp((f << 2) - 1, e - 2) : DiyFp((f << 1) - 1, e - 1);
219                mi.f <<= mi.e - pl.e;
220                mi.e = pl.e;
221                *plus = pl;
222                *minus = mi;
223            }
224            */
225            #[cfg_attr(feature = "no-panic", no_panic)]
226            fn normalized_boundaries(self) -> (DiyFp, DiyFp) {
227                let pl = DiyFp::new((self.f << 1) + 1, self.e - 1).normalize_boundary();
228                let mut mi = if self.f == $hidden_bit {
229                    DiyFp::new((self.f << 2) - 1, self.e - 2)
230                } else {
231                    DiyFp::new((self.f << 1) - 1, self.e - 1)
232                };
233                mi.f <<= mi.e - pl.e;
234                mi.e = pl.e;
235                (mi, pl)
236            }
237        }
238
239        /*
240        inline DiyFp GetCachedPower(int e, int* K) {
241            //int k = static_cast<int>(ceil((-61 - e) * 0.30102999566398114)) + 374;
242            double dk = (-61 - e) * 0.30102999566398114 + 347;  // dk must be positive, so can do ceiling in positive
243            int k = static_cast<int>(dk);
244            if (dk - k > 0.0)
245                k++;
246
247            unsigned index = static_cast<unsigned>((k >> 3) + 1);
248            *K = -(-348 + static_cast<int>(index << 3));    // decimal exponent no need lookup table
249
250            return GetCachedPowerByIndex(index);
251        }
252        */
253        #[inline]
254        #[cfg_attr(feature = "no-panic", no_panic)]
255        fn get_cached_power(e: $expty) -> (DiyFp, isize) {
256            let dk = (3 - $diy_significand_size - e) as f64 * 0.30102999566398114f64
257                - ($min_power + 1) as f64;
258            let mut k = dk as isize;
259            if dk - k as f64 > 0.0 {
260                k += 1;
261            }
262
263            let index = ((k >> 3) + 1) as usize;
264            let k = -($min_power + (index << 3) as isize);
265
266            (
267                DiyFp::new(*unsafe { $cached_powers_f.get_unchecked(index) }, *unsafe {
268                    $cached_powers_e.get_unchecked(index)
269                }
270                    as $expty),
271                k,
272            )
273        }
274    };
275}