malachite_q/arithmetic/
log_base.rs

1// Copyright © 2025 Mikhail Hogrefe
2//
3// Uses code adopted from the FLINT Library.
4//
5//      Copyright (C) 2011 Fredrik Johansson
6//
7// This file is part of Malachite.
8//
9// Malachite is free software: you can redistribute it and/or modify it under the terms of the GNU
10// Lesser General Public License (LGPL) as published by the Free Software Foundation; either version
11// 3 of the License, or (at your option) any later version. See <https://www.gnu.org/licenses/>.
12
13use crate::Rational;
14use core::cmp::Ordering::*;
15use malachite_base::num::arithmetic::traits::{
16    CeilingLogBase, CeilingLogBasePowerOf2, CheckedLogBase, CheckedLogBase2,
17    CheckedLogBasePowerOf2, FloorLogBase, FloorLogBasePowerOf2, Pow,
18};
19use malachite_base::num::comparison::traits::OrdAbs;
20use malachite_base::num::conversion::traits::{RoundingFrom, SciMantissaAndExponent};
21use malachite_base::rounding_modes::RoundingMode::*;
22
23fn approx_log_helper(x: &Rational) -> f64 {
24    let (mantissa, exponent): (f64, i64) = x.sci_mantissa_and_exponent();
25    mantissa.ln() + (exponent as f64) * core::f64::consts::LN_2
26}
27
28impl Rational {
29    /// Calculates the approximate natural logarithm of a positive [`Rational`].
30    ///
31    /// $f(x) = (1+\varepsilon)(\log x)$, where $|\varepsilon| < 2^{-52}.$
32    ///
33    /// # Worst-case complexity
34    /// $T(n) = O(n)$
35    ///
36    /// $M(n) = O(1)$
37    ///
38    /// where $T$ is time, $M$ is additional memory, and $n$ is `self.significant_bits()`.
39    ///
40    /// # Examples
41    /// ```
42    /// use malachite_base::num::arithmetic::traits::{Pow, PowerOf2};
43    /// use malachite_base::num::float::NiceFloat;
44    /// use malachite_q::Rational;
45    ///
46    /// assert_eq!(
47    ///     NiceFloat(Rational::from(10i32).approx_log()),
48    ///     NiceFloat(2.3025850929940455)
49    /// );
50    /// assert_eq!(
51    ///     NiceFloat(Rational::from(10i32).pow(100u64).approx_log()),
52    ///     NiceFloat(230.25850929940455)
53    /// );
54    /// assert_eq!(
55    ///     NiceFloat(Rational::power_of_2(1000000u64).approx_log()),
56    ///     NiceFloat(693147.1805599453)
57    /// );
58    /// assert_eq!(
59    ///     NiceFloat(Rational::power_of_2(-1000000i64).approx_log()),
60    ///     NiceFloat(-693147.1805599453)
61    /// );
62    /// ```
63    ///
64    /// This is equivalent to `fmpz_dlog` from `fmpz/dlog.c`, FLINT 2.7.1.
65    #[inline]
66    pub fn approx_log(&self) -> f64 {
67        assert!(*self > 0u32);
68        approx_log_helper(self)
69    }
70}
71
72// # Worst-case complexity
73// $T(n, m) = O(nm \log (nm) \log\log (nm))$
74//
75// $M(n, m) = O(nm \log (nm))$
76//
77// where $T$ is time, $M$ is additional memory, $n$ is `base.significant_bits()`, and $m$ is
78// $|\log_b x|$, where $b$ is `base` and $x$ is `x`.
79pub(crate) fn log_base_helper(x: &Rational, base: &Rational) -> (i64, bool) {
80    assert!(*base > 0u32);
81    assert_ne!(*base, 1u32);
82    if *x == 1u32 {
83        return (0, true);
84    }
85    let mut log = i64::rounding_from(approx_log_helper(x) / approx_log_helper(base), Floor).0;
86    let mut power = base.pow(log);
87    if *base > 1u32 {
88        match power.cmp_abs(x) {
89            Equal => (log, true),
90            Less => loop {
91                power *= base;
92                match power.cmp_abs(x) {
93                    Equal => {
94                        return (log + 1, true);
95                    }
96                    Less => {
97                        log += 1;
98                    }
99                    Greater => {
100                        return (log, false);
101                    }
102                }
103            },
104            Greater => loop {
105                power /= base;
106                match power.cmp_abs(x) {
107                    Equal => {
108                        return (log - 1, true);
109                    }
110                    Less => {
111                        return (log - 1, false);
112                    }
113                    Greater => {
114                        log -= 1;
115                    }
116                }
117            },
118        }
119    } else {
120        match power.cmp_abs(x) {
121            Equal => (log, true),
122            Less => loop {
123                power /= base;
124                match power.cmp_abs(x) {
125                    Equal => {
126                        return (log - 1, true);
127                    }
128                    Less => {
129                        log -= 1;
130                    }
131                    Greater => {
132                        return (log - 1, false);
133                    }
134                }
135            },
136            Greater => loop {
137                power *= base;
138                match power.cmp_abs(x) {
139                    Equal => {
140                        return (log + 1, true);
141                    }
142                    Less => {
143                        return (log, false);
144                    }
145                    Greater => {
146                        log += 1;
147                    }
148                }
149            },
150        }
151    }
152}
153
154impl FloorLogBase<&Rational> for &Rational {
155    type Output = i64;
156
157    /// Returns the floor of the base-$b$ logarithm of a positive [`Rational`].
158    ///
159    /// Note that this function may be slow if the base is very close to 1.
160    ///
161    /// $f(x, b) = \lfloor\log_b x\rfloor$.
162    ///
163    /// # Worst-case complexity
164    /// $T(n, m) = O(nm \log (nm) \log\log (nm))$
165    ///
166    /// $M(n, m) = O(nm \log (nm))$
167    ///
168    /// where $T$ is time, $M$ is additional memory, $n$ is `base.significant_bits()`, and $m$ is
169    /// $|\log_b x|$, where $b$ is `base` and $x$ is `x`.
170    ///
171    /// # Panics
172    /// Panics if `self` less than or equal to zero or `base` is 1.
173    ///
174    /// # Examples
175    /// ```
176    /// use malachite_base::num::arithmetic::traits::FloorLogBase;
177    /// use malachite_q::Rational;
178    ///
179    /// assert_eq!(
180    ///     Rational::from(80u32).floor_log_base(&Rational::from(3u32)),
181    ///     3
182    /// );
183    /// assert_eq!(
184    ///     Rational::from(81u32).floor_log_base(&Rational::from(3u32)),
185    ///     4
186    /// );
187    /// assert_eq!(
188    ///     Rational::from(82u32).floor_log_base(&Rational::from(3u32)),
189    ///     4
190    /// );
191    /// assert_eq!(
192    ///     Rational::from(4294967296u64).floor_log_base(&Rational::from(10u32)),
193    ///     9
194    /// );
195    /// assert_eq!(
196    ///     Rational::from_signeds(936851431250i64, 1397).floor_log_base(&Rational::from(10u32)),
197    ///     8
198    /// );
199    /// assert_eq!(
200    ///     Rational::from_signeds(5153632, 16807).floor_log_base(&Rational::from_signeds(22, 7)),
201    ///     5
202    /// );
203    /// ```
204    fn floor_log_base(self, base: &Rational) -> i64 {
205        assert!(*self > 0u32);
206        if let Some(log_base) = base.checked_log_base_2() {
207            return self.floor_log_base_power_of_2(log_base);
208        }
209        log_base_helper(self, base).0
210    }
211}
212
213impl CeilingLogBase<&Rational> for &Rational {
214    type Output = i64;
215
216    /// Returns the ceiling of the base-$b$ logarithm of a positive [`Rational`].
217    ///
218    /// Note that this function may be slow if the base is very close to 1.
219    ///
220    /// $f(x, b) = \lceil\log_b x\rceil$.
221    ///
222    /// # Worst-case complexity
223    /// $T(n, m) = O(nm \log (nm) \log\log (nm))$
224    ///
225    /// $M(n, m) = O(nm \log (nm))$
226    ///
227    /// where $T$ is time, $M$ is additional memory, $n$ is `base.significant_bits()`, and $m$ is
228    /// $|\log_b x|$, where $b$ is `base` and $x$ is `x`.
229    ///
230    /// # Panics
231    /// Panics if `self` less than or equal to zero or `base` is 1.
232    ///
233    /// # Examples
234    /// ```
235    /// use malachite_base::num::arithmetic::traits::CeilingLogBase;
236    /// use malachite_q::Rational;
237    ///
238    /// assert_eq!(
239    ///     Rational::from(80u32).ceiling_log_base(&Rational::from(3u32)),
240    ///     4
241    /// );
242    /// assert_eq!(
243    ///     Rational::from(81u32).ceiling_log_base(&Rational::from(3u32)),
244    ///     4
245    /// );
246    /// assert_eq!(
247    ///     Rational::from(82u32).ceiling_log_base(&Rational::from(3u32)),
248    ///     5
249    /// );
250    /// assert_eq!(
251    ///     Rational::from(4294967296u64).ceiling_log_base(&Rational::from(10u32)),
252    ///     10
253    /// );
254    /// assert_eq!(
255    ///     Rational::from_signeds(936851431250i64, 1397).ceiling_log_base(&Rational::from(10u32)),
256    ///     9
257    /// );
258    /// assert_eq!(
259    ///     Rational::from_signeds(5153632, 16807).ceiling_log_base(&Rational::from_signeds(22, 7)),
260    ///     5
261    /// );
262    /// ```
263    fn ceiling_log_base(self, base: &Rational) -> i64 {
264        assert!(*self > 0u32);
265        if let Some(log_base) = base.checked_log_base_2() {
266            return self.ceiling_log_base_power_of_2(log_base);
267        }
268        let (log, exact) = log_base_helper(self, base);
269        if exact { log } else { log + 1 }
270    }
271}
272
273impl CheckedLogBase<&Rational> for &Rational {
274    type Output = i64;
275
276    /// Returns the base-$b$ logarithm of a positive [`Rational`]. If the [`Rational`] is not a
277    /// power of $b$, then `None` is returned.
278    ///
279    /// Note that this function may be slow if the base is very close to 1.
280    ///
281    /// $$
282    /// f(x, b) = \\begin{cases}
283    ///     \operatorname{Some}(\log_b x) & \text{if} \\quad \log_b x \in \Z, \\\\
284    ///     \operatorname{None} & \textrm{otherwise}.
285    /// \\end{cases}
286    /// $$
287    ///
288    /// # Worst-case complexity
289    /// $T(n, m) = O(nm \log (nm) \log\log (nm))$
290    ///
291    /// $M(n, m) = O(nm \log (nm))$
292    ///
293    /// where $T$ is time, $M$ is additional memory, $n$ is `base.significant_bits()`, and $m$ is
294    /// $|\log_b x|$, where $b$ is `base` and $x$ is `x`.
295    ///
296    /// # Panics
297    /// Panics if `self` less than or equal to zero or `base` is 1.
298    ///
299    /// # Examples
300    /// ```
301    /// use malachite_base::num::arithmetic::traits::CheckedLogBase;
302    /// use malachite_q::Rational;
303    ///
304    /// assert_eq!(
305    ///     Rational::from(80u32).checked_log_base(&Rational::from(3u32)),
306    ///     None
307    /// );
308    /// assert_eq!(
309    ///     Rational::from(81u32).checked_log_base(&Rational::from(3u32)),
310    ///     Some(4)
311    /// );
312    /// assert_eq!(
313    ///     Rational::from(82u32).checked_log_base(&Rational::from(3u32)),
314    ///     None
315    /// );
316    /// assert_eq!(
317    ///     Rational::from(4294967296u64).checked_log_base(&Rational::from(10u32)),
318    ///     None
319    /// );
320    /// assert_eq!(
321    ///     Rational::from_signeds(936851431250i64, 1397).checked_log_base(&Rational::from(10u32)),
322    ///     None
323    /// );
324    /// assert_eq!(
325    ///     Rational::from_signeds(5153632, 16807).checked_log_base(&Rational::from_signeds(22, 7)),
326    ///     Some(5)
327    /// );
328    /// ```
329    fn checked_log_base(self, base: &Rational) -> Option<i64> {
330        assert!(*self > 0u32);
331        if let Some(log_base) = base.checked_log_base_2() {
332            return self.checked_log_base_power_of_2(log_base);
333        }
334        let (log, exact) = log_base_helper(self, base);
335        if exact { Some(log) } else { None }
336    }
337}