malachite_float/float/arithmetic/log_base.rs
1// Copyright © 2026 Mikhail Hogrefe
2//
3// Uses code adopted from the GNU MPFR Library.
4//
5// Copyright 2001-2026 Free Software Foundation, Inc.
6//
7// Contributed by the Pascaline and Caramba projects, INRIA.
8//
9// This file is part of Malachite.
10//
11// Malachite is free software: you can redistribute it and/or modify it under the terms of the GNU
12// Lesser General Public License (LGPL) as published by the Free Software Foundation; either version
13// 3 of the License, or (at your option) any later version. See <https://www.gnu.org/licenses/>.
14
15use crate::InnerFloat::{Finite, Infinity, NaN, Zero};
16use crate::float::arithmetic::ln::{SliverOfOne, sliver_of_one};
17use crate::float::arithmetic::log_base_2::extended_log_base_2_of_rational;
18use crate::float::basic::extended::ExtendedFloat;
19use crate::{
20 Float, emulate_float_to_float_fn, emulate_rational_to_float_fn, float_either_zero,
21 float_infinity, float_nan, float_negative_infinity,
22};
23use alloc::vec::Vec;
24use core::cmp::Ordering::{self, *};
25use malachite_base::num::arithmetic::traits::{
26 CeilingLogBase2, CheckedLogBase, DivisibleBy, Gcd, IsPowerOf2, LogBase, LogBaseAssign, Mod,
27 ModAdd, ModMul, ModPow, Pow, Sign,
28};
29use malachite_base::num::basic::floats::PrimitiveFloat;
30use malachite_base::num::basic::integers::PrimitiveInt;
31use malachite_base::num::basic::traits::{One, Two, Zero as ZeroTrait};
32use malachite_base::num::conversion::traits::{ExactFrom, RoundingFrom};
33use malachite_base::num::factorization::traits::ExpressAsPower;
34use malachite_base::num::logic::traits::SignificantBits;
35use malachite_base::rounding_modes::RoundingMode::{self, *};
36use malachite_nz::integer::Integer;
37use malachite_nz::natural::Natural;
38use malachite_nz::natural::arithmetic::float::round::float_can_round;
39use malachite_nz::platform::Limb;
40use malachite_q::Rational;
41
42// Decomposes a positive finite nonzero `Float` into `(s, t)` with value `s * 2^t` and `s` odd. `s`
43// has at most `prec` bits, so everything downstream of this decomposition works on small numbers no
44// matter how extreme `x`'s exponent is.
45pub(crate) fn odd_significand_and_exponent(x: &Float) -> (Natural, i64) {
46 let sig = x.significand_ref().unwrap();
47 let strip = sig.trailing_zeros().unwrap();
48 let s = sig >> strip;
49 let t = i64::from(x.get_exponent().unwrap()) - i64::exact_from(sig.significant_bits())
50 + i64::exact_from(strip);
51 (s, t)
52}
53
54// Given a positive dyadic value `s * 2^t` (`s` odd) and a root `2^z * h` (`h` an odd `Natural`, the
55// root positive and different from 1), returns the integer `m` with `s * 2^t = (2^z * h)^m`, or
56// `None` if no such integer exists.
57//
58// All arithmetic is on `s`, `h`, and `i64` exponents: nothing the size of `2^|t|` is ever
59// materialized, so the check is cheap for arbitrary exponents. `m` may be negative only when `h =
60// 1` (an odd `h >= 3` in the root makes every negative power non-dyadic).
61pub(crate) fn dyadic_log_of_root(s: &Natural, t: i64, z: i64, h: &Natural) -> Option<i64> {
62 if *h == 1u32 {
63 // The root is 2^z, so s * 2^t = 2^(z * m) requires s = 1 and z | t.
64 return if *s == 1u32 && z != 0 && t.divisible_by(z) {
65 Some(t / z)
66 } else {
67 None
68 };
69 }
70 // h >= 3: (2^z * h)^m = 2^(z * m) * h^m with h^m odd, so s = h^m and t = z * m. A negative m
71 // would put h in the denominator, which the dyadic s * 2^t cannot cancel, so m >= 0; and m = 0
72 // means the value is 1, which callers handle separately.
73 if *s == 1u32 {
74 return None;
75 }
76 let m = i64::exact_from(s.checked_log_base(h)?);
77 if z.checked_mul(m)? == t {
78 Some(m)
79 } else {
80 None
81 }
82}
83
84// Given a positive dyadic base `s_b * 2^t_b` (`s_b` odd) with base not 1, returns `(z, h, e_base)`
85// such that `base = (2^z * h)^e_base` with `h` odd and `e_base` maximal (the primitive root). All
86// work is on `s_b` (at most the base's precision in bits) and `i64` exponents, so the base is never
87// materialized as an integer or `Rational`, no matter how extreme its exponent.
88pub(crate) fn dyadic_primitive_root(s_b: &Natural, t_b: i64) -> (i64, Natural, u64) {
89 if *s_b == 1u32 {
90 // base = 2^t_b with t_b != 0: the primitive root is 2 (or 1/2 for a base below 1).
91 return if t_b > 0 {
92 (1, Natural::ONE, u64::exact_from(t_b))
93 } else {
94 (-1, Natural::ONE, u64::exact_from(-t_b))
95 };
96 }
97 // s_b = h0^k with k maximal (k = 1 when s_b is not a perfect power).
98 let (h0, k) = s_b.express_as_power().unwrap_or_else(|| (s_b.clone(), 1));
99 // base = h0^k * 2^t_b = (2^(t_b / e) * h0^(k / e))^e for any common divisor e of k and t_b; the
100 // primitive root takes e maximal.
101 let e = if t_b == 0 {
102 k
103 } else {
104 k.gcd(t_b.unsigned_abs())
105 };
106 (t_b / i64::exact_from(e), h0.pow(k / e), e)
107}
108
109// Given a positive `Rational` `g` other than 1, decomposes it as `g = 2^z * hn / hd` with `hn` and
110// `hd` odd (coprime) `Natural`s. At most one of the numerator and denominator is even, since they
111// are coprime.
112pub(crate) fn rational_root_parts(g: &Rational) -> (i64, Natural, Natural) {
113 let num = g.numerator_ref();
114 let den = g.denominator_ref();
115 let num_z = num.trailing_zeros().unwrap();
116 let den_z = den.trailing_zeros().unwrap();
117 (
118 i64::exact_from(num_z) - i64::exact_from(den_z),
119 num >> num_z,
120 den >> den_z,
121 )
122}
123
124// Given a positive dyadic value `x = s * 2^t` (`s` odd) and a `Rational` root `g != 1` (`g > 0`),
125// returns the integer `m` with `x = g^m`, or `None` if no such integer exists. Writing `g = 2^z *
126// hn / hd` with `hn`, `hd` odd and coprime: a positive `m` requires `hd = 1` (an odd denominator
127// could never cancel against the dyadic `x`), and a negative `m` symmetrically requires `hn = 1`,
128// in which case `x = (2^(-z) * hd)^(-m)`.
129pub(crate) fn dyadic_log_of_rational_root(s: &Natural, t: i64, g: &Rational) -> Option<i64> {
130 let (z, hn, hd) = rational_root_parts(g);
131 if hd == 1u32 {
132 dyadic_log_of_root(s, t, z, &hn)
133 } else if hn == 1u32 {
134 dyadic_log_of_root(s, t, -z, &hd).map(|m| -m)
135 } else {
136 // Both an odd numerator and an odd denominator: no nonzero power is dyadic.
137 None
138 }
139}
140
141// Given a positive `Rational` `x` and a root `2^z * h` (`h` an odd `Natural`, the root positive and
142// different from 1), returns the integer `m` with `x = (2^z * h)^m`, or `None` if no such integer
143// exists. All big-number work is on `x`'s odd numerator and denominator parts and on `h`; the
144// power-of-2 parts stay as `i64` exponents, so the root's 2-power is never materialized.
145pub(crate) fn rational_value_log_of_dyadic_root(x: &Rational, z: i64, h: &Natural) -> Option<i64> {
146 let num = x.numerator_ref();
147 let den = x.denominator_ref();
148 let num_z = num.trailing_zeros().unwrap();
149 let den_z = den.trailing_zeros().unwrap();
150 let v2 = i64::exact_from(num_z) - i64::exact_from(den_z);
151 let on = num >> num_z;
152 let od = den >> den_z;
153 if *h == 1u32 {
154 // The root is 2^z: x = 2^(z * m) requires odd parts 1 and z | v2.
155 return if on == 1u32 && od == 1u32 && z != 0 && v2.divisible_by(z) {
156 Some(v2 / z)
157 } else {
158 None
159 };
160 }
161 // h >= 3: a positive power puts h^m in the numerator, a negative one in the denominator.
162 if on != 1u32 && od == 1u32 {
163 let m = i64::exact_from(on.checked_log_base(h)?);
164 if z.checked_mul(m)? == v2 {
165 Some(m)
166 } else {
167 None
168 }
169 } else if on == 1u32 && od != 1u32 {
170 let m = i64::exact_from(od.checked_log_base(h)?);
171 if z.checked_mul(m)? == -v2 {
172 Some(-m)
173 } else {
174 None
175 }
176 } else {
177 // on = od = 1 means x is a power of 2, requiring m = 0 and hence x = 1 (callers handle);
178 // on, od > 1 cannot both come from a single power of h.
179 None
180 }
181}
182
183// A large prime modulus for the congruence filter in `dyadic_1p_log_of_root`: 2^64 - 59, the
184// largest 64-bit prime.
185const FILTER_PRIME: u64 = 0xFFFFFFFFFFFFFFC5;
186
187// Whether `n = h^m`, where `n` is only available as the implicit sum `high * 2^shift + low` (with
188// `shift` possibly enormous). A congruence modulo `FILTER_PRIME` proves inequality cheaply -- the
189// power and the shift reduce via modular exponentiation -- and only a match (in practice a genuine
190// power) is verified exactly, at cost proportional to `shift`.
191fn implicit_sum_is_pow(high: &Natural, shift: u64, low: &Integer, h: &Natural, m: u64) -> bool {
192 let p = Natural::from(FILTER_PRIME);
193 let lhs = (high % &p)
194 .mod_mul(Natural::TWO.mod_pow(Natural::from(shift), &p), &p)
195 .mod_add(Natural::exact_from(low.mod_op(Integer::from(&p))), &p);
196 if (h % &p).mod_pow(Natural::from(m), &p) != lhs {
197 return false;
198 }
199 // The filter passed; verify exactly.
200 Integer::from(high << shift) + low == h.pow(m)
201}
202
203// Given a `Float` `x` (finite, nonzero, greater than -1) and a root `2^z * h` (`h` an odd
204// `Natural`, the root positive and different from 1), returns the integer `m` with `1 + x = (2^z *
205// h)^m`, or `None` if no such integer exists.
206//
207// `1 + x` is never materialized up front: its bit length can be as large as `|EXP(x)|` (up to
208// ~2^30) even when `x` itself has few bits. Instead, the structure of `1 + x` -- an implicit sum
209// with an odd significand pinned by `x`'s own odd significand `s` and exponent `t` -- determines at
210// most a couple of candidate exponents `m`, each checked by a congruence filter that proves
211// non-powers unequal; the expensive exact verification runs only for a match, whose cost is
212// proportional to the true size of `1 + x`.
213pub(crate) fn dyadic_1p_log_of_root(x: &Float, z: i64, h: &Natural) -> Option<i64> {
214 let (s, t) = odd_significand_and_exponent(x);
215 let neg = *x < 0u32;
216 if t >= 0 {
217 // x is an integer; x > -1 and x != 0, so x >= 1 and 1 + x = s * 2^t + 1.
218 debug_assert!(!neg);
219 if t == 0 {
220 // 1 + x = s + 1 is small: decompose it and match directly.
221 let n = &s + Natural::ONE;
222 let r = n.trailing_zeros().unwrap();
223 return dyadic_log_of_root(&(n >> r), i64::exact_from(r), z, h);
224 }
225 // t > 0: 1 + x = s * 2^t + 1 is odd and greater than 1, so the root's power-of-2 part must
226 // vanish and h contributes a positive power.
227 if z != 0 || *h == 1u32 {
228 return None;
229 }
230 // 1 + x has bit length exactly t + bits(s). Each candidate m must reproduce it.
231 let l = u64::exact_from(t) + s.significant_bits();
232 for m in pow_bit_length_candidates(h, l) {
233 if implicit_sum_is_pow(&s, u64::exact_from(t), &Integer::ONE, h, m) {
234 return Some(i64::exact_from(m));
235 }
236 }
237 return None;
238 }
239 // t < 0: 1 + x = (2^|t| ± s) * 2^t, with an odd numerator n = 2^|t| - s (x < 0) or 2^|t| + s
240 // (x > 0). The 2-adic valuation of 1 + x is exactly t, so z * m = t.
241 let t_abs = u64::exact_from(-t);
242 if *h == 1u32 {
243 // Root 2^z: 1 + x = 2^(z * m) requires n = 1, i.e. x < 0 and s = 2^|t| - 1 (an all-ones odd
244 // number, checked without materializing 2^|t|).
245 return if neg
246 && z != 0
247 && t.divisible_by(z)
248 && s.significant_bits() == t_abs
249 && (&s + Natural::ONE).is_power_of_2()
250 {
251 Some(t / z)
252 } else {
253 None
254 };
255 }
256 if z == 0 || !t.divisible_by(z) {
257 // An odd root can't produce the nonzero 2-adic valuation t, and a mixed root needs z | t.
258 return None;
259 }
260 // The candidate exponent is pinned by the 2-adic valuation alone: n = h^m needs m >= 1, since n
261 // is a positive integer and n = 1 (forcing m = 0) was the h = 1 case above. The congruence
262 // filter then proves or refutes n = h^m without materializing n.
263 let m = t / z;
264 if m <= 0 {
265 return None;
266 }
267 let low = if neg {
268 -Integer::from(&s)
269 } else {
270 Integer::from(&s)
271 };
272 if implicit_sum_is_pow(&Natural::ONE, t_abs, &low, h, u64::exact_from(m)) {
273 Some(m)
274 } else {
275 None
276 }
277}
278
279// The integers `m` for which `h^m` (with `h` an odd `Natural` at least 3) can have bit length
280// exactly `l`: at most a couple of values, pinned by 80-bit directed bounds on `log_2(h)`.
281fn pow_bit_length_candidates(h: &Natural, l: u64) -> Vec<u64> {
282 // log_2(h) is irrational (h is odd and at least 3), so directed rounding gives strict bounds.
283 let h_float = Float::exact_from(h.clone());
284 let lo = Rational::exact_from(h_float.log_base_2_prec_round_ref(80, Floor).0);
285 let hi = Rational::exact_from(h_float.log_base_2_prec_round_ref(80, Ceiling).0);
286 // h^m has bit length l iff l - 1 <= m * log_2(h) < l.
287 let m_min = Integer::rounding_from(Rational::from(l - 1) / hi, Ceiling).0;
288 let m_max = Integer::rounding_from(Rational::from(l) / lo, Floor).0;
289 let mut candidates = Vec::new();
290 let mut m = m_min;
291 while m <= m_max && candidates.len() < 4 {
292 if m > 0u32 {
293 candidates.push(u64::exact_from(&m));
294 }
295 m += Integer::ONE;
296 }
297 candidates
298}
299
300// `log_base(x)` is rational exactly when `x` and `base` are both powers of a common root `g`, say
301// `x = g ^ m` and `base = g ^ e_base`; then `log_base(x) = m / e_base`. Taking `g` to be the
302// smallest integer of which `base` is a power (obtained by stripping `base` of perfect-power
303// factors via `express_as_power`) and writing `g = 2^z * h` with `h` odd, this holds iff `x`'s odd
304// significand is the corresponding power of `h` and its exponent matches (see
305// `dyadic_log_of_root`).
306//
307// Detecting these rational results up front is essential, not just an optimization: when the result
308// is exactly representable (for example `log_9(3) = 1/2`), the Ziv loop in
309// `log_base_prec_round_normal` would never terminate, because the rounding test can never certify a
310// value that sits exactly on a representable point (or exactly on a tie). This generalizes the
311// `10^n` exactness check in mpfr_log10, which only catches integer results.
312//
313// The check is complete and cheap for any input: representable results can have enormous exponents
314// with few significant bits (`log_4` of the smallest positive `Float` is `-2^29`, exact at
315// precision 1), so no size cutoff on `x`'s exponent is sound; instead the decomposition keeps all
316// big-number work on `x`'s odd significand, which has at most `prec(x)` bits.
317pub(crate) fn rational_log_base(x: &Float, base: u64) -> Option<Rational> {
318 // `express_as_power` returns `None` when `base` is not a perfect power, in which case `base`
319 // itself is `g` (with exponent 1).
320 let (g, e_base) = base.express_as_power().unwrap_or((base, 1));
321 let z = i64::exact_from(g.trailing_zeros());
322 let h = Natural::from(g >> g.trailing_zeros());
323 let (s, t) = odd_significand_and_exponent(x);
324 let m = dyadic_log_of_root(&s, t, z, &h)?;
325 Some(Rational::from_signeds(m, i64::exact_from(e_base)))
326}
327
328// Returns `Some(m / e_base)` -- the value of `log_base(x)` -- when the positive `Rational` `x`
329// equals `g ^ m` for the root `g` of `base` (so `base = g ^ e_base` and `log_base(x)` is rational),
330// and `None` when `log_base(x)` is irrational. `x` must be positive and `base > 1`.
331//
332// `m` (signed) is found by `Rational::checked_log_base`, which also covers `x < 1` (negative `m`).
333// Detecting these rational results up front is essential: the Ziv loop could never certify an
334// exactly-representable one (see `rational_log_base` for the `Float` analog).
335pub(crate) fn rational_log_base_of_rational(x: &Rational, base: u64) -> Option<Rational> {
336 let (g, e_base) = base.express_as_power().unwrap_or((base, 1));
337 x.checked_log_base(g)
338 .map(|m| Rational::from_signeds(m, i64::exact_from(e_base)))
339}
340
341// The computation of log_base(x, base) is done by log_base(x) = ln(x) / ln(base). When `base` is a
342// power of 2 the caller delegates to `log_base_power_of_2`, so here `base` is not a power of 2.
343//
344// This is mpfr_log10 from log10.c, MPFR 4.3.0, generalized from base 10 to an arbitrary non-power-
345// of-2 `base`. The input is finite, nonzero, and positive.
346fn log_base_prec_round_normal(
347 x: &Float,
348 base: u64,
349 prec: u64,
350 rm: RoundingMode,
351) -> (Float, Ordering) {
352 // If x is 1, the result is 0.
353 if *x == 1u32 {
354 return (Float::ZERO, Equal);
355 }
356 // If log_base(x) is rational -- x and base are both powers of a common integer -- compute it
357 // directly. This includes the exactly-representable results (integers like log_8(64) = 2 and
358 // dyadics like log_9(3) = 1/2), which the Ziv loop below could never certify, as well as
359 // non-representable rationals like log_27(9) = 2/3, which it could but for which the direct
360 // computation is cheaper and exact.
361 if let Some(q) = rational_log_base(x, base) {
362 return Float::from_rational_prec_round(q, prec, rm);
363 }
364 // log_base(x) for x in a sliver of 1 can fall below the smallest positive Float; the 1-plus-x
365 // form handles that underflow region.
366 match sliver_of_one(x) {
367 SliverOfOne::Representable(d) => return d.log_base_1_plus_x_prec_round(base, prec, rm),
368 SliverOfOne::Underflow => {
369 return Float::log_base_rational_prec_round(Rational::exact_from(x), base, prec, rm);
370 }
371 SliverOfOne::No => {}
372 }
373 // The result is irrational, so it is never exactly representable.
374 assert_ne!(rm, Exact, "Inexact log_base");
375 let base_float = Float::from(base);
376 // Compute the precision of the intermediary variable: the optimal number of bits, see
377 // algorithms.tex.
378 let mut working_prec = prec + 4 + prec.ceiling_log_base_2();
379 let mut increment = Limb::WIDTH;
380 loop {
381 // ln(x) / ln(base). ln(x), ln(base), and the division are each correctly rounded (at most
382 // 1/2 ulp), so the relative error is below 2^(2 - working_prec) and working_prec - 4
383 // correct bits suffice for rounding (mpfr_log10 uses Nt - 4).
384 let t = x
385 .ln_prec_ref(working_prec)
386 .0
387 .div_prec(base_float.ln_prec_ref(working_prec).0, working_prec)
388 .0;
389 if float_can_round(t.significand_ref().unwrap(), working_prec - 4, prec, rm) {
390 return Float::from_float_prec_round(t, prec, rm);
391 }
392 // Increase the precision.
393 working_prec += increment;
394 increment = working_prec >> 1;
395 }
396}
397
398// Computes log_base(x) for a positive `Rational` x whose logarithm is irrational, in a Ziv loop.
399// `base > 1` is not a power of 2.
400//
401// log_base(x) = log_2(x) / log_2(base). Routing through `log_base_2_rational` (rather than
402// computing `ln(x) / ln(base)` directly) reuses its handling of x near a power of 2 -- in
403// particular x near 1, where the result is near 0 and a direct computation would need a working
404// precision proportional to how close x is to 1. log_2(x), log_2(base), and the division are each
405// correctly rounded (at most 1/2 ulp), so the relative error is below 2^(2 - working_prec) and
406// working_prec - 4 correct bits suffice for rounding.
407fn log_base_rational_prec_round_helper(
408 x: &Rational,
409 base: u64,
410 prec: u64,
411 rm: RoundingMode,
412) -> (Float, Ordering) {
413 let base_float = Float::from(base);
414 // The initial slack keeps working_prec at least 7, so the working_prec - 6 below stays
415 // positive.
416 let mut working_prec = prec + 6 + prec.ceiling_log_base_2();
417 let mut increment = Limb::WIDTH;
418 loop {
419 // log_2(x) in the extended exponent range: for x within a sliver of 1 the ordinary Float
420 // form would flush to zero or clamp, and the rounding test below could never resolve it.
421 let num = extended_log_base_2_of_rational(x, working_prec);
422 let den = ExtendedFloat::from(base_float.log_base_2_prec_ref(working_prec).0);
423 let quotient = num.div_prec_val_ref(&den, working_prec).0;
424 // log_2(x) is within 2 ulps, log_2(base) within 1/2, and the division adds 1/2 more, so
425 // working_prec - 6 correct bits comfortably suffice.
426 if float_can_round(
427 quotient.x.significand_ref().unwrap(),
428 working_prec - 6,
429 prec,
430 rm,
431 ) {
432 let (rounded, o) = Float::from_float_prec_round(quotient.x, prec, rm);
433 let mut result = ExtendedFloat::from(rounded);
434 result.exp = result.exp.checked_add(quotient.exp).unwrap();
435 return result.into_float_helper(prec, rm, o);
436 }
437 // Increase the precision.
438 working_prec += increment;
439 increment = working_prec >> 1;
440 }
441}
442
443impl Float {
444 /// Computes $\log_b x$, where $x$ is a [`Float`] and $b$ is a `u64` greater than 1, rounding
445 /// the result to the specified precision and with the specified rounding mode. The [`Float`] is
446 /// taken by value. An [`Ordering`] is also returned, indicating whether the rounded value is
447 /// less than, equal to, or greater than the exact value. Although `NaN`s are not comparable to
448 /// any [`Float`], whenever this function returns a `NaN` it also returns `Equal`.
449 ///
450 /// The base-$b$ logarithm of any nonzero negative number is `NaN`.
451 ///
452 /// When `base` is a power of 2, this function delegates to
453 /// [`Float::log_base_power_of_2_prec_round`]; otherwise it computes $\ln x / \ln b$.
454 ///
455 /// See [`RoundingMode`] for a description of the possible rounding modes.
456 ///
457 /// $$
458 /// f(x,b,p,m) = \log_b x+\varepsilon.
459 /// $$
460 /// - If $\log_b x$ is infinite, zero, or `NaN`, $\varepsilon$ may be ignored or assumed to be
461 /// 0.
462 /// - If $\log_b x$ is finite and nonzero, and $m$ is not `Nearest`, then $|\varepsilon| <
463 /// 2^{\lfloor\log_2 |\log_b x|\rfloor-p+1}$.
464 /// - If $\log_b x$ is finite and nonzero, and $m$ is `Nearest`, then $|\varepsilon| \leq
465 /// 2^{\lfloor\log_2 |\log_b x|\rfloor-p}$.
466 ///
467 /// If the output has a precision, it is `prec`.
468 ///
469 /// Special cases:
470 /// - $f(\text{NaN},b,p,m)=\text{NaN}$
471 /// - $f(\infty,b,p,m)=\infty$
472 /// - $f(-\infty,b,p,m)=\text{NaN}$
473 /// - $f(\pm0.0,b,p,m)=-\infty$
474 /// - $f(1.0,b,p,m)=0.0$, and the result is exact
475 /// - $f(b^n,b,p,m)=n$, rounded to precision $p$; the result is exact if and only if $n$ is
476 /// representable with precision $p$
477 /// - $f(x,b,p,m)=\text{NaN}$ for $x<0$
478 ///
479 /// Neither overflow nor underflow is possible.
480 ///
481 /// If you know you'll be using `Nearest`, consider using [`Float::log_base_prec`] instead. If
482 /// you know that your target precision is the precision of the input, consider using
483 /// [`Float::log_base_round`] instead. If both of these things are true, consider using
484 /// [`Float::log_base`] instead.
485 ///
486 /// # Worst-case complexity
487 /// $T(n, m) = O(n (\log n)^2 \log\log n + m \log m \log\log m)$
488 ///
489 /// $M(n, m) = O(n \log n + m \log m)$
490 ///
491 /// where $T$ is time, $M$ is additional memory, $n$ is `prec`, and $m$ is
492 /// `self.significant_bits()`.
493 ///
494 /// # Panics
495 /// Panics if `prec` is zero, if `base` is less than 2, or if `rm` is `Exact` but the result
496 /// cannot be represented exactly with the given precision.
497 ///
498 /// # Examples
499 /// ```
500 /// use malachite_base::rounding_modes::RoundingMode::*;
501 /// use malachite_float::Float;
502 /// use std::cmp::Ordering::*;
503 ///
504 /// let (log, o) = Float::from(1000).log_base_prec_round(10, 10, Nearest);
505 /// assert_eq!(log.to_string(), "3.0000");
506 /// assert_eq!(o, Equal);
507 ///
508 /// let (log, o) = Float::from(50).log_base_prec_round(10, 10, Floor);
509 /// assert_eq!(log.to_string(), "1.6973");
510 /// assert_eq!(o, Less);
511 ///
512 /// let (log, o) = Float::from(50).log_base_prec_round(10, 10, Ceiling);
513 /// assert_eq!(log.to_string(), "1.6992");
514 /// assert_eq!(o, Greater);
515 /// ```
516 #[inline]
517 pub fn log_base_prec_round(self, base: u64, prec: u64, rm: RoundingMode) -> (Self, Ordering) {
518 assert_ne!(prec, 0);
519 assert!(base > 1, "Logarithm base must be greater than 1");
520 if base.is_power_of_2() {
521 return self.log_base_power_of_2_prec_round(i64::from(base.trailing_zeros()), prec, rm);
522 }
523 match self {
524 Self(NaN | Infinity { sign: false } | Finite { sign: false, .. }) => {
525 (float_nan!(), Equal)
526 }
527 float_either_zero!() => (float_negative_infinity!(), Equal),
528 float_infinity!() => (float_infinity!(), Equal),
529 _ => log_base_prec_round_normal(&self, base, prec, rm),
530 }
531 }
532
533 /// Computes $\log_b x$, where $x$ is a [`Float`] and $b$ is a `u64` greater than 1, rounding
534 /// the result to the specified precision and with the specified rounding mode. The [`Float`] is
535 /// taken by reference. An [`Ordering`] is also returned, indicating whether the rounded value
536 /// is less than, equal to, or greater than the exact value. Although `NaN`s are not comparable
537 /// to any [`Float`], whenever this function returns a `NaN` it also returns `Equal`.
538 ///
539 /// See [`Float::log_base_prec_round`] for details, special cases, and a description of the
540 /// rounding behavior.
541 ///
542 /// # Worst-case complexity
543 /// $T(n, m) = O(n (\log n)^2 \log\log n + m \log m \log\log m)$
544 ///
545 /// $M(n, m) = O(n \log n + m \log m)$
546 ///
547 /// where $T$ is time, $M$ is additional memory, $n$ is `prec`, and $m$ is
548 /// `self.significant_bits()`.
549 ///
550 /// # Panics
551 /// Panics if `prec` is zero, if `base` is less than 2, or if `rm` is `Exact` but the result
552 /// cannot be represented exactly with the given precision.
553 ///
554 /// # Examples
555 /// ```
556 /// use malachite_base::rounding_modes::RoundingMode::*;
557 /// use malachite_float::Float;
558 /// use std::cmp::Ordering::*;
559 ///
560 /// let (log, o) = Float::from(1000).log_base_prec_round_ref(10, 10, Nearest);
561 /// assert_eq!(log.to_string(), "3.0000");
562 /// assert_eq!(o, Equal);
563 /// ```
564 #[inline]
565 pub fn log_base_prec_round_ref(
566 &self,
567 base: u64,
568 prec: u64,
569 rm: RoundingMode,
570 ) -> (Self, Ordering) {
571 assert_ne!(prec, 0);
572 assert!(base > 1, "Logarithm base must be greater than 1");
573 if base.is_power_of_2() {
574 return self.log_base_power_of_2_prec_round_ref(
575 i64::from(base.trailing_zeros()),
576 prec,
577 rm,
578 );
579 }
580 match self {
581 Self(NaN | Infinity { sign: false } | Finite { sign: false, .. }) => {
582 (float_nan!(), Equal)
583 }
584 float_either_zero!() => (float_negative_infinity!(), Equal),
585 float_infinity!() => (float_infinity!(), Equal),
586 _ => log_base_prec_round_normal(self, base, prec, rm),
587 }
588 }
589
590 /// Computes $\log_b x$, where $x$ is a [`Float`] and $b$ is a `u64` greater than 1, rounding
591 /// the result to the nearest value of the specified precision. The [`Float`] is taken by value.
592 /// An [`Ordering`] is also returned, indicating whether the rounded value is less than, equal
593 /// to, or greater than the exact value.
594 ///
595 /// See [`Float::log_base_prec_round`] for details and special cases.
596 ///
597 /// # Worst-case complexity
598 /// $T(n, m) = O(n (\log n)^2 \log\log n + m \log m \log\log m)$
599 ///
600 /// $M(n, m) = O(n \log n + m \log m)$
601 ///
602 /// where $T$ is time, $M$ is additional memory, $n$ is `prec`, and $m$ is
603 /// `self.significant_bits()`.
604 ///
605 /// # Panics
606 /// Panics if `prec` is zero or if `base` is less than 2.
607 ///
608 /// # Examples
609 /// ```
610 /// use malachite_float::Float;
611 /// use std::cmp::Ordering::*;
612 ///
613 /// let (log, o) = Float::from(50).log_base_prec(10, 10);
614 /// assert_eq!(log.to_string(), "1.6992");
615 /// assert_eq!(o, Greater);
616 /// ```
617 #[inline]
618 pub fn log_base_prec(self, base: u64, prec: u64) -> (Self, Ordering) {
619 self.log_base_prec_round(base, prec, Nearest)
620 }
621
622 /// Computes $\log_b x$, where $x$ is a [`Float`] and $b$ is a `u64` greater than 1, rounding
623 /// the result to the nearest value of the specified precision. The [`Float`] is taken by
624 /// reference. An [`Ordering`] is also returned, indicating whether the rounded value is less
625 /// than, equal to, or greater than the exact value.
626 ///
627 /// See [`Float::log_base_prec_round`] for details and special cases.
628 ///
629 /// # Worst-case complexity
630 /// $T(n, m) = O(n (\log n)^2 \log\log n + m \log m \log\log m)$
631 ///
632 /// $M(n, m) = O(n \log n + m \log m)$
633 ///
634 /// where $T$ is time, $M$ is additional memory, $n$ is `prec`, and $m$ is
635 /// `self.significant_bits()`.
636 ///
637 /// # Panics
638 /// Panics if `prec` is zero or if `base` is less than 2.
639 ///
640 /// # Examples
641 /// ```
642 /// use malachite_float::Float;
643 /// use std::cmp::Ordering::*;
644 ///
645 /// let (log, o) = Float::from(50).log_base_prec_ref(10, 10);
646 /// assert_eq!(log.to_string(), "1.6992");
647 /// assert_eq!(o, Greater);
648 /// ```
649 #[inline]
650 pub fn log_base_prec_ref(&self, base: u64, prec: u64) -> (Self, Ordering) {
651 self.log_base_prec_round_ref(base, prec, Nearest)
652 }
653
654 /// Computes $\log_b x$, where $x$ is a [`Float`] and $b$ is a `u64` greater than 1, rounding
655 /// the result to the precision of the input and with the specified rounding mode. The [`Float`]
656 /// is taken by value. An [`Ordering`] is also returned, indicating whether the rounded value is
657 /// less than, equal to, or greater than the exact value.
658 ///
659 /// See [`Float::log_base_prec_round`] for details and special cases.
660 ///
661 /// # Worst-case complexity
662 /// $T(n) = O(n (\log n)^2 \log\log n)$
663 ///
664 /// $M(n) = O(n \log n)$
665 ///
666 /// where $T$ is time, $M$ is additional memory, and $n$ is the precision of the input.
667 ///
668 /// # Panics
669 /// Panics if `base` is less than 2, or if `rm` is `Exact` but the result cannot be represented
670 /// exactly with the input's precision.
671 ///
672 /// # Examples
673 /// ```
674 /// use malachite_base::rounding_modes::RoundingMode::*;
675 /// use malachite_float::Float;
676 /// use std::cmp::Ordering::*;
677 ///
678 /// let (log, o) = Float::from(1000).log_base_round(10, Floor);
679 /// assert_eq!(log.to_string(), "3.000");
680 /// assert_eq!(o, Equal);
681 /// ```
682 #[inline]
683 pub fn log_base_round(self, base: u64, rm: RoundingMode) -> (Self, Ordering) {
684 let prec = self.significant_bits();
685 self.log_base_prec_round(base, prec, rm)
686 }
687
688 /// Computes $\log_b x$, where $x$ is a [`Float`] and $b$ is a `u64` greater than 1, rounding
689 /// the result to the precision of the input and with the specified rounding mode. The [`Float`]
690 /// is taken by reference. An [`Ordering`] is also returned, indicating whether the rounded
691 /// value is less than, equal to, or greater than the exact value.
692 ///
693 /// See [`Float::log_base_prec_round`] for details and special cases.
694 ///
695 /// # Worst-case complexity
696 /// $T(n) = O(n (\log n)^2 \log\log n)$
697 ///
698 /// $M(n) = O(n \log n)$
699 ///
700 /// where $T$ is time, $M$ is additional memory, and $n$ is the precision of the input.
701 ///
702 /// # Panics
703 /// Panics if `base` is less than 2, or if `rm` is `Exact` but the result cannot be represented
704 /// exactly with the input's precision.
705 ///
706 /// # Examples
707 /// ```
708 /// use malachite_base::rounding_modes::RoundingMode::*;
709 /// use malachite_float::Float;
710 /// use std::cmp::Ordering::*;
711 ///
712 /// let (log, o) = Float::from(81).log_base_round_ref(3, Ceiling);
713 /// assert_eq!(log.to_string(), "4.000");
714 /// assert_eq!(o, Equal);
715 /// ```
716 #[inline]
717 pub fn log_base_round_ref(&self, base: u64, rm: RoundingMode) -> (Self, Ordering) {
718 self.log_base_prec_round_ref(base, self.significant_bits(), rm)
719 }
720
721 /// Computes $\log_b x$, where $x$ is a [`Float`] and $b$ is a `u64` greater than 1, in place,
722 /// rounding the result to the specified precision and with the specified rounding mode. An
723 /// [`Ordering`] is returned, indicating whether the rounded value is less than, equal to, or
724 /// greater than the exact value.
725 ///
726 /// See [`Float::log_base_prec_round`] for details and special cases.
727 ///
728 /// # Worst-case complexity
729 /// $T(n, m) = O(n (\log n)^2 \log\log n + m \log m \log\log m)$
730 ///
731 /// $M(n, m) = O(n \log n + m \log m)$
732 ///
733 /// where $T$ is time, $M$ is additional memory, $n$ is `prec`, and $m$ is
734 /// `self.significant_bits()`.
735 ///
736 /// # Panics
737 /// Panics if `prec` is zero, if `base` is less than 2, or if `rm` is `Exact` but the result
738 /// cannot be represented exactly with the given precision.
739 ///
740 /// # Examples
741 /// ```
742 /// use malachite_base::rounding_modes::RoundingMode::*;
743 /// use malachite_float::Float;
744 /// use std::cmp::Ordering::*;
745 ///
746 /// let mut x = Float::from(50);
747 /// let o = x.log_base_prec_round_assign(10, 10, Floor);
748 /// assert_eq!(x.to_string(), "1.6973");
749 /// assert_eq!(o, Less);
750 /// ```
751 #[inline]
752 pub fn log_base_prec_round_assign(
753 &mut self,
754 base: u64,
755 prec: u64,
756 rm: RoundingMode,
757 ) -> Ordering {
758 let (result, o) = core::mem::take(self).log_base_prec_round(base, prec, rm);
759 *self = result;
760 o
761 }
762
763 /// Computes $\log_b x$, where $x$ is a [`Float`] and $b$ is a `u64` greater than 1, in place,
764 /// rounding the result to the nearest value of the specified precision. An [`Ordering`] is
765 /// returned, indicating whether the rounded value is less than, equal to, or greater than the
766 /// exact value.
767 ///
768 /// See [`Float::log_base_prec_round`] for details and special cases.
769 ///
770 /// # Worst-case complexity
771 /// $T(n, m) = O(n (\log n)^2 \log\log n + m \log m \log\log m)$
772 ///
773 /// $M(n, m) = O(n \log n + m \log m)$
774 ///
775 /// where $T$ is time, $M$ is additional memory, $n$ is `prec`, and $m$ is
776 /// `self.significant_bits()`.
777 ///
778 /// # Panics
779 /// Panics if `prec` is zero or if `base` is less than 2.
780 ///
781 /// # Examples
782 /// ```
783 /// use malachite_float::Float;
784 /// use std::cmp::Ordering::*;
785 ///
786 /// let mut x = Float::from(1000);
787 /// let o = x.log_base_prec_assign(10, 10);
788 /// assert_eq!(x.to_string(), "3.0000");
789 /// assert_eq!(o, Equal);
790 /// ```
791 #[inline]
792 pub fn log_base_prec_assign(&mut self, base: u64, prec: u64) -> Ordering {
793 self.log_base_prec_round_assign(base, prec, Nearest)
794 }
795
796 /// Computes $\log_b x$, where $x$ is a [`Float`] and $b$ is a `u64` greater than 1, in place,
797 /// rounding the result to the precision of the input and with the specified rounding mode. An
798 /// [`Ordering`] is returned, indicating whether the rounded value is less than, equal to, or
799 /// greater than the exact value.
800 ///
801 /// See [`Float::log_base_prec_round`] for details and special cases.
802 ///
803 /// # Worst-case complexity
804 /// $T(n) = O(n (\log n)^2 \log\log n)$
805 ///
806 /// $M(n) = O(n \log n)$
807 ///
808 /// where $T$ is time, $M$ is additional memory, and $n$ is the precision of the input.
809 ///
810 /// # Panics
811 /// Panics if `base` is less than 2, or if `rm` is `Exact` but the result cannot be represented
812 /// exactly with the input's precision.
813 ///
814 /// # Examples
815 /// ```
816 /// use malachite_base::rounding_modes::RoundingMode::*;
817 /// use malachite_float::Float;
818 /// use std::cmp::Ordering::*;
819 ///
820 /// let mut x = Float::from(81);
821 /// let o = x.log_base_round_assign(3, Nearest);
822 /// assert_eq!(x.to_string(), "4.000");
823 /// assert_eq!(o, Equal);
824 /// ```
825 #[inline]
826 pub fn log_base_round_assign(&mut self, base: u64, rm: RoundingMode) -> Ordering {
827 let prec = self.significant_bits();
828 self.log_base_prec_round_assign(base, prec, rm)
829 }
830
831 /// Computes $\log_b x$, where $x$ is a [`Rational`] and $b$ is a `u64` greater than 1, rounding
832 /// the result to the specified precision and with the specified rounding mode and returning the
833 /// result as a [`Float`]. The [`Rational`] is taken by value. An [`Ordering`] is also returned,
834 /// indicating whether the rounded value is less than, equal to, or greater than the exact
835 /// value. Although `NaN`s are not comparable to any [`Float`], whenever this function returns a
836 /// `NaN` it also returns `Equal`.
837 ///
838 /// The base-$b$ logarithm of any negative number is `NaN`.
839 ///
840 /// Inputs of any magnitude are handled, including [`Rational`]s whose magnitudes are too large
841 /// or too small to be representable as [`Float`]s. Neither overflow nor underflow of the output
842 /// is possible.
843 ///
844 /// When `base` is a power of 2, this function delegates to
845 /// [`Float::log_base_power_of_2_rational_prec_round`].
846 ///
847 /// See [`Float::log_base_prec_round`] for details and a description of the rounding behavior.
848 ///
849 /// Special cases:
850 /// - $f(0,b,p,m)=-\infty$
851 /// - $f(x,b,p,m)=\text{NaN}$ for $x<0$
852 /// - $f(1,b,p,m)=0.0$, and the result is exact
853 ///
854 /// # Worst-case complexity
855 /// $T(n, m) = O(n (\log n)^2 \log\log n + m \log m \log\log m)$
856 ///
857 /// $M(n, m) = O(n \log n + m \log m)$
858 ///
859 /// where $T$ is time, $M$ is additional memory, $n$ is `prec`, and $m$ is
860 /// `x.significant_bits()`.
861 ///
862 /// # Panics
863 /// Panics if `prec` is zero, if `base` is less than 2, or if `rm` is `Exact` but the result
864 /// cannot be represented exactly with the given precision. (The result is exactly representable
865 /// if and only if $x \leq 0$ or $\log_b x$ is rational and representable with the given
866 /// precision.)
867 ///
868 /// # Examples
869 /// ```
870 /// use malachite_base::num::basic::traits::Two;
871 /// use malachite_base::rounding_modes::RoundingMode::*;
872 /// use malachite_float::Float;
873 /// use malachite_q::Rational;
874 /// use std::cmp::Ordering::*;
875 ///
876 /// let (log, o) = Float::log_base_rational_prec_round(Rational::from(3), 9, 10, Exact);
877 /// assert_eq!(log.to_string(), "0.50000"); // log_9(3) = 1/2
878 /// assert_eq!(o, Equal);
879 ///
880 /// let (log, o) = Float::log_base_rational_prec_round(Rational::TWO, 3, 20, Nearest);
881 /// assert_eq!(log.to_string(), "0.63092995");
882 /// assert_eq!(o, Greater);
883 /// ```
884 #[allow(clippy::needless_pass_by_value)]
885 #[inline]
886 pub fn log_base_rational_prec_round(
887 x: Rational,
888 base: u64,
889 prec: u64,
890 rm: RoundingMode,
891 ) -> (Self, Ordering) {
892 Self::log_base_rational_prec_round_ref(&x, base, prec, rm)
893 }
894
895 /// Computes $\log_b x$, where $x$ is a [`Rational`] and $b$ is a `u64` greater than 1, rounding
896 /// the result to the specified precision and with the specified rounding mode and returning the
897 /// result as a [`Float`]. The [`Rational`] is taken by reference. An [`Ordering`] is also
898 /// returned, indicating whether the rounded value is less than, equal to, or greater than the
899 /// exact value. Although `NaN`s are not comparable to any [`Float`], whenever this function
900 /// returns a `NaN` it also returns `Equal`.
901 ///
902 /// See [`Float::log_base_rational_prec_round`] for details, special cases, and a description of
903 /// the rounding behavior.
904 ///
905 /// # Worst-case complexity
906 /// $T(n, m) = O(n (\log n)^2 \log\log n + m \log m \log\log m)$
907 ///
908 /// $M(n, m) = O(n \log n + m \log m)$
909 ///
910 /// where $T$ is time, $M$ is additional memory, $n$ is `prec`, and $m$ is
911 /// `x.significant_bits()`.
912 ///
913 /// # Panics
914 /// Panics if `prec` is zero, if `base` is less than 2, or if `rm` is `Exact` but the result
915 /// cannot be represented exactly with the given precision.
916 ///
917 /// # Examples
918 /// ```
919 /// use malachite_base::rounding_modes::RoundingMode::*;
920 /// use malachite_float::Float;
921 /// use malachite_q::Rational;
922 /// use std::cmp::Ordering::*;
923 ///
924 /// let (log, o) =
925 /// Float::log_base_rational_prec_round_ref(&Rational::from_signeds(1, 9), 3, 10, Exact);
926 /// assert_eq!(log.to_string(), "-2.0000"); // log_3(1/9) = -2
927 /// assert_eq!(o, Equal);
928 /// ```
929 pub fn log_base_rational_prec_round_ref(
930 x: &Rational,
931 base: u64,
932 prec: u64,
933 rm: RoundingMode,
934 ) -> (Self, Ordering) {
935 assert_ne!(prec, 0);
936 assert!(base > 1, "Logarithm base must be greater than 1");
937 if base.is_power_of_2() {
938 return Self::log_base_power_of_2_rational_prec_round_ref(
939 x,
940 i64::from(base.trailing_zeros()),
941 prec,
942 rm,
943 );
944 }
945 match x.sign() {
946 Equal => return (float_negative_infinity!(), Equal),
947 Less => return (float_nan!(), Equal),
948 Greater => {}
949 }
950 // If x = g^m for the base's root g (so base = g^e_base), then log_base(x) = m / e_base is
951 // rational, and exact -- the Ziv loop could never certify it (see rational_log_base).
952 if let Some(q) = rational_log_base_of_rational(x, base) {
953 return Self::from_rational_prec_round(q, prec, rm);
954 }
955 // The result is irrational, so it is never exactly representable.
956 assert_ne!(rm, Exact, "Inexact log_base");
957 log_base_rational_prec_round_helper(x, base, prec, rm)
958 }
959
960 /// Computes $\log_b x$, where $x$ is a [`Rational`] and $b$ is a `u64` greater than 1, rounding
961 /// the result to the nearest value of the specified precision and returning the result as a
962 /// [`Float`]. The [`Rational`] is taken by value. An [`Ordering`] is also returned, indicating
963 /// whether the rounded value is less than, equal to, or greater than the exact value.
964 ///
965 /// See [`Float::log_base_rational_prec_round`] for details and special cases.
966 ///
967 /// # Worst-case complexity
968 /// $T(n, m) = O(n (\log n)^2 \log\log n + m \log m \log\log m)$
969 ///
970 /// $M(n, m) = O(n \log n + m \log m)$
971 ///
972 /// where $T$ is time, $M$ is additional memory, $n$ is `prec`, and $m$ is
973 /// `x.significant_bits()`.
974 ///
975 /// # Panics
976 /// Panics if `prec` is zero or if `base` is less than 2.
977 ///
978 /// # Examples
979 /// ```
980 /// use malachite_float::Float;
981 /// use malachite_q::Rational;
982 /// use std::cmp::Ordering::*;
983 ///
984 /// let (log, o) = Float::log_base_rational_prec(Rational::from_signeds(1, 9), 3, 10);
985 /// assert_eq!(log.to_string(), "-2.0000");
986 /// assert_eq!(o, Equal);
987 /// ```
988 #[inline]
989 pub fn log_base_rational_prec(x: Rational, base: u64, prec: u64) -> (Self, Ordering) {
990 Self::log_base_rational_prec_round(x, base, prec, Nearest)
991 }
992
993 /// Computes $\log_b x$, where $x$ is a [`Rational`] and $b$ is a `u64` greater than 1, rounding
994 /// the result to the nearest value of the specified precision and returning the result as a
995 /// [`Float`]. The [`Rational`] is taken by reference. An [`Ordering`] is also returned,
996 /// indicating whether the rounded value is less than, equal to, or greater than the exact
997 /// value.
998 ///
999 /// See [`Float::log_base_rational_prec_round`] for details and special cases.
1000 ///
1001 /// # Worst-case complexity
1002 /// $T(n, m) = O(n (\log n)^2 \log\log n + m \log m \log\log m)$
1003 ///
1004 /// $M(n, m) = O(n \log n + m \log m)$
1005 ///
1006 /// where $T$ is time, $M$ is additional memory, $n$ is `prec`, and $m$ is
1007 /// `x.significant_bits()`.
1008 ///
1009 /// # Panics
1010 /// Panics if `prec` is zero or if `base` is less than 2.
1011 ///
1012 /// # Examples
1013 /// ```
1014 /// use malachite_base::num::basic::traits::Two;
1015 /// use malachite_float::Float;
1016 /// use malachite_q::Rational;
1017 /// use std::cmp::Ordering::*;
1018 ///
1019 /// let (log, o) = Float::log_base_rational_prec_ref(&Rational::TWO, 3, 20);
1020 /// assert_eq!(log.to_string(), "0.63092995");
1021 /// assert_eq!(o, Greater);
1022 /// ```
1023 #[inline]
1024 pub fn log_base_rational_prec_ref(x: &Rational, base: u64, prec: u64) -> (Self, Ordering) {
1025 Self::log_base_rational_prec_round_ref(x, base, prec, Nearest)
1026 }
1027}
1028
1029impl LogBase<u64> for Float {
1030 type Output = Self;
1031
1032 /// Computes $\log_b x$, where $x$ is a [`Float`] and $b$ is a `u64` greater than 1, rounding
1033 /// the result to the nearest value of the input's precision. The [`Float`] is taken by value.
1034 ///
1035 /// The base-$b$ logarithm of any nonzero negative number is `NaN`. See
1036 /// [`Float::log_base_prec_round`] for the special cases.
1037 ///
1038 /// $$
1039 /// f(x,b) = \log_b x+\varepsilon,
1040 /// $$
1041 /// where $|\varepsilon| \leq 2^{\lfloor\log_2 |\log_b x|\rfloor-p}$ and $p$ is the precision of
1042 /// the input.
1043 ///
1044 /// # Worst-case complexity
1045 /// $T(n) = O(n (\log n)^2 \log\log n)$
1046 ///
1047 /// $M(n) = O(n \log n)$
1048 ///
1049 /// where $T$ is time, $M$ is additional memory, and $n$ is the precision of the input.
1050 ///
1051 /// # Panics
1052 /// Panics if `base` is less than 2.
1053 ///
1054 /// # Examples
1055 /// ```
1056 /// use malachite_base::num::arithmetic::traits::LogBase;
1057 /// use malachite_float::Float;
1058 ///
1059 /// assert_eq!(Float::from(1000).log_base(10).to_string(), "3.000");
1060 /// assert_eq!(Float::from(81).log_base(3).to_string(), "4.000");
1061 /// ```
1062 #[inline]
1063 fn log_base(self, base: u64) -> Self {
1064 let prec = self.significant_bits();
1065 self.log_base_prec_round(base, prec, Nearest).0
1066 }
1067}
1068
1069impl LogBase<u64> for &Float {
1070 type Output = Float;
1071
1072 /// Computes $\log_b x$, where $x$ is a [`Float`] and $b$ is a `u64` greater than 1, rounding
1073 /// the result to the nearest value of the input's precision. The [`Float`] is taken by
1074 /// reference.
1075 ///
1076 /// The base-$b$ logarithm of any nonzero negative number is `NaN`. See
1077 /// [`Float::log_base_prec_round`] for the special cases.
1078 ///
1079 /// $$
1080 /// f(x,b) = \log_b x+\varepsilon,
1081 /// $$
1082 /// where $|\varepsilon| \leq 2^{\lfloor\log_2 |\log_b x|\rfloor-p}$ and $p$ is the precision of
1083 /// the input.
1084 ///
1085 /// # Worst-case complexity
1086 /// $T(n) = O(n (\log n)^2 \log\log n)$
1087 ///
1088 /// $M(n) = O(n \log n)$
1089 ///
1090 /// where $T$ is time, $M$ is additional memory, and $n$ is the precision of the input.
1091 ///
1092 /// # Panics
1093 /// Panics if `base` is less than 2.
1094 ///
1095 /// # Examples
1096 /// ```
1097 /// use malachite_base::num::arithmetic::traits::LogBase;
1098 /// use malachite_float::Float;
1099 ///
1100 /// assert_eq!((&Float::from(1000)).log_base(10).to_string(), "3.000");
1101 /// ```
1102 #[inline]
1103 fn log_base(self, base: u64) -> Float {
1104 self.log_base_prec_round_ref(base, self.significant_bits(), Nearest)
1105 .0
1106 }
1107}
1108
1109impl LogBaseAssign<u64> for Float {
1110 /// Replaces a [`Float`] $x$ with $\log_b x$, where $b$ is a `u64` greater than 1, rounding the
1111 /// result to the nearest value of the input's precision.
1112 ///
1113 /// The base-$b$ logarithm of any nonzero negative number is `NaN`. See
1114 /// [`Float::log_base_prec_round`] for the special cases.
1115 ///
1116 /// # Worst-case complexity
1117 /// $T(n) = O(n (\log n)^2 \log\log n)$
1118 ///
1119 /// $M(n) = O(n \log n)$
1120 ///
1121 /// where $T$ is time, $M$ is additional memory, and $n$ is the precision of the input.
1122 ///
1123 /// # Panics
1124 /// Panics if `base` is less than 2.
1125 ///
1126 /// # Examples
1127 /// ```
1128 /// use malachite_base::num::arithmetic::traits::LogBaseAssign;
1129 /// use malachite_float::Float;
1130 ///
1131 /// let mut x = Float::from(1000);
1132 /// x.log_base_assign(10);
1133 /// assert_eq!(x.to_string(), "3.000");
1134 /// ```
1135 #[inline]
1136 fn log_base_assign(&mut self, base: u64) {
1137 let prec = self.significant_bits();
1138 self.log_base_prec_round_assign(base, prec, Nearest);
1139 }
1140}
1141
1142/// Computes $\log_b x$, the base-$b$ logarithm of a primitive float, where $b$ is a `u64` greater
1143/// than 1. Using this function is more accurate than computing the logarithm using the standard
1144/// library, whose `log` is not always correctly rounded.
1145///
1146/// The base-$b$ logarithm of any negative number is `NaN`.
1147///
1148/// $$
1149/// f(x,b) = \log_b x+\varepsilon.
1150/// $$
1151/// - If $\log_b x$ is infinite, zero, or `NaN`, $\varepsilon$ may be ignored or assumed to be 0.
1152/// - If $\log_b x$ is finite and nonzero, then $|\varepsilon| < 2^{\lfloor\log_2 |\log_b
1153/// x|\rfloor-p}$, where $p$ is precision of the output (typically 24 if `T` is a [`f32`] and 53
1154/// if `T` is a [`f64`], but less if the output is subnormal).
1155///
1156/// Special cases:
1157/// - $f(\text{NaN},b)=\text{NaN}$
1158/// - $f(\infty,b)=\infty$
1159/// - $f(-\infty,b)=\text{NaN}$
1160/// - $f(\pm0.0,b)=-\infty$
1161/// - $f(1.0,b)=0.0$
1162/// - $f(x,b)=\text{NaN}$ for $x<0$
1163///
1164/// Neither overflow nor underflow is possible.
1165///
1166/// # Worst-case complexity
1167/// Constant time and additional memory.
1168///
1169/// # Panics
1170/// Panics if `base` is less than 2.
1171///
1172/// # Examples
1173/// ```
1174/// use malachite_base::num::basic::traits::NegativeInfinity;
1175/// use malachite_base::num::float::NiceFloat;
1176/// use malachite_float::float::arithmetic::log_base::primitive_float_log_base;
1177///
1178/// assert!(primitive_float_log_base(f32::NAN, 10).is_nan());
1179/// assert_eq!(
1180/// NiceFloat(primitive_float_log_base(f32::INFINITY, 10)),
1181/// NiceFloat(f32::INFINITY)
1182/// );
1183/// assert_eq!(
1184/// NiceFloat(primitive_float_log_base(0.0f32, 10)),
1185/// NiceFloat(f32::NEGATIVE_INFINITY)
1186/// );
1187/// // log_10(1000) = 3
1188/// assert_eq!(
1189/// NiceFloat(primitive_float_log_base(1000.0f32, 10)),
1190/// NiceFloat(3.0)
1191/// );
1192/// // log_3(9) = 2
1193/// assert_eq!(
1194/// NiceFloat(primitive_float_log_base(9.0f32, 3)),
1195/// NiceFloat(2.0)
1196/// );
1197/// // log_10(50)
1198/// assert_eq!(
1199/// NiceFloat(primitive_float_log_base(50.0f32, 10)),
1200/// NiceFloat(1.69897)
1201/// );
1202/// assert!(primitive_float_log_base(-1.0f32, 10).is_nan());
1203/// ```
1204#[inline]
1205#[allow(clippy::type_repetition_in_bounds)]
1206pub fn primitive_float_log_base<T: PrimitiveFloat>(x: T, base: u64) -> T
1207where
1208 Float: From<T> + PartialOrd<T>,
1209 for<'a> T: ExactFrom<&'a Float> + RoundingFrom<&'a Float>,
1210{
1211 emulate_float_to_float_fn(|x, prec| Float::log_base_prec(x, base, prec), x)
1212}
1213
1214/// Computes $\log_b x$, the base-$b$ logarithm of a [`Rational`], where $b$ is a `u64` greater than
1215/// 1, returning a primitive float result.
1216///
1217/// If the logarithm is equidistant from two primitive floats, the primitive float with fewer 1s in
1218/// its binary expansion is chosen. See [`RoundingMode`] for a description of the `Nearest` rounding
1219/// mode.
1220///
1221/// The base-$b$ logarithm of any negative number is `NaN`.
1222///
1223/// $$
1224/// f(x,b) = \log_b x+\varepsilon.
1225/// $$
1226/// - If $\log_b x$ is infinite, zero, or `NaN`, $\varepsilon$ may be ignored or assumed to be 0.
1227/// - If $\log_b x$ is finite and nonzero, then $|\varepsilon| < 2^{\lfloor\log_2 |\log_b
1228/// x|\rfloor-p}$, where $p$ is precision of the output (typically 24 if `T` is a [`f32`] and 53
1229/// if `T` is a [`f64`], but less if the output is subnormal).
1230///
1231/// Special cases:
1232/// - $f(0,b)=-\infty$
1233/// - $f(x,b)=\text{NaN}$ for $x<0$
1234/// - $f(1,b)=0.0$
1235///
1236/// Neither overflow nor underflow is possible.
1237///
1238/// # Worst-case complexity
1239/// $T(m) = O(m \log m \log\log m)$
1240///
1241/// $M(m) = O(m \log m)$
1242///
1243/// where $T$ is time, $M$ is additional memory, and $m$ is `x.significant_bits()`.
1244///
1245/// # Panics
1246/// Panics if `base` is less than 2.
1247///
1248/// # Examples
1249/// ```
1250/// use malachite_base::num::basic::traits::{NegativeInfinity, Zero};
1251/// use malachite_base::num::float::NiceFloat;
1252/// use malachite_float::float::arithmetic::log_base::primitive_float_log_base_rational;
1253/// use malachite_q::Rational;
1254///
1255/// assert_eq!(
1256/// NiceFloat(primitive_float_log_base_rational::<f64>(
1257/// &Rational::ZERO,
1258/// 10
1259/// )),
1260/// NiceFloat(f64::NEGATIVE_INFINITY)
1261/// );
1262/// // log_10(1000) = 3
1263/// assert_eq!(
1264/// NiceFloat(primitive_float_log_base_rational::<f64>(
1265/// &Rational::from(1000),
1266/// 10
1267/// )),
1268/// NiceFloat(3.0)
1269/// );
1270/// // log_3(1/9) = -2
1271/// assert_eq!(
1272/// NiceFloat(primitive_float_log_base_rational::<f64>(
1273/// &Rational::from_unsigneds(1u8, 9),
1274/// 3
1275/// )),
1276/// NiceFloat(-2.0)
1277/// );
1278/// // log_10(1/3)
1279/// assert_eq!(
1280/// NiceFloat(primitive_float_log_base_rational::<f64>(
1281/// &Rational::from_unsigneds(1u8, 3),
1282/// 10
1283/// )),
1284/// NiceFloat(-0.47712125471966244)
1285/// );
1286/// assert_eq!(
1287/// NiceFloat(primitive_float_log_base_rational::<f64>(
1288/// &Rational::from(-1000),
1289/// 10
1290/// )),
1291/// NiceFloat(f64::NAN)
1292/// );
1293/// ```
1294#[inline]
1295#[allow(clippy::type_repetition_in_bounds)]
1296pub fn primitive_float_log_base_rational<T: PrimitiveFloat>(x: &Rational, base: u64) -> T
1297where
1298 Float: PartialOrd<T>,
1299 for<'a> T: ExactFrom<&'a Float> + RoundingFrom<&'a Float>,
1300{
1301 emulate_rational_to_float_fn(
1302 |x, prec| Float::log_base_rational_prec_ref(x, base, prec),
1303 x,
1304 )
1305}