malachite_float/conversion/from_bits.rs
1// Copyright © 2025 Mikhail Hogrefe
2//
3// This file is part of Malachite.
4//
5// Malachite is free software: you can redistribute it and/or modify it under the terms of the GNU
6// Lesser General Public License (LGPL) as published by the Free Software Foundation; either version
7// 3 of the License, or (at your option) any later version. See <https://www.gnu.org/licenses/>.
8
9use crate::Float;
10use crate::InnerFloat::Finite;
11use alloc::vec;
12use core::cmp::Ordering::{self, *};
13use malachite_base::num::arithmetic::traits::{
14 DivisibleByPowerOf2, NegModPowerOf2, PowerOf2, ShrRound,
15};
16use malachite_base::num::basic::integers::PrimitiveInt;
17use malachite_base::num::conversion::traits::ExactFrom;
18use malachite_base::num::logic::traits::SignificantBits;
19use malachite_base::rounding_modes::RoundingMode::{self, *};
20use malachite_nz::natural::Natural;
21use malachite_nz::platform::Limb;
22
23const HIGH_BIT: Limb = 1 << (Limb::WIDTH - 1);
24
25impl Float {
26 /// Returns an approximation of a real number, given the number's bits. To avoid troublesome
27 /// edge cases, the number should not be a dyadic rational (and the iterator of bits should
28 /// therefore be infinite, and not eventually 0 or 1). Given this assumption, the rounding mode
29 /// `Exact` should never be passed.
30 ///
31 /// The approximation has precision `prec` and is rounded according to the provided rounding
32 /// mode.
33 ///
34 /// This function reads `prec + z` bits, or `prec + z + 1` bits if `rm` is `Nearest`, where `z`
35 /// is the number of leading false bits in `bits`.
36 ///
37 /// This function always produces a value in the interval $[1/2,1]$. In particular, it never
38 /// overflows or underflows.
39 ///
40 /// $$
41 /// f((x_k),p,m) = C+\varepsilon,
42 /// $$
43 /// where
44 /// $$
45 /// C=\sum_{k=0}^\infty x_k 2^{-(k+1)}.
46 /// $$
47 /// - If $m$ is not `Nearest`, then $|\varepsilon| < 2^{\lfloor\log_2 C\rfloor-p+1}$.
48 /// - If $m$ is `Nearest`, then $|\varepsilon| < 2^{\lfloor\log_2 C\rfloor-p}$.
49 ///
50 /// The output has precision `prec`.
51 ///
52 /// # Worst-case complexity
53 /// $T(n) = O(n)$
54 ///
55 /// $M(n) = O(n)$
56 ///
57 /// where $T$ is time, $M$ is additional memory, and $n$ is `prec`.
58 ///
59 /// # Panics
60 /// Panics if `prec` is zero or `rm` is `Exact`.
61 ///
62 /// # Examples
63 /// ```
64 /// use malachite_base::rounding_modes::RoundingMode::*;
65 /// use malachite_float::Float;
66 /// use std::cmp::Ordering::*;
67 ///
68 /// // Produces 10100100010000100000...
69 /// struct Bits {
70 /// b: bool,
71 /// k: usize,
72 /// j: usize,
73 /// }
74 ///
75 /// impl Iterator for Bits {
76 /// type Item = bool;
77 ///
78 /// fn next(&mut self) -> Option<bool> {
79 /// Some(if self.b {
80 /// self.b = false;
81 /// self.j = self.k;
82 /// true
83 /// } else {
84 /// self.j -= 1;
85 /// if self.j == 0 {
86 /// self.k += 1;
87 /// self.b = true;
88 /// }
89 /// false
90 /// })
91 /// }
92 /// }
93 ///
94 /// impl Bits {
95 /// fn new() -> Bits {
96 /// Bits {
97 /// b: true,
98 /// k: 1,
99 /// j: 1,
100 /// }
101 /// }
102 /// }
103 ///
104 /// let (c, o) = Float::non_dyadic_from_bits_prec_round(Bits::new(), 100, Floor);
105 /// assert_eq!(c.to_string(), "0.6416325606551538662938427702254");
106 /// assert_eq!(o, Less);
107 ///
108 /// let (c, o) = Float::non_dyadic_from_bits_prec_round(Bits::new(), 100, Ceiling);
109 /// assert_eq!(c.to_string(), "0.641632560655153866293842770226");
110 /// assert_eq!(o, Greater);
111 /// ```
112 pub fn non_dyadic_from_bits_prec_round<I: Iterator<Item = bool>>(
113 mut bits: I,
114 prec: u64,
115 rm: RoundingMode,
116 ) -> (Self, Ordering) {
117 assert_ne!(prec, 0);
118 assert_ne!(rm, Exact);
119 let len = usize::exact_from(prec.shr_round(Limb::LOG_WIDTH, Ceiling).0);
120 let mut limbs = vec![0; len];
121 let mut limbs_it = limbs.iter_mut().rev();
122 let mut x = limbs_it.next().unwrap();
123 let mut mask = HIGH_BIT;
124 let mut seen_one = false;
125 let mut exponent: i32 = 0;
126 let mut remaining = prec;
127 for b in &mut bits {
128 if !seen_one {
129 if b {
130 seen_one = true;
131 } else {
132 exponent = exponent.checked_sub(1).unwrap();
133 continue;
134 }
135 }
136 if b {
137 *x |= mask;
138 }
139 remaining -= 1;
140 if remaining == 0 {
141 break;
142 }
143 if mask == 1 {
144 x = limbs_it.next().unwrap();
145 mask = HIGH_BIT;
146 } else {
147 mask >>= 1;
148 }
149 }
150 let mut significand = Natural::from_owned_limbs_asc(limbs);
151 let increment = rm == Up || rm == Ceiling || (rm == Nearest && bits.next() == Some(true));
152 if increment {
153 significand +=
154 Natural::from(Limb::power_of_2(prec.neg_mod_power_of_2(Limb::LOG_WIDTH)));
155 if !significand
156 .significant_bits()
157 .divisible_by_power_of_2(Limb::LOG_WIDTH)
158 {
159 significand >>= 1;
160 exponent += 1;
161 }
162 }
163 (
164 Self(Finite {
165 sign: true,
166 exponent,
167 precision: prec,
168 significand,
169 }),
170 if increment { Greater } else { Less },
171 )
172 }
173
174 /// Returns an approximation of a real number, given the number's bits. To avoid troublesome
175 /// edge cases, the number should not be a dyadic rational (and the iterator of bits should
176 /// therefore be infinite, and not eventually 0 or 1).
177 ///
178 /// The approximation has precision `prec` and is rounded according to the `Nearest` rounding
179 /// mode.
180 ///
181 /// This function reads `prec + z + 1` bits, where `z` is the number of leading false bits in
182 /// `bits`.
183 ///
184 /// This function always produces a value in the interval $[1/2,1]$. In particular, it never
185 /// overflows or underflows.
186 ///
187 /// $$
188 /// f((x_k),p,m) = C+\varepsilon,
189 /// $$
190 /// where
191 /// $$
192 /// C=\sum_{k=0}^\infty x_k 2^{-(k+1)}
193 /// $$
194 /// and $|\varepsilon| < 2^{\lfloor\log_2 C\rfloor-p}$.
195 ///
196 /// The output has precision `prec`.
197 ///
198 /// # Worst-case complexity
199 /// $T(n) = O(n)$
200 ///
201 /// $M(n) = O(n)$
202 ///
203 /// where $T$ is time, $M$ is additional memory, and $n$ is `prec`.
204 ///
205 /// # Panics
206 /// Panics if `prec` is zero.
207 ///
208 /// # Examples
209 /// ```
210 /// use malachite_float::Float;
211 /// use std::cmp::Ordering::*;
212 ///
213 /// // Produces 10100100010000100000...
214 /// struct Bits {
215 /// b: bool,
216 /// k: usize,
217 /// j: usize,
218 /// }
219 ///
220 /// impl Iterator for Bits {
221 /// type Item = bool;
222 ///
223 /// fn next(&mut self) -> Option<bool> {
224 /// Some(if self.b {
225 /// self.b = false;
226 /// self.j = self.k;
227 /// true
228 /// } else {
229 /// self.j -= 1;
230 /// if self.j == 0 {
231 /// self.k += 1;
232 /// self.b = true;
233 /// }
234 /// false
235 /// })
236 /// }
237 /// }
238 ///
239 /// impl Bits {
240 /// fn new() -> Bits {
241 /// Bits {
242 /// b: true,
243 /// k: 1,
244 /// j: 1,
245 /// }
246 /// }
247 /// }
248 ///
249 /// let (c, o) = Float::non_dyadic_from_bits_prec(Bits::new(), 1);
250 /// assert_eq!(c.to_string(), "0.5");
251 /// assert_eq!(o, Less);
252 ///
253 /// let (c, o) = Float::non_dyadic_from_bits_prec(Bits::new(), 10);
254 /// assert_eq!(c.to_string(), "0.642");
255 /// assert_eq!(o, Less);
256 ///
257 /// let (c, o) = Float::non_dyadic_from_bits_prec(Bits::new(), 100);
258 /// assert_eq!(c.to_string(), "0.6416325606551538662938427702254");
259 /// assert_eq!(o, Less);
260 /// ```
261 #[inline]
262 pub fn non_dyadic_from_bits_prec<I: Iterator<Item = bool>>(
263 bits: I,
264 prec: u64,
265 ) -> (Self, Ordering) {
266 Self::non_dyadic_from_bits_prec_round(bits, prec, Nearest)
267 }
268}