malachite_float/float/arithmetic/log_base_float_base_1_plus_x.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::InnerFloat::{Infinity, NaN};
10use crate::float::arithmetic::log_base::{
11 dyadic_1p_log_of_root, dyadic_primitive_root, odd_significand_and_exponent,
12};
13use crate::float::basic::extended::ExtendedFloat;
14use crate::{
15 Float, emulate_float_float_to_float_fn, float_infinity, float_nan, float_negative_infinity,
16};
17use core::cmp::Ordering::{self, *};
18use malachite_base::num::arithmetic::traits::{
19 CeilingLogBase2, LogBaseOf1PlusX, LogBaseOf1PlusXAssign,
20};
21use malachite_base::num::basic::floats::PrimitiveFloat;
22use malachite_base::num::basic::integers::PrimitiveInt;
23use malachite_base::num::basic::traits::{NegativeZero, Zero as ZeroTrait};
24use malachite_base::num::conversion::traits::{ExactFrom, RoundingFrom};
25use malachite_base::num::logic::traits::SignificantBits;
26use malachite_base::rounding_modes::RoundingMode::{self, *};
27use malachite_nz::natural::arithmetic::float_extras::float_can_round;
28use malachite_nz::platform::Limb;
29use malachite_q::Rational;
30
31// Returns `Some(log_base(1 + x))` when it is rational, and `None` when it is irrational. `x` must
32// be a finite [`Float`] in (-1, 0) or positive (not 0), and `base` a finite positive [`Float`] not
33// equal to 1.
34//
35// `log_base(1 + x)` is rational exactly when `1 + x` and `base` are commensurable. `1 + x` is
36// formed exactly as a `Rational`, and `base` is dyadic, so this reuses
37// `rational_log_base_rational_rational_base`; for a base in (0, 1) -- where
38// `Rational::checked_log_base` requires a base above 1 -- the identity `log_b(y) = -log_{1/b}(y)`
39// reduces to a base above 1. Balloon-safe via the `64 * prec` size bound on `x`'s exponent and the
40// operands' precisions (an `x` near -1 has a near-zero exponent and is materialized, but `1 + x` is
41// then bounded by `x`'s precision).
42pub(crate) fn log_base_float_base_1_plus_x_rational(x: &Float, base: &Float) -> Option<Rational> {
43 if *x == 0u32 {
44 return Some(Rational::ZERO);
45 }
46 // The base's primitive root comes from its odd significand and exponent without materializing
47 // it, and `1 + x` is matched against that root implicitly (see `dyadic_1p_log_of_root`): its
48 // integer form can be enormous even when `x` has few bits, so it is only materialized to verify
49 // a match that a cheap congruence filter has already endorsed. No size cutoff: skipping the
50 // check when the result is exactly representable would leave the Ziv loop unable to terminate.
51 let (s_b, t_b) = odd_significand_and_exponent(base);
52 let (z, h, e_base) = dyadic_primitive_root(&s_b, t_b);
53 let m = dyadic_1p_log_of_root(x, z, &h)?;
54 Some(Rational::from_signeds(m, i64::exact_from(e_base)))
55}
56
57// The computation of log_base(1 + x) for a `Float` base is done by log_base(1 + x) = log_2(1 + x) /
58// log_2(base). The inputs are a finite `Float` `x` in (-1, 0) or positive (not 0), and a finite
59// positive `Float` base not equal to 1.
60//
61// Both logs are ordinary native `Float`s. `log_2(1 + x)` is routed through `log_base_2_1_plus_x` to
62// preserve accuracy for x near 0; it cannot underflow (|x| is at least the smallest positive
63// `Float`). Their quotient can overflow (base near 1) or underflow (x near 0), so the operands are
64// wrapped as `ExtendedFloat`s, divided in the extended range, and converted back with a single
65// `into_float_helper` clamp. A base in (0, 1) gives a negative `log_2(base)`, so the division
66// yields the (sign-flipped) result for free.
67fn log_base_float_base_1_plus_x_normal(
68 x: &Float,
69 base: &Float,
70 prec: u64,
71 rm: RoundingMode,
72) -> (Float, Ordering) {
73 // If log_base(1 + x) is rational -- 1 + x and base commensurable -- compute it directly.
74 if let Some(q) = log_base_float_base_1_plus_x_rational(x, base) {
75 return Float::from_rational_prec_round(q, prec, rm);
76 }
77 // The result is irrational, so it is never exactly representable.
78 assert_ne!(rm, Exact, "Inexact log_base_float_base_1_plus_x");
79 // The initial slack keeps working_prec at least 7, so the working_prec - 6 below stays
80 // positive.
81 let mut working_prec = prec + 6 + prec.ceiling_log_base_2();
82 let mut increment = Limb::WIDTH;
83 loop {
84 // log_2(1 + x) and log_2(base), correctly rounded and wrapped; both finite and nonzero (x
85 // is not 0 and base is not 1), neither underflowing.
86 let num = ExtendedFloat::from(x.log_base_2_1_plus_x_prec_ref(working_prec).0);
87 let den = ExtendedFloat::from(base.log_base_2_prec_ref(working_prec).0);
88 // log_2(1 + x) / log_2(base) in the extended range; cannot overflow or underflow here.
89 let quotient = num.div_prec_val_ref(&den, working_prec).0;
90 // Two correctly-rounded logs (<= 1/2 ulp each) and the division (<= 1/2 ulp) give under 2
91 // ulps total; working_prec - 6 correct bits comfortably suffice for the rounding test.
92 if float_can_round(
93 quotient.x.significand_ref().unwrap(),
94 working_prec - 6,
95 prec,
96 rm,
97 ) {
98 // Round the mantissa to prec, then place the extended exponent, clamping once to the
99 // Float range as the rounding mode dictates.
100 let (rounded, o) = Float::from_float_prec_round(quotient.x, prec, rm);
101 let mut result = ExtendedFloat::from(rounded);
102 result.exp = result.exp.checked_add(quotient.exp).unwrap();
103 return result.into_float_helper(prec, rm, o);
104 }
105 // Increase the precision.
106 working_prec += increment;
107 increment = working_prec >> 1;
108 }
109}
110
111// Computes log_base(1 + x) = ln(1 + x) / ln(base) for `Float` `x` and `base`, following IEEE
112// division of the natural logs for every special case (so the function is total: no input value
113// panics). `ln(1 + x)` uses the sign-preserving `ln_1p` convention, so `ln_1p(x)` has the sign of
114// `x` for `x` in (-1, infinity].
115fn log_base_float_base_1_plus_x_helper(
116 x: &Float,
117 base: &Float,
118 prec: u64,
119 rm: RoundingMode,
120) -> (Float, Ordering) {
121 // ln(1 + x) or ln(base) is NaN: x or base is NaN, base is negative, or 1 + x < 0.
122 if x.is_nan() || base.is_nan() {
123 return (float_nan!(), Equal);
124 }
125 if *base < 0u32 {
126 return (float_nan!(), Equal); // base negative finite or -infinity
127 }
128 if *x < -1i32 {
129 return (float_nan!(), Equal); // 1 + x < 0 (including x = -infinity)
130 }
131 // x is in [-1, infinity] and not NaN; base is +infinity, zero, or positive finite. ln_1p(x) has
132 // the sign of x (negative for x in [-1, 0) including -0.0).
133 let x_neg = x.is_sign_negative();
134 if base.is_infinite() {
135 // ln(base) = +infinity. ln_1p(x) / +infinity = 0 for finite ln_1p(x) (NaN when it is
136 // +-infinity, i.e. x = +infinity or x = -1), with the sign of ln_1p(x).
137 if x.is_infinite() || *x == -1i32 {
138 return (float_nan!(), Equal);
139 }
140 return if x_neg {
141 (Float::NEGATIVE_ZERO, Equal)
142 } else {
143 (Float::ZERO, Equal)
144 };
145 }
146 if *base == 0u32 {
147 // ln(base) = -infinity. Sign-flipped from the +infinity case.
148 if x.is_infinite() || *x == -1i32 {
149 return (float_nan!(), Equal);
150 }
151 return if x_neg {
152 (Float::ZERO, Equal)
153 } else {
154 (Float::NEGATIVE_ZERO, Equal)
155 };
156 }
157 if *base == 1u32 {
158 // ln(base) = +0. ln_1p(x) / +0 = +-infinity by the sign of ln_1p(x), or NaN for ln_1p(x) =
159 // +-0 (x = +-0).
160 if *x == 0u32 {
161 return (float_nan!(), Equal);
162 }
163 return if x_neg {
164 (float_negative_infinity!(), Equal)
165 } else {
166 (float_infinity!(), Equal)
167 };
168 }
169 // base is positive finite and not 1.
170 if x.is_infinite() {
171 // ln_1p(+infinity) = +infinity. +infinity / ln(base): +infinity for base > 1, -infinity for
172 // base < 1.
173 return if *base < 1u32 {
174 (float_negative_infinity!(), Equal)
175 } else {
176 (float_infinity!(), Equal)
177 };
178 }
179 if *x == -1i32 {
180 // ln_1p(-1) = ln(0) = -infinity. -infinity / ln(base): -infinity for base > 1, +infinity
181 // for base < 1.
182 return if *base < 1u32 {
183 (float_infinity!(), Equal)
184 } else {
185 (float_negative_infinity!(), Equal)
186 };
187 }
188 if *x == 0u32 {
189 // ln_1p(+-0) = +-0. +-0 / ln(base) = 0, with the sign of ln_1p(x) times the sign of
190 // ln(base) (positive for base > 1, negative for base < 1).
191 return if x_neg == (*base < 1u32) {
192 (Float::ZERO, Equal)
193 } else {
194 (Float::NEGATIVE_ZERO, Equal)
195 };
196 }
197 // x is finite in (-1, 0) or positive (not 0), and base is positive finite and not 1.
198 log_base_float_base_1_plus_x_normal(x, base, prec, rm)
199}
200
201impl Float {
202 /// Computes $\log_b(1+x)$, where $x$ and the base $b$ are both [`Float`]s, rounding the result
203 /// to the specified precision and with the specified rounding mode. The [`Float`] is taken by
204 /// value and the base by reference. An [`Ordering`] is also returned, indicating whether the
205 /// rounded value is less than, equal to, or greater than the exact value. Although `NaN`s are
206 /// not comparable to any [`Float`], whenever this function returns a `NaN` it also returns
207 /// `Equal`.
208 ///
209 /// $\log_b(1+x)$ is undefined for $x<-1$, so whenever $x<-1$, `NaN` is returned. Otherwise the
210 /// base may be any [`Float`]: the function is defined as $\ln(1+x) / \ln b$ for every pair,
211 /// applying IEEE division to the natural logs, and never panics on an input value. In
212 /// particular a base in $(0,1)$ gives a (sign-flipped) logarithm, and the non-normal and
213 /// degenerate bases follow the limits below.
214 ///
215 /// This computes $\log_2(1+x) / \log_2 b$, routing through
216 /// [`Float::log_base_2_1_plus_x_prec_ref`] to preserve accuracy for $x$ near 0, and wrapping
217 /// the quotient so it may overflow (base near 1) or underflow (x near 0) and be clamped exactly
218 /// once.
219 ///
220 /// See [`RoundingMode`] for a description of the possible rounding modes.
221 ///
222 /// $$
223 /// f(x,b,p,m) = \log_b(1+x)+\varepsilon.
224 /// $$
225 /// - If $\log_b(1+x)$ is infinite, zero, or `NaN`, $\varepsilon$ may be ignored or assumed to
226 /// be 0.
227 /// - If $\log_b(1+x)$ is finite and nonzero, and $m$ is not `Nearest`, then $|\varepsilon| <
228 /// 2^{\lfloor\log_2 |\log_b(1+x)|\rfloor-p+1}$.
229 /// - If $\log_b(1+x)$ is finite and nonzero, and $m$ is `Nearest`, then $|\varepsilon| \leq
230 /// 2^{\lfloor\log_2 |\log_b(1+x)|\rfloor-p}$.
231 ///
232 /// If the output has a precision, it is `prec`.
233 ///
234 /// Special cases (with $b$ the base):
235 /// - $f(\text{NaN},b,p,m)=\text{NaN}$, and $f(x,\text{NaN},p,m)=\text{NaN}$
236 /// - $f(x,b,p,m)=\text{NaN}$ for $x<-1$ or $b<0$
237 /// - $f(\infty,b,p,m)=\infty$ for $b>1$, and $-\infty$ for $0\leq b<1$
238 /// - $f(-1.0,b,p,m)=-\infty$ for $b>1$, and $\infty$ for $0<b<1$
239 /// - $f(\pm0.0,b,p,m)=0$ (the sign of $\pm0.0$ times the sign of $1/\ln b$)
240 /// - $f(x,\infty,p,m)=0$ for finite $x>-1$ with $x\neq0$ (and $\text{NaN}$ for
241 /// $x\in\{\infty,-1\}$)
242 /// - $f(x,\pm0.0,p,m)=0$ for finite $x>-1$ with $x\neq0$ (and $\text{NaN}$ for
243 /// $x\in\{\infty,-1\}$)
244 /// - $f(x,1.0,p,m)=\infty$ for $x>0$ or $x=\infty$, $-\infty$ for $-1\leq x<0$, and
245 /// $\text{NaN}$ for $x=\pm0.0$
246 /// - $f(g^a-1,g^e,p,m)=a/e$ for a common rational $g$, rounded to precision $p$; the result is
247 /// exact if and only if $a/e$ is representable with precision $p$ (for example
248 /// $\log_4(1+1)=1/2$)
249 ///
250 /// This function can both overflow (for a base near 1) and underflow (for an $x$ near 0).
251 ///
252 /// # Worst-case complexity
253 /// $T(n) = O(n (\log n)^2 \log\log n)$
254 ///
255 /// $M(n) = O(n (\log n)^2)$
256 ///
257 /// where $T$ is time, $M$ is additional memory, and $n$ is `prec`.
258 ///
259 /// # Panics
260 /// Panics if `prec` is zero, or if `rm` is `Exact` but the result cannot be represented exactly
261 /// with the given precision.
262 ///
263 /// # Examples
264 /// ```
265 /// use malachite_base::rounding_modes::RoundingMode::*;
266 /// use malachite_float::Float;
267 /// use std::cmp::Ordering::*;
268 ///
269 /// let (log, o) =
270 /// Float::from(8).log_base_float_base_1_plus_x_prec_round(&Float::from(3), 10, Exact);
271 /// assert_eq!(log.to_string(), "2.0000"); // log_3(1 + 8) = log_3(9) = 2
272 /// assert_eq!(o, Equal);
273 ///
274 /// let (log, o) =
275 /// Float::from(3).log_base_float_base_1_plus_x_prec_round(&Float::from(0.5), 10, Exact);
276 /// assert_eq!(log.to_string(), "-2.0000"); // log_{1/2}(1 + 3) = log_{1/2}(4) = -2
277 /// assert_eq!(o, Equal);
278 /// ```
279 #[inline]
280 pub fn log_base_float_base_1_plus_x_prec_round(
281 self,
282 base: &Self,
283 prec: u64,
284 rm: RoundingMode,
285 ) -> (Self, Ordering) {
286 assert_ne!(prec, 0);
287 log_base_float_base_1_plus_x_helper(&self, base, prec, rm)
288 }
289
290 /// Computes $\log_b(1+x)$, where $x$ and the base $b$ are both [`Float`]s, rounding the result
291 /// to the specified precision and with the specified rounding mode. Both are taken by
292 /// reference. An [`Ordering`] is also returned, indicating whether the rounded value is less
293 /// than, equal to, or greater than the exact value.
294 ///
295 /// See [`Float::log_base_float_base_1_plus_x_prec_round`] for details, special cases, and a
296 /// description of the rounding behavior.
297 ///
298 /// # Worst-case complexity
299 /// $T(n) = O(n (\log n)^2 \log\log n)$
300 ///
301 /// $M(n) = O(n (\log n)^2)$
302 ///
303 /// where $T$ is time, $M$ is additional memory, and $n$ is `prec`.
304 ///
305 /// # Panics
306 /// Panics if `prec` is zero, or if `rm` is `Exact` but the result cannot be represented exactly
307 /// with the given precision.
308 ///
309 /// # Examples
310 /// ```
311 /// use malachite_base::rounding_modes::RoundingMode::*;
312 /// use malachite_float::Float;
313 /// use std::cmp::Ordering::*;
314 ///
315 /// let x = Float::from(7);
316 /// let (log, o) = x.log_base_float_base_1_plus_x_prec_round_ref(&Float::from(2), 10, Exact);
317 /// assert_eq!(log.to_string(), "3.0000"); // log_2(1 + 7) = log_2(8) = 3
318 /// assert_eq!(o, Equal);
319 ///
320 /// let x = Float::from(1);
321 /// let (log, o) = x.log_base_float_base_1_plus_x_prec_round_ref(&Float::from(3), 20, Floor);
322 /// assert_eq!(log.to_string(), "0.63092899"); // log_3(2), rounded down
323 /// assert_eq!(o, Less);
324 /// ```
325 pub fn log_base_float_base_1_plus_x_prec_round_ref(
326 &self,
327 base: &Self,
328 prec: u64,
329 rm: RoundingMode,
330 ) -> (Self, Ordering) {
331 assert_ne!(prec, 0);
332 log_base_float_base_1_plus_x_helper(self, base, prec, rm)
333 }
334
335 /// Computes $\log_b(1+x)$, where $x$ and the base $b$ are both [`Float`]s, rounding the result
336 /// to the nearest value of the specified precision. The [`Float`] is taken by value and the
337 /// base by reference. An [`Ordering`] is also returned.
338 ///
339 /// See [`Float::log_base_float_base_1_plus_x_prec_round`] for details and special cases.
340 ///
341 /// # Worst-case complexity
342 /// $T(n) = O(n (\log n)^2 \log\log n)$
343 ///
344 /// $M(n) = O(n (\log n)^2)$
345 ///
346 /// where $T$ is time, $M$ is additional memory, and $n$ is `prec`.
347 ///
348 /// # Panics
349 /// Panics if `prec` is zero.
350 ///
351 /// # Examples
352 /// ```
353 /// use malachite_float::Float;
354 /// use std::cmp::Ordering::*;
355 ///
356 /// let (log, o) = Float::from(8).log_base_float_base_1_plus_x_prec(&Float::from(3), 10);
357 /// assert_eq!(log.to_string(), "2.0000"); // log_3(1 + 8) = log_3(9) = 2
358 /// assert_eq!(o, Equal);
359 /// ```
360 #[inline]
361 pub fn log_base_float_base_1_plus_x_prec(self, base: &Self, prec: u64) -> (Self, Ordering) {
362 self.log_base_float_base_1_plus_x_prec_round(base, prec, Nearest)
363 }
364
365 /// Computes $\log_b(1+x)$, where $x$ and the base $b$ are both [`Float`]s, rounding the result
366 /// to the nearest value of the specified precision. Both are taken by reference. An
367 /// [`Ordering`] is also returned.
368 ///
369 /// See [`Float::log_base_float_base_1_plus_x_prec_round`] for details and special cases.
370 ///
371 /// # Worst-case complexity
372 /// $T(n) = O(n (\log n)^2 \log\log n)$
373 ///
374 /// $M(n) = O(n (\log n)^2)$
375 ///
376 /// where $T$ is time, $M$ is additional memory, and $n$ is `prec`.
377 ///
378 /// # Panics
379 /// Panics if `prec` is zero.
380 ///
381 /// # Examples
382 /// ```
383 /// use malachite_float::Float;
384 /// use std::cmp::Ordering::*;
385 ///
386 /// let (log, o) = (&Float::from(8)).log_base_float_base_1_plus_x_prec_ref(&Float::from(3), 10);
387 /// assert_eq!(log.to_string(), "2.0000"); // log_3(1 + 8) = log_3(9) = 2
388 /// assert_eq!(o, Equal);
389 /// ```
390 #[inline]
391 pub fn log_base_float_base_1_plus_x_prec_ref(
392 &self,
393 base: &Self,
394 prec: u64,
395 ) -> (Self, Ordering) {
396 self.log_base_float_base_1_plus_x_prec_round_ref(base, prec, Nearest)
397 }
398
399 /// Computes $\log_b(1+x)$, where $x$ and the base $b$ are both [`Float`]s, rounding the result
400 /// to the precision of the input and with the specified rounding mode. The [`Float`] is taken
401 /// by value and the base by reference. An [`Ordering`] is also returned.
402 ///
403 /// See [`Float::log_base_float_base_1_plus_x_prec_round`] for details and special cases.
404 ///
405 /// # Worst-case complexity
406 /// $T(n) = O(n (\log n)^2 \log\log n)$
407 ///
408 /// $M(n) = O(n (\log n)^2)$
409 ///
410 /// where $T$ is time, $M$ is additional memory, and $n$ is the precision of the input.
411 ///
412 /// # Panics
413 /// Panics if `rm` is `Exact` but the result cannot be represented exactly with the input's
414 /// precision.
415 ///
416 /// # Examples
417 /// ```
418 /// use malachite_base::rounding_modes::RoundingMode::*;
419 /// use malachite_float::Float;
420 /// use std::cmp::Ordering::*;
421 ///
422 /// let (log, o) = Float::from(8).log_base_float_base_1_plus_x_round(&Float::from(3), Exact);
423 /// assert_eq!(log.to_string(), "2.0"); // log_3(1 + 8) = log_3(9) = 2
424 /// assert_eq!(o, Equal);
425 /// ```
426 #[inline]
427 pub fn log_base_float_base_1_plus_x_round(
428 self,
429 base: &Self,
430 rm: RoundingMode,
431 ) -> (Self, Ordering) {
432 let prec = self.significant_bits();
433 self.log_base_float_base_1_plus_x_prec_round(base, prec, rm)
434 }
435
436 /// Computes $\log_b(1+x)$, where $x$ and the base $b$ are both [`Float`]s, rounding the result
437 /// to the precision of the input and with the specified rounding mode. Both are taken by
438 /// reference. An [`Ordering`] is also returned.
439 ///
440 /// See [`Float::log_base_float_base_1_plus_x_prec_round`] for details and special cases.
441 ///
442 /// # Worst-case complexity
443 /// $T(n) = O(n (\log n)^2 \log\log n)$
444 ///
445 /// $M(n) = O(n (\log n)^2)$
446 ///
447 /// where $T$ is time, $M$ is additional memory, and $n$ is the precision of the input.
448 ///
449 /// # Panics
450 /// Panics if `rm` is `Exact` but the result cannot be represented exactly with the input's
451 /// precision.
452 ///
453 /// # Examples
454 /// ```
455 /// use malachite_base::rounding_modes::RoundingMode::*;
456 /// use malachite_float::Float;
457 /// use std::cmp::Ordering::*;
458 ///
459 /// let (log, o) =
460 /// (&Float::from(8)).log_base_float_base_1_plus_x_round_ref(&Float::from(3), Exact);
461 /// assert_eq!(log.to_string(), "2.0"); // log_3(1 + 8) = log_3(9) = 2
462 /// assert_eq!(o, Equal);
463 /// ```
464 #[inline]
465 pub fn log_base_float_base_1_plus_x_round_ref(
466 &self,
467 base: &Self,
468 rm: RoundingMode,
469 ) -> (Self, Ordering) {
470 self.log_base_float_base_1_plus_x_prec_round_ref(base, self.significant_bits(), rm)
471 }
472
473 /// Computes $\log_b(1+x)$, where $x$ and the base $b$ are both [`Float`]s, in place, rounding
474 /// the result to the specified precision and with the specified rounding mode. The base is
475 /// taken by reference. An [`Ordering`] is returned.
476 ///
477 /// See [`Float::log_base_float_base_1_plus_x_prec_round`] for details and special cases.
478 ///
479 /// # Worst-case complexity
480 /// $T(n) = O(n (\log n)^2 \log\log n)$
481 ///
482 /// $M(n) = O(n (\log n)^2)$
483 ///
484 /// where $T$ is time, $M$ is additional memory, and $n$ is `prec`.
485 ///
486 /// # Panics
487 /// Panics if `prec` is zero, or if `rm` is `Exact` but the result cannot be represented exactly
488 /// with the given precision.
489 ///
490 /// # Examples
491 /// ```
492 /// use malachite_base::rounding_modes::RoundingMode::*;
493 /// use malachite_float::Float;
494 /// use std::cmp::Ordering::*;
495 ///
496 /// let mut x = Float::from(8);
497 /// assert_eq!(
498 /// x.log_base_float_base_1_plus_x_prec_round_assign(&Float::from(3), 10, Exact),
499 /// Equal
500 /// );
501 /// assert_eq!(x.to_string(), "2.0000"); // log_3(1 + 8) = log_3(9) = 2
502 /// ```
503 #[inline]
504 pub fn log_base_float_base_1_plus_x_prec_round_assign(
505 &mut self,
506 base: &Self,
507 prec: u64,
508 rm: RoundingMode,
509 ) -> Ordering {
510 let (result, o) =
511 core::mem::take(self).log_base_float_base_1_plus_x_prec_round(base, prec, rm);
512 *self = result;
513 o
514 }
515
516 /// Computes $\log_b(1+x)$, where $x$ and the base $b$ are both [`Float`]s, in place, rounding
517 /// the result to the nearest value of the specified precision. The base is taken by reference.
518 /// An [`Ordering`] is returned.
519 ///
520 /// See [`Float::log_base_float_base_1_plus_x_prec_round`] for details and special cases.
521 ///
522 /// # Worst-case complexity
523 /// $T(n) = O(n (\log n)^2 \log\log n)$
524 ///
525 /// $M(n) = O(n (\log n)^2)$
526 ///
527 /// where $T$ is time, $M$ is additional memory, and $n$ is `prec`.
528 ///
529 /// # Panics
530 /// Panics if `prec` is zero.
531 ///
532 /// # Examples
533 /// ```
534 /// use malachite_float::Float;
535 ///
536 /// let mut x = Float::from(8);
537 /// x.log_base_float_base_1_plus_x_prec_assign(&Float::from(3), 10);
538 /// assert_eq!(x.to_string(), "2.0000"); // log_3(1 + 8) = log_3(9) = 2
539 /// ```
540 #[inline]
541 pub fn log_base_float_base_1_plus_x_prec_assign(&mut self, base: &Self, prec: u64) -> Ordering {
542 self.log_base_float_base_1_plus_x_prec_round_assign(base, prec, Nearest)
543 }
544
545 /// Computes $\log_b(1+x)$, where $x$ and the base $b$ are both [`Float`]s, in place, rounding
546 /// the result to the precision of the input and with the specified rounding mode. The base is
547 /// taken by reference. An [`Ordering`] is returned.
548 ///
549 /// See [`Float::log_base_float_base_1_plus_x_prec_round`] for details and special cases.
550 ///
551 /// # Worst-case complexity
552 /// $T(n) = O(n (\log n)^2 \log\log n)$
553 ///
554 /// $M(n) = O(n (\log n)^2)$
555 ///
556 /// where $T$ is time, $M$ is additional memory, and $n$ is the precision of the input.
557 ///
558 /// # Panics
559 /// Panics if `rm` is `Exact` but the result cannot be represented exactly with the input's
560 /// precision.
561 ///
562 /// # Examples
563 /// ```
564 /// use malachite_base::rounding_modes::RoundingMode::*;
565 /// use malachite_float::Float;
566 ///
567 /// let mut x = Float::from(8);
568 /// x.log_base_float_base_1_plus_x_round_assign(&Float::from(3), Exact);
569 /// assert_eq!(x.to_string(), "2.0"); // log_3(1 + 8) = log_3(9) = 2
570 /// ```
571 #[inline]
572 pub fn log_base_float_base_1_plus_x_round_assign(
573 &mut self,
574 base: &Self,
575 rm: RoundingMode,
576 ) -> Ordering {
577 let prec = self.significant_bits();
578 self.log_base_float_base_1_plus_x_prec_round_assign(base, prec, rm)
579 }
580}
581
582impl LogBaseOf1PlusX<Self> for Float {
583 type Output = Self;
584
585 /// Computes $\log_b(1+x)$, where $x$ and the base $b$ are both [`Float`]s, rounding the result
586 /// to the nearest value of the input's precision. Both are taken by value.
587 ///
588 /// See [`Float::log_base_float_base_1_plus_x_prec_round`] for special cases.
589 ///
590 /// # Worst-case complexity
591 /// $T(n) = O(n (\log n)^2 \log\log n)$
592 ///
593 /// $M(n) = O(n (\log n)^2)$
594 ///
595 /// where $T$ is time, $M$ is additional memory, and $n$ is the precision of the input.
596 ///
597 /// # Examples
598 /// ```
599 /// use malachite_base::num::arithmetic::traits::LogBaseOf1PlusX;
600 /// use malachite_float::Float;
601 ///
602 /// // log_3(1 + 8) = log_3(9) = 2
603 /// assert_eq!(
604 /// Float::from(8).log_base_1_plus_x(Float::from(3)).to_string(),
605 /// "2.0"
606 /// );
607 /// ```
608 #[inline]
609 fn log_base_1_plus_x(self, base: Self) -> Self {
610 let prec = self.significant_bits();
611 self.log_base_float_base_1_plus_x_prec_round(&base, prec, Nearest)
612 .0
613 }
614}
615
616impl LogBaseOf1PlusX<&Float> for &Float {
617 type Output = Float;
618
619 /// Computes $\log_b(1+x)$, where $x$ and the base $b$ are both [`Float`]s, rounding the result
620 /// to the nearest value of the input's precision. Both are taken by reference.
621 ///
622 /// See [`Float::log_base_float_base_1_plus_x_prec_round`] for special cases.
623 ///
624 /// # Worst-case complexity
625 /// $T(n) = O(n (\log n)^2 \log\log n)$
626 ///
627 /// $M(n) = O(n (\log n)^2)$
628 ///
629 /// where $T$ is time, $M$ is additional memory, and $n$ is the precision of the input.
630 ///
631 /// # Examples
632 /// ```
633 /// use malachite_base::num::arithmetic::traits::LogBaseOf1PlusX;
634 /// use malachite_float::Float;
635 ///
636 /// // log_3(1 + 8) = log_3(9) = 2
637 /// assert_eq!(
638 /// (&Float::from(8))
639 /// .log_base_1_plus_x(&Float::from(3))
640 /// .to_string(),
641 /// "2.0"
642 /// );
643 /// ```
644 #[inline]
645 fn log_base_1_plus_x(self, base: &Float) -> Float {
646 self.log_base_float_base_1_plus_x_prec_round_ref(base, self.significant_bits(), Nearest)
647 .0
648 }
649}
650
651impl LogBaseOf1PlusXAssign<&Self> for Float {
652 /// Replaces a [`Float`] $x$ with $\log_b(1+x)$, where the base $b$ is a [`Float`], rounding the
653 /// result to the nearest value of the input's precision. The base is taken by reference.
654 ///
655 /// See [`Float::log_base_float_base_1_plus_x_prec_round`] for special cases.
656 ///
657 /// # Worst-case complexity
658 /// $T(n) = O(n (\log n)^2 \log\log n)$
659 ///
660 /// $M(n) = O(n (\log n)^2)$
661 ///
662 /// where $T$ is time, $M$ is additional memory, and $n$ is the precision of the input.
663 ///
664 /// # Examples
665 /// ```
666 /// use malachite_base::num::arithmetic::traits::LogBaseOf1PlusXAssign;
667 /// use malachite_float::Float;
668 ///
669 /// let mut x = Float::from(8);
670 /// x.log_base_1_plus_x_assign(&Float::from(3));
671 /// assert_eq!(x.to_string(), "2.0"); // log_3(1 + 8) = log_3(9) = 2
672 /// ```
673 #[inline]
674 fn log_base_1_plus_x_assign(&mut self, base: &Self) {
675 let prec = self.significant_bits();
676 self.log_base_float_base_1_plus_x_prec_round_assign(base, prec, Nearest);
677 }
678}
679
680/// Computes $\log_b(1+x)$, the base-$b$ logarithm of one plus a primitive float, where the base $b$
681/// is also a primitive float, returning a primitive float result. Using this function is more
682/// accurate than computing the logarithm using the standard library, both because $1+x$ may not be
683/// representable as a primitive float and because the standard library's `log` is not always
684/// correctly rounded.
685///
686/// $\log_b(1+x)$ is undefined for $x<-1$, so whenever $x<-1$, `NaN` is returned. Otherwise the base
687/// may be any primitive float: the function is defined as $\ln(1+x) / \ln b$ and never panics on an
688/// input value. A base in $(0,1)$ gives a (sign-flipped) logarithm, and the non-normal and
689/// degenerate bases follow the limits below.
690///
691/// $$
692/// f(x,b) = \log_b(1+x)+\varepsilon.
693/// $$
694/// - If $\log_b(1+x)$ is infinite, zero, or `NaN`, $\varepsilon$ may be ignored or assumed to be 0.
695/// - If $\log_b(1+x)$ is finite and nonzero, then $|\varepsilon| < 2^{\lfloor\log_2
696/// |\log_b(1+x)|\rfloor-p}$, where $p$ is precision of the output (typically 24 if `T` is a
697/// [`f32`] and 53 if `T` is a [`f64`], but less if the output is subnormal).
698///
699/// Special cases (with $b$ the base):
700/// - $f(\text{NaN},b)=\text{NaN}$, and $f(x,\text{NaN})=\text{NaN}$
701/// - $f(x,b)=\text{NaN}$ for $x<-1$ or $b<0$
702/// - $f(\infty,b)=\infty$ for $b>1$, and $-\infty$ for $0\leq b<1$
703/// - $f(-1.0,b)=-\infty$ for $b>1$, and $\infty$ for $0<b<1$
704/// - $f(\pm0.0,b)=0$ (the sign of $\pm0.0$ times the sign of $1/\ln b$)
705/// - $f(x,\infty)=0$ for finite $x>-1$ with $x\neq0$ (and $\text{NaN}$ for $x\in\{\infty,-1\}$)
706/// - $f(x,\pm0.0)=0$ for finite $x>-1$ with $x\neq0$ (and $\text{NaN}$ for $x\in\{\infty,-1\}$)
707/// - $f(x,1.0)=\infty$ for $x>0$ or $x=\infty$, $-\infty$ for $-1\leq x<0$, and $\text{NaN}$ for
708/// $x=\pm0.0$
709///
710/// This function can both overflow (for a base near 1) and underflow (for an $x$ near 0).
711///
712/// # Worst-case complexity
713/// Constant time and additional memory.
714///
715/// # Examples
716/// ```
717/// use malachite_base::num::basic::traits::NegativeInfinity;
718/// use malachite_base::num::float::NiceFloat;
719/// use malachite_float::float::arithmetic::log_base_float_base_1_plus_x::*;
720///
721/// // log_4(1 + 3) = log_4(4) = 1
722/// assert_eq!(
723/// NiceFloat(primitive_float_log_base_float_base_1_plus_x(3.0f32, 4.0)),
724/// NiceFloat(1.0)
725/// );
726/// // log_4(1 + 1) = log_4(2) = 1/2
727/// assert_eq!(
728/// NiceFloat(primitive_float_log_base_float_base_1_plus_x(1.0f32, 4.0)),
729/// NiceFloat(0.5)
730/// );
731/// // log_(1/2)(1 + 3) = log_(1/2)(4) = -2
732/// assert_eq!(
733/// NiceFloat(primitive_float_log_base_float_base_1_plus_x(3.0f32, 0.5)),
734/// NiceFloat(-2.0)
735/// );
736/// assert_eq!(
737/// NiceFloat(primitive_float_log_base_float_base_1_plus_x(-1.0f32, 10.0)),
738/// NiceFloat(f32::NEGATIVE_INFINITY)
739/// );
740/// assert!(primitive_float_log_base_float_base_1_plus_x(-2.0f32, 10.0).is_nan());
741/// assert!(primitive_float_log_base_float_base_1_plus_x(3.0f32, f32::NAN).is_nan());
742/// ```
743#[inline]
744#[allow(clippy::type_repetition_in_bounds)]
745pub fn primitive_float_log_base_float_base_1_plus_x<T: PrimitiveFloat>(x: T, base: T) -> T
746where
747 Float: From<T> + PartialOrd<T>,
748 for<'a> T: ExactFrom<&'a Float> + RoundingFrom<&'a Float>,
749{
750 emulate_float_float_to_float_fn(
751 |x, base, prec| x.log_base_float_base_1_plus_x_prec(&base, prec),
752 x,
753 base,
754 )
755}