malachite_float/float/conversion/from_digits.rs
1// Copyright © 2026 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 alloc::vec::Vec;
11use core::cmp::Ordering::{self, *};
12use core::cmp::max;
13use malachite_base::num::arithmetic::traits::{CheckedLogBase2, FloorLogBase2, Pow, PowerOf2};
14use malachite_base::num::basic::integers::PrimitiveInt;
15use malachite_base::num::basic::traits::One;
16use malachite_base::num::conversion::traits::{Digits, ExactFrom};
17use malachite_base::num::logic::traits::{BitAccess, SignificantBits};
18use malachite_base::rounding_modes::RoundingMode::{self, *};
19use malachite_nz::natural::Natural;
20use malachite_q::Rational;
21
22// Rounds the exact quotient of two `Natural`s, the second nonzero, to a `Float` of the given
23// precision. The conversions are exact, so the division sees the exact quotient and the returned
24// `Ordering` describes the quotient itself.
25fn quotient_prec_round(n: Natural, d: Natural, prec: u64, rm: RoundingMode) -> (Float, Ordering) {
26 // A `Natural` needs an exponent of `significant_bits` to be a `Float`, so past `MAX_EXPONENT`
27 // it has no `Float` of its own even when the quotient is perfectly ordinary -- which is the
28 // usual case here, since the quotient is a digit expansion less than 1 while the numerator and
29 // denominator both grow with the precision. Building the `Rational` first never forms either
30 // endpoint as a `Float`. It costs a gcd, so it is worth avoiding until it is needed.
31 //
32 // Both arguments are taken by value: at these sizes the conversions would otherwise copy them,
33 // and one caller-side clone is cheaper than two internal ones.
34 if max(n.significant_bits(), d.significant_bits()) > Float::MAX_EXPONENT_U64 {
35 Float::from_rational_prec_round(Rational::from_naturals(n, d), prec, rm)
36 } else {
37 Float::exact_from(n).div_prec_round(Float::exact_from(d), prec, rm)
38 }
39}
40
41impl Float {
42 /// Returns an approximation of a real number, given the number's digits in a base that is a
43 /// power of 2.
44 ///
45 /// Each digit contributes exactly $\log_2 b$ bits, so this is
46 /// [`non_dyadic_from_bits_prec_round`](Float::non_dyadic_from_bits_prec_round) with the digits
47 /// expanded, and it reads the same number of digits that that function reads bits, rounded up
48 /// to a whole digit.
49 ///
50 /// $$
51 /// f((x_k),b,p,m) = C+\varepsilon, \quad C=\sum_{k=0}^\infty x_k b^{-(k+1)}.
52 /// $$
53 /// - If $m$ is not `Nearest`, then $|\varepsilon| < 2^{\lfloor\log_2 C\rfloor-p+1}$.
54 /// - If $m$ is `Nearest`, then $|\varepsilon| < 2^{\lfloor\log_2 C\rfloor-p}$.
55 ///
56 /// The output has precision `prec`.
57 ///
58 /// # Preconditions
59 /// $C$ must not be a dyadic rational: the digits must be infinite and not eventually all zero
60 /// or all $b-1$. Given that, the rounded value never equals $C$, so the returned [`Ordering`]
61 /// is never `Equal` and `Exact` is never a sensible rounding mode. $C$ must also be less than
62 /// 1, which holds whenever the digits are read as lying wholly after the point.
63 ///
64 /// # Worst-case complexity
65 /// $T(n) = O(n)$
66 ///
67 /// $M(n) = O(n)$
68 ///
69 /// where $T$ is time, $M$ is additional memory, and $n$ is `prec`.
70 ///
71 /// # Panics
72 /// Panics if `log_base` is zero or greater than 64, if a digit is greater than or equal to
73 /// $2^{\ell}$, if `prec` is zero, or if `rm` is `Exact`.
74 ///
75 /// # Examples
76 /// ```
77 /// use malachite_base::rounding_modes::RoundingMode::*;
78 /// use malachite_float::Float;
79 /// use std::cmp::Ordering::*;
80 ///
81 /// // 0.4444... in base 16 is 4/15
82 /// let (x, o) = Float::non_dyadic_from_power_of_2_digits_prec_round(
83 /// core::iter::repeat(4),
84 /// 4,
85 /// 20,
86 /// Floor,
87 /// );
88 /// assert_eq!(x.to_string(), "0.26666641");
89 /// assert_eq!(o, Less);
90 /// ```
91 pub fn non_dyadic_from_power_of_2_digits_prec_round<I: Iterator<Item = u64>>(
92 digits: I,
93 log_base: u64,
94 prec: u64,
95 rm: RoundingMode,
96 ) -> (Self, Ordering) {
97 assert_ne!(log_base, 0);
98 assert!(log_base <= u64::WIDTH);
99 Self::non_dyadic_from_bits_prec_round(
100 digits.flat_map(move |d| {
101 assert!(
102 log_base == u64::WIDTH || d < u64::power_of_2(log_base),
103 "digit out of range"
104 );
105 (0..log_base).rev().map(move |i| d.get_bit(i))
106 }),
107 prec,
108 rm,
109 )
110 }
111
112 /// Returns an approximation of a real number, given the number's digits in a base that is a
113 /// power of 2, rounding to nearest.
114 ///
115 /// See [`non_dyadic_from_power_of_2_digits_prec_round`](
116 /// Float::non_dyadic_from_power_of_2_digits_prec_round) for details and preconditions.
117 ///
118 /// # Worst-case complexity
119 /// $T(n) = O(n)$
120 ///
121 /// $M(n) = O(n)$
122 ///
123 /// where $T$ is time, $M$ is additional memory, and $n$ is `prec`.
124 ///
125 /// # Panics
126 /// Panics if `log_base` is zero or greater than 64, if a digit is greater than or equal to
127 /// $2^{\ell}$, or if `prec` is zero.
128 ///
129 /// # Examples
130 /// ```
131 /// use malachite_float::Float;
132 /// use std::cmp::Ordering::*;
133 ///
134 /// let (x, o) = Float::non_dyadic_from_power_of_2_digits_prec(core::iter::repeat(4), 4, 20);
135 /// assert_eq!(x.to_string(), "0.26666689");
136 /// assert_eq!(o, Greater);
137 /// ```
138 #[inline]
139 pub fn non_dyadic_from_power_of_2_digits_prec<I: Iterator<Item = u64>>(
140 digits: I,
141 log_base: u64,
142 prec: u64,
143 ) -> (Self, Ordering) {
144 Self::non_dyadic_from_power_of_2_digits_prec_round(digits, log_base, prec, Nearest)
145 }
146
147 /// Returns an approximation of a real number, given the number's digits in an arbitrary base.
148 ///
149 /// A digit in a base that is not a power of 2 does not correspond to a whole number of bits, so
150 /// this reads a batch of digits, brackets $C$ between the values those digits allow, and reads
151 /// more if the bracket is not yet narrow enough to determine both the rounded value and its
152 /// position relative to $C$. When the base is a power of 2 it defers to
153 /// [`non_dyadic_from_power_of_2_digits_prec_round`](
154 /// Float::non_dyadic_from_power_of_2_digits_prec_round), which needs no such loop.
155 ///
156 /// $$
157 /// f((x_k),b,p,m) = C+\varepsilon, \quad C=\sum_{k=0}^\infty x_k b^{-(k+1)}.
158 /// $$
159 /// - If $m$ is not `Nearest`, then $|\varepsilon| < 2^{\lfloor\log_2 C\rfloor-p+1}$.
160 /// - If $m$ is `Nearest`, then $|\varepsilon| < 2^{\lfloor\log_2 C\rfloor-p}$.
161 ///
162 /// The output has precision `prec`.
163 ///
164 /// # Preconditions
165 /// $C$ must not be a dyadic rational. If it is, the bracket can never separate it from the
166 /// [`Float`] that equals it and this function does not terminate. Note that this is a condition
167 /// on $C$, not on the digits: in a base that is not a power of 2 a dyadic rational has a
168 /// non-terminating expansion, as $1/2$ does in base 3. Given the precondition, the rounded
169 /// value never equals $C$, so the returned [`Ordering`] is never `Equal` and `Exact` is never a
170 /// sensible rounding mode. $C$ must also be less than 1, and the iterator must be infinite.
171 ///
172 /// # Worst-case complexity
173 /// $T(n) = O(n (\log n)^2 \log\log n)$
174 ///
175 /// $M(n) = O(n \log n)$
176 ///
177 /// where $T$ is time, $M$ is additional memory, and $n$ is `prec`.
178 ///
179 /// # Panics
180 /// Panics if `base` is less than 2, if a digit is greater than or equal to `base`, if `prec` is
181 /// zero, or if `rm` is `Exact`.
182 ///
183 /// # Examples
184 /// ```
185 /// use malachite_base::rounding_modes::RoundingMode::*;
186 /// use malachite_float::Float;
187 /// use std::cmp::Ordering::*;
188 ///
189 /// // 0.3333... in base 10 is 1/3
190 /// let (x, o) = Float::non_dyadic_from_digits_prec_round(core::iter::repeat(3), 10, 20, Floor);
191 /// assert_eq!(x.to_string(), "0.33333302");
192 /// assert_eq!(o, Less);
193 /// ```
194 pub fn non_dyadic_from_digits_prec_round<I: Iterator<Item = u64>>(
195 mut digits: I,
196 base: u64,
197 prec: u64,
198 rm: RoundingMode,
199 ) -> (Self, Ordering) {
200 assert!(base >= 2, "base out of range");
201 assert_ne!(prec, 0);
202 assert_ne!(rm, Exact);
203 if let Some(log_base) = base.checked_log_base_2() {
204 return Self::non_dyadic_from_power_of_2_digits_prec_round(digits, log_base, prec, rm);
205 }
206 // Each digit carries at least `floor(log2(base))` bits, so this many digits is enough to
207 // cover the precision, with a few to spare for the rounding decision.
208 let per_digit = base.floor_log_base_2();
209 let mut target = (prec + 8).div_ceil(per_digit);
210 let base_n = Natural::from(base);
211 let mut buf: Vec<u64> = Vec::new();
212 loop {
213 while (buf.len() as u64) < target {
214 let d = digits
215 .next()
216 .expect("the digit iterator must not run out; see the preconditions");
217 assert!(d < base, "digit out of range");
218 buf.push(d);
219 }
220 let n = Natural::from_digits_desc(&base, buf.iter().copied()).unwrap();
221 let power = (&base_n).pow(target);
222 // C lies strictly between n / power and (n + 1) / power. The second call consumes both
223 // values, so only the first has to clone.
224 let (lo, o_lo) = quotient_prec_round(n.clone(), power.clone(), prec, rm);
225 let (hi, o_hi) = quotient_prec_round(n + Natural::ONE, power, prec, rm);
226 if lo == hi {
227 // Rounding is monotonic, so `lo` is the rounding of everything in the bracket, and
228 // hence of C. It only remains to place it relative to C, which the orderings of the
229 // endpoints settle whenever `lo` falls outside the bracket.
230 if o_lo != Greater {
231 return (lo, Less);
232 }
233 if o_hi != Less {
234 return (hi, Greater);
235 }
236 }
237 target += max(4, target >> 1);
238 }
239 }
240
241 /// Returns an approximation of a real number, given the number's digits in an arbitrary base,
242 /// rounding to nearest.
243 ///
244 /// See [`non_dyadic_from_digits_prec_round`](Float::non_dyadic_from_digits_prec_round) for
245 /// details and preconditions.
246 ///
247 /// # Worst-case complexity
248 /// $T(n) = O(n (\log n)^2 \log\log n)$
249 ///
250 /// $M(n) = O(n \log n)$
251 ///
252 /// where $T$ is time, $M$ is additional memory, and $n$ is `prec`.
253 ///
254 /// # Panics
255 /// Panics if `base` is less than 2, if a digit is greater than or equal to `base`, or if `prec`
256 /// is zero.
257 ///
258 /// # Examples
259 /// ```
260 /// use malachite_float::Float;
261 /// use std::cmp::Ordering::*;
262 ///
263 /// let (x, o) = Float::non_dyadic_from_digits_prec(core::iter::repeat(3), 10, 20);
264 /// assert_eq!(x.to_string(), "0.33333349");
265 /// assert_eq!(o, Greater);
266 /// ```
267 #[inline]
268 pub fn non_dyadic_from_digits_prec<I: Iterator<Item = u64>>(
269 digits: I,
270 base: u64,
271 prec: u64,
272 ) -> (Self, Ordering) {
273 Self::non_dyadic_from_digits_prec_round(digits, base, prec, Nearest)
274 }
275}