Skip to main content

malachite_float/float/arithmetic/
sum.rs

1// Copyright © 2026 Mikhail Hogrefe
2//
3// Uses code adopted from the GNU MPFR Library.
4//
5//      Copyright 2014-2025 Free Software Foundation, Inc.
6//
7//      Contributed by the AriC 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::emulate_float_slice_to_float_fn;
17use crate::float::{MAX_EXPONENT_I64, MIN_EXPONENT_I64};
18use crate::{
19    Float, float_infinity, float_nan, float_negative_infinity, float_negative_zero, float_zero,
20};
21use alloc::vec::Vec;
22use core::cmp::Ordering::{self, *};
23use core::iter::Sum;
24use malachite_base::num::arithmetic::traits::ShlRoundAssign;
25use malachite_base::num::basic::floats::PrimitiveFloat;
26use malachite_base::num::conversion::traits::ExactFrom;
27use malachite_base::num::logic::traits::SignificantBits;
28use malachite_base::rounding_modes::RoundingMode::{self, *};
29use malachite_nz::natural::arithmetic::float::sum::{
30    FloatSumInput, FloatSumResult, sum_float_significands,
31};
32
33// Update the sticky sign of an all-zero result with a new zero term of sign `s` (1 or -1): if all
34// the zeros seen so far have the same sign, the result keeps that sign; otherwise the sign of the
35// zero result depends only on the rounding mode.
36pub(crate) fn update_zero_sign(sign_zero: &mut i8, s: i8, rm: RoundingMode) {
37    if *sign_zero == 0 {
38        *sign_zero = s;
39    } else if *sign_zero != s {
40        *sign_zero = if rm == Floor { -1 } else { 1 };
41    }
42}
43
44// This is mpfr_sum from sum.c, MPFR 4.2.2, split into the singular scan below and the
45// significand-level sum_aux in malachite-nz. The `Exact` rounding mode is handled by computing with
46// `Nearest` and panicking if the result is inexact.
47fn sum_prec_round_helper(xs: &[&Float], prec: u64, rm: RoundingMode) -> (Float, Ordering) {
48    assert_ne!(prec, 0);
49    let n = xs.len();
50    if n == 0 {
51        return (float_zero!(), Equal);
52    } else if n == 1 {
53        return Float::from_float_prec_round_ref(xs[0], prec, rm);
54    } else if n == 2 {
55        return xs[0].add_prec_round_ref_ref(xs[1], prec, rm);
56    }
57    // Check for special inputs, and determine the sign of an infinite result, the sign of an
58    // all-zero result, and the regular inputs.
59    let mut sign_inf = 0i8;
60    let mut sign_zero = 0i8;
61    let mut regulars: Vec<&Float> = Vec::new();
62    for x in xs {
63        match x {
64            float_nan!() => return (float_nan!(), Equal),
65            float_infinity!() => {
66                if sign_inf == 0 {
67                    sign_inf = 1;
68                } else if sign_inf < 0 {
69                    return (float_nan!(), Equal);
70                }
71            }
72            float_negative_infinity!() => {
73                if sign_inf == 0 {
74                    sign_inf = -1;
75                } else if sign_inf > 0 {
76                    return (float_nan!(), Equal);
77                }
78            }
79            float_zero!() => {
80                if regulars.is_empty() {
81                    // This choice is sticky when new zeros are considered.
82                    update_zero_sign(&mut sign_zero, 1, rm);
83                }
84            }
85            float_negative_zero!() => {
86                if regulars.is_empty() {
87                    update_zero_sign(&mut sign_zero, -1, rm);
88                }
89            }
90            _ => regulars.push(x),
91        }
92    }
93    // At this point the result cannot be NaN.
94    if sign_inf != 0 {
95        // At least one infinity, and all of them have the same sign. The sum is the infinity of
96        // this sign.
97        return if sign_inf > 0 {
98            (float_infinity!(), Equal)
99        } else {
100            (float_negative_infinity!(), Equal)
101        };
102    }
103    // At this point, all the inputs are finite numbers.
104    if regulars.is_empty() {
105        // All the numbers were zeros (and there is at least one). The sum is zero with sign
106        // sign_zero.
107        assert_ne!(sign_zero, 0);
108        return if sign_zero > 0 {
109            (float_zero!(), Equal)
110        } else {
111            (float_negative_zero!(), Equal)
112        };
113    }
114    // Optimize the case where there are only one or two regular numbers.
115    if regulars.len() == 1 {
116        return Float::from_float_prec_round_ref(regulars[0], prec, rm);
117    } else if regulars.len() == 2 {
118        return regulars[0].add_prec_round_ref_ref(regulars[1], prec, rm);
119    }
120    let (kernel_rm, exact) = if rm == Exact {
121        (Nearest, true)
122    } else {
123        (rm, false)
124    };
125    let inputs: Vec<FloatSumInput> = regulars
126        .iter()
127        .map(|x| {
128            let Float(Finite {
129                sign,
130                exponent,
131                precision,
132                significand,
133            }) = x
134            else {
135                unreachable!()
136            };
137            FloatSumInput {
138                sign: *sign,
139                exp: i64::from(*exponent),
140                prec: *precision,
141                significand,
142            }
143        })
144        .collect();
145    complete_sum_result(
146        sum_float_significands(&inputs, prec, kernel_rm),
147        prec,
148        rm,
149        exact,
150        "Inexact Float sum",
151    )
152}
153
154// Convert a summation-kernel result into a `Float`, applying the cancellation-zero sign rule and
155// the exponent range check. `exact_message` is the panic message demanded when the caller's
156// rounding mode was `Exact` (indicated by `exact`) but the result is inexact.
157pub(crate) fn complete_sum_result(
158    result: FloatSumResult,
159    prec: u64,
160    rm: RoundingMode,
161    exact: bool,
162    exact_message: &str,
163) -> (Float, Ordering) {
164    match result {
165        FloatSumResult::Zero => {
166            // The exact sum of nonzero values is zero, which is +0 except in the Floor rounding
167            // mode, as specified according to the IEEE 754 rules for the addition of two numbers.
168            if rm == Floor {
169                (float_negative_zero!(), Equal)
170            } else {
171                (float_zero!(), Equal)
172            }
173        }
174        FloatSumResult::Regular {
175            sign,
176            exp,
177            significand,
178            o,
179        } => {
180            if exact {
181                assert_eq!(o, Equal, "{exact_message}");
182            }
183            if (MIN_EXPONENT_I64..=MAX_EXPONENT_I64).contains(&exp) {
184                (
185                    Float(Finite {
186                        sign,
187                        exponent: i32::exact_from(exp),
188                        precision: prec,
189                        significand,
190                    }),
191                    o,
192                )
193            } else {
194                // The exponent is out of range; construct the value at a safe exponent and use a
195                // saturating shift to apply the overflow or underflow rules, in the same
196                // round-then-check-range order as MPFR.
197                assert!(!exact, "{exact_message}");
198                let mut f = Float(Finite {
199                    sign,
200                    exponent: 1,
201                    precision: prec,
202                    significand,
203                });
204                let o_shift = f.shl_round_assign(exp - 1, rm);
205                (f, if o_shift == Equal { o } else { o_shift })
206            }
207        }
208    }
209}
210
211// The precision used by `sum_round` and the `Sum` implementations: the maximum precision of the
212// inputs, or 1 if there are none.
213pub(crate) fn max_prec<'a, I: Iterator<Item = &'a Float>>(xs: I) -> u64 {
214    xs.map(SignificantBits::significant_bits).max().unwrap_or(1)
215}
216
217impl Float {
218    /// Computes the sum of a slice of [`Float`]s, rounding the result to the specified precision
219    /// and with the specified rounding mode. An [`Ordering`] is also returned, indicating whether
220    /// the rounded sum is less than, equal to, or greater than the exact sum. Although `NaN`s are
221    /// not comparable to any [`Float`], whenever this function returns a `NaN` it also returns
222    /// `Equal`.
223    ///
224    /// Only a single rounding is performed, no matter how many inputs there are: the result is the
225    /// correctly-rounded exact sum, with no intermediate rounding, overflow, or underflow.
226    ///
227    /// See [`RoundingMode`] for a description of the possible rounding modes.
228    ///
229    /// $$
230    /// f((x_i)_ {i=0}^{n-1}, p, m) = \sum_ {i=0}^{n-1} x_i + \varepsilon.
231    /// $$
232    /// - If $\sum_ {i=0}^{n-1} x_i$ is infinite, zero, or `NaN`, $\varepsilon$ may be ignored or
233    ///   assumed to be 0.
234    /// - If $\sum_ {i=0}^{n-1} x_i$ is finite and nonzero, and $m$ is not `Nearest`, then
235    ///   $|\varepsilon| < 2^{\lfloor\log_2 |\sum_ {i=0}^{n-1} x_i|\rfloor-p+1}$.
236    /// - If $\sum_ {i=0}^{n-1} x_i$ is finite and nonzero, and $m$ is `Nearest`, then
237    ///   $|\varepsilon| \leq 2^{\lfloor\log_2 |\sum_ {i=0}^{n-1} x_i|\rfloor-p}$.
238    ///
239    /// If the output has a precision, it is `prec`.
240    ///
241    /// Special cases:
242    /// - The sum of no [`Float`]s is $0.0$.
243    /// - If any input is `NaN`, or if the inputs include both $\infty$ and $-\infty$, the sum is
244    ///   `NaN`.
245    /// - Otherwise, if any input is $\infty$, the sum is $\infty$, and if any input is $-\infty$,
246    ///   the sum is $-\infty$.
247    /// - If every input is a zero and all of them have the same sign, the sum is a zero of that
248    ///   sign.
249    /// - If every input is a zero and they do not all have the same sign, the sum is $0.0$, unless
250    ///   $m$ is `Floor`, in which case it is $-0.0$.
251    /// - If the inputs include a nonzero value but sum to zero exactly, the sum is $0.0$, unless
252    ///   $m$ is `Floor`, in which case it is $-0.0$.
253    ///
254    /// Overflow and underflow:
255    /// - If $f((x_i)_ {i=0}^{n-1},p,m)\geq 2^{2^{30}-1}$ and $m$ is `Ceiling`, `Up`, or `Nearest`,
256    ///   $\infty$ is returned instead.
257    /// - If $f((x_i)_ {i=0}^{n-1},p,m)\geq 2^{2^{30}-1}$ and $m$ is `Floor` or `Down`,
258    ///   $(1-(1/2)^p)2^{2^{30}-1}$ is returned instead.
259    /// - If $f((x_i)_ {i=0}^{n-1},p,m)\leq -2^{2^{30}-1}$ and $m$ is `Floor`, `Up`, or `Nearest`,
260    ///   $-\infty$ is returned instead.
261    /// - If $f((x_i)_ {i=0}^{n-1},p,m)\leq -2^{2^{30}-1}$ and $m$ is `Ceiling` or `Down`,
262    ///   $-(1-(1/2)^p)2^{2^{30}-1}$ is returned instead.
263    /// - If $0<f((x_i)_ {i=0}^{n-1},p,m)<2^{-2^{30}}$, and $m$ is `Floor` or `Down`, $0.0$ is
264    ///   returned instead.
265    /// - If $0<f((x_i)_ {i=0}^{n-1},p,m)<2^{-2^{30}}$, and $m$ is `Ceiling` or `Up`, $2^{-2^{30}}$
266    ///   is returned instead.
267    /// - If $0<f((x_i)_ {i=0}^{n-1},p,m)\leq2^{-2^{30}-1}$, and $m$ is `Nearest`, $0.0$ is returned
268    ///   instead.
269    /// - If $2^{-2^{30}-1}<f((x_i)_ {i=0}^{n-1},p,m)<2^{-2^{30}}$, and $m$ is `Nearest`,
270    ///   $2^{-2^{30}}$ is returned instead.
271    /// - If $-2^{-2^{30}}<f((x_i)_ {i=0}^{n-1},p,m)<0$, and $m$ is `Ceiling` or `Down`, $-0.0$ is
272    ///   returned instead.
273    /// - If $-2^{-2^{30}}<f((x_i)_ {i=0}^{n-1},p,m)<0$, and $m$ is `Floor` or `Up`, $-2^{-2^{30}}$
274    ///   is returned instead.
275    /// - If $-2^{-2^{30}-1}\leq f((x_i)_ {i=0}^{n-1},p,m)<0$, and $m$ is `Nearest`, $-0.0$ is
276    ///   returned instead.
277    /// - If $-2^{-2^{30}}<f((x_i)_ {i=0}^{n-1},p,m)<-2^{-2^{30}-1}$, and $m$ is `Nearest`,
278    ///   $-2^{-2^{30}}$ is returned instead.
279    ///
280    /// If you know you'll be using `Nearest`, consider using [`Float::sum_prec`] instead. If you
281    /// know that your target precision is the maximum of the precisions of the inputs, consider
282    /// using [`Float::sum_round`] instead. If both of these things are true, consider summing an
283    /// iterator with [`Sum`] instead.
284    ///
285    /// # Worst-case complexity
286    /// $T(n, m, p) = O(n + m (n + p))$
287    ///
288    /// $M(n, p) = O(n + p)$
289    ///
290    /// where $T$ is time, $M$ is additional memory, $n$ is `xs.len()`, $m$ is
291    /// `u64::sum(xs.map(Float::significant_bits))`, and $p$ is `prec`: each significand bit enters
292    /// the accumulator at most once, but adversarially-placed cancelling clusters can force a pass
293    /// over all $n$ inputs and the $O(p + \log n)$-bit accumulator for every such bit.
294    ///
295    /// # Panics
296    /// Panics if `prec` is zero, or if `rm` is `Exact` and the exact sum is not exactly
297    /// representable with `prec` bits.
298    ///
299    /// # Examples
300    /// ```
301    /// use malachite_base::num::arithmetic::traits::PowerOf2;
302    /// use malachite_base::num::basic::traits::One;
303    /// use malachite_base::rounding_modes::RoundingMode::*;
304    /// use malachite_float::Float;
305    /// use std::cmp::Ordering::*;
306    ///
307    /// let xs = [Float::ONE, Float::power_of_2(-20i64), Float::power_of_2(-40i64)];
308    ///
309    /// let (sum, o) = Float::sum_prec_round(&xs, 10, Floor);
310    /// assert_eq!(sum.to_string(), "1.0000");
311    /// assert_eq!(o, Less);
312    ///
313    /// let (sum, o) = Float::sum_prec_round(&xs, 10, Ceiling);
314    /// assert_eq!(sum.to_string(), "1.0020");
315    /// assert_eq!(o, Greater);
316    ///
317    /// let (sum, o) = Float::sum_prec_round(&xs, 10, Nearest);
318    /// assert_eq!(sum.to_string(), "1.0000");
319    /// assert_eq!(o, Less);
320    ///
321    /// let (sum, o) = Float::sum_prec_round(&xs, 30, Floor);
322    /// assert_eq!(sum.to_string(), "1.0000009537");
323    /// assert_eq!(o, Less);
324    ///
325    /// let (sum, o) = Float::sum_prec_round(&xs, 30, Ceiling);
326    /// assert_eq!(sum.to_string(), "1.0000009555");
327    /// assert_eq!(o, Greater);
328    ///
329    /// let (sum, o) = Float::sum_prec_round(&xs, 30, Nearest);
330    /// assert_eq!(sum.to_string(), "1.0000009537");
331    /// assert_eq!(o, Less);
332    /// ```
333    pub fn sum_prec_round(xs: &[Self], prec: u64, rm: RoundingMode) -> (Self, Ordering) {
334        let refs: Vec<&Self> = xs.iter().collect();
335        sum_prec_round_helper(&refs, prec, rm)
336    }
337
338    /// Computes the sum of a slice of [`Float`]s, rounding the result to the nearest value of the
339    /// specified precision. An [`Ordering`] is also returned, indicating whether the rounded sum is
340    /// less than, equal to, or greater than the exact sum. Although `NaN`s are not comparable to
341    /// any [`Float`], whenever this function returns a `NaN` it also returns `Equal`.
342    ///
343    /// Only a single rounding is performed, no matter how many inputs there are: the result is the
344    /// correctly-rounded exact sum, with no intermediate rounding, overflow, or underflow.
345    ///
346    /// If the sum is equidistant from two [`Float`]s with the specified precision, the [`Float`]
347    /// with fewer 1s in its binary expansion is chosen. See [`RoundingMode`] for a description of
348    /// the `Nearest` rounding mode.
349    ///
350    /// $$
351    /// f((x_i)_ {i=0}^{n-1}, p) = \sum_ {i=0}^{n-1} x_i + \varepsilon.
352    /// $$
353    /// - If $\sum_ {i=0}^{n-1} x_i$ is infinite, zero, or `NaN`, $\varepsilon$ may be ignored or
354    ///   assumed to be 0.
355    /// - If $\sum_ {i=0}^{n-1} x_i$ is finite and nonzero, then $|\varepsilon| \leq
356    ///   2^{\lfloor\log_2 |\sum_ {i=0}^{n-1} x_i|\rfloor-p}$.
357    ///
358    /// If the output has a precision, it is `prec`.
359    ///
360    /// See [`Float::sum_prec_round`] for a description of the special cases and of overflow and
361    /// underflow behavior.
362    ///
363    /// If you know that your target precision is the maximum of the precisions of the inputs,
364    /// consider summing an iterator with [`Sum`] instead.
365    ///
366    /// # Worst-case complexity
367    /// $T(n, m, p) = O(n + m (n + p))$
368    ///
369    /// $M(n, p) = O(n + p)$
370    ///
371    /// where $T$ is time, $M$ is additional memory, $n$ is `xs.len()`, $m$ is
372    /// `u64::sum(xs.map(Float::significant_bits))`, and $p$ is `prec`.
373    ///
374    /// # Panics
375    /// Panics if `prec` is zero.
376    ///
377    /// # Examples
378    /// ```
379    /// use malachite_base::num::arithmetic::traits::PowerOf2;
380    /// use malachite_base::num::basic::traits::{One, Two};
381    /// use malachite_float::Float;
382    /// use std::cmp::Ordering::*;
383    ///
384    /// let (sum, o) = Float::sum_prec(&[Float::ONE, Float::TWO, Float::from(3)], 10);
385    /// assert_eq!(sum.to_string(), "6.0000");
386    /// assert_eq!(o, Equal);
387    ///
388    /// let (sum, o) = Float::sum_prec(
389    ///     &[Float::ONE, Float::power_of_2(-20i64), Float::power_of_2(-40i64)],
390    ///     30,
391    /// );
392    /// assert_eq!(sum.to_string(), "1.0000009537");
393    /// assert_eq!(o, Less);
394    /// ```
395    #[inline]
396    pub fn sum_prec(xs: &[Self], prec: u64) -> (Self, Ordering) {
397        Self::sum_prec_round(xs, prec, Nearest)
398    }
399
400    /// Computes the sum of a slice of [`Float`]s, rounding the result with the specified rounding
401    /// mode. The precision of the result is the maximum of the precisions of the inputs (or 1 if
402    /// there are no inputs). An [`Ordering`] is also returned, indicating whether the rounded sum
403    /// is less than, equal to, or greater than the exact sum. Although `NaN`s are not comparable to
404    /// any [`Float`], whenever this function returns a `NaN` it also returns `Equal`.
405    ///
406    /// Only a single rounding is performed, no matter how many inputs there are: the result is the
407    /// correctly-rounded exact sum, with no intermediate rounding, overflow, or underflow.
408    ///
409    /// See [`RoundingMode`] for a description of the possible rounding modes.
410    ///
411    /// $$
412    /// f((x_i)_ {i=0}^{n-1}, m) = \sum_ {i=0}^{n-1} x_i + \varepsilon.
413    /// $$
414    /// - If $\sum_ {i=0}^{n-1} x_i$ is infinite, zero, or `NaN`, $\varepsilon$ may be ignored or
415    ///   assumed to be 0.
416    /// - If $\sum_ {i=0}^{n-1} x_i$ is finite and nonzero, and $m$ is not `Nearest`, then
417    ///   $|\varepsilon| < 2^{\lfloor\log_2 |\sum_ {i=0}^{n-1} x_i|\rfloor-p+1}$, where $p$ is the
418    ///   maximum precision of the inputs.
419    /// - If $\sum_ {i=0}^{n-1} x_i$ is finite and nonzero, and $m$ is `Nearest`, then
420    ///   $|\varepsilon| \leq 2^{\lfloor\log_2 |\sum_ {i=0}^{n-1} x_i|\rfloor-p}$, where $p$ is the
421    ///   maximum precision of the inputs.
422    ///
423    /// See [`Float::sum_prec_round`] for a description of the special cases and of overflow and
424    /// underflow behavior.
425    ///
426    /// If you know you'll be using `Nearest`, consider summing an iterator with [`Sum`] instead.
427    ///
428    /// # Worst-case complexity
429    /// $T(n, m) = O(m (n + m))$
430    ///
431    /// $M(n, m) = O(n + m)$
432    ///
433    /// where $T$ is time, $M$ is additional memory, $n$ is `xs.len()`, and $m$ is
434    /// `u64::sum(xs.map(Float::significant_bits))`.
435    ///
436    /// # Panics
437    /// Panics if `rm` is `Exact` and the exact sum is not exactly representable with the maximum of
438    /// the precisions of the inputs.
439    ///
440    /// # Examples
441    /// ```
442    /// use malachite_base::num::arithmetic::traits::PowerOf2;
443    /// use malachite_base::rounding_modes::RoundingMode::*;
444    /// use malachite_float::Float;
445    /// use std::cmp::Ordering::*;
446    ///
447    /// let xs = [Float::one_prec(10), Float::power_of_2(-20i64), Float::power_of_2(-40i64)];
448    ///
449    /// let (sum, o) = Float::sum_round(&xs, Floor);
450    /// assert_eq!(sum.to_string(), "1.0000");
451    /// assert_eq!(o, Less);
452    ///
453    /// let (sum, o) = Float::sum_round(&xs, Ceiling);
454    /// assert_eq!(sum.to_string(), "1.0020");
455    /// assert_eq!(o, Greater);
456    /// ```
457    #[inline]
458    pub fn sum_round(xs: &[Self], rm: RoundingMode) -> (Self, Ordering) {
459        Self::sum_prec_round(xs, max_prec(xs.iter()), rm)
460    }
461}
462
463/// Computes the sum of a slice of primitive floats, with a single rounding.
464///
465/// The result is correctly rounded to the nearest value: the sum is computed as if in infinite
466/// precision and rounded only once, at the end, no matter how many inputs there are. This includes
467/// gradual underflow: results in the subnormal range are correctly rounded to their reduced
468/// precisions.
469///
470/// $$
471/// f((x_i)_ {i=0}^{n-1}) = \sum_ {i=0}^{n-1} x_i + \varepsilon.
472/// $$
473/// - If $\sum_ {i=0}^{n-1} x_i$ is infinite, zero, or `NaN`, $\varepsilon$ may be ignored or
474///   assumed to be 0.
475/// - If $\sum_ {i=0}^{n-1} x_i$ is finite and nonzero, then $|\varepsilon| \leq 2^{\lfloor\log_2
476///   |\sum_ {i=0}^{n-1} x_i|\rfloor-p}$, where $p$ is the precision of the output (typically 24 if
477///   `T` is a [`f32`] and 53 if `T` is a [`f64`], but less if the output is subnormal).
478///
479/// Special cases:
480/// - The sum of no floats is $0.0$.
481/// - If any input is `NaN`, or if the inputs include both $\infty$ and $-\infty$, the sum is `NaN`.
482/// - Otherwise, if any input is $\infty$, the sum is $\infty$, and if any input is $-\infty$, the
483///   sum is $-\infty$.
484/// - If every input is a zero and all of them have the same sign, the sum is a zero of that sign.
485///   If they do not all have the same sign, or if the inputs include a nonzero value but sum to
486///   zero exactly, the sum is $0.0$.
487///
488/// If the result overflows, $\pm\infty$ is returned, and if it underflows, $\pm0.0$ is returned.
489///
490/// # Worst-case complexity
491/// $T(n) = O(n)$
492///
493/// $M(n) = O(n)$
494///
495/// where $T$ is time, $M$ is additional memory, and $n$ is `xs.len()`: a primitive float's exponent
496/// range is bounded, so the summation window is repositioned only a constant number of times.
497///
498/// # Examples
499/// ```
500/// use malachite_base::num::float::NiceFloat;
501/// use malachite_float::float::arithmetic::sum::primitive_float_sum;
502///
503/// // Each addition of 0.1 in a naive fold rounds, but here only one rounding is performed.
504/// assert_eq!(
505///     NiceFloat(primitive_float_sum(&[0.1f64; 10])),
506///     NiceFloat(1.0)
507/// );
508/// assert_eq!(
509///     NiceFloat([0.1f64; 10].iter().sum::<f64>()),
510///     NiceFloat(0.9999999999999999)
511/// );
512/// ```
513#[allow(clippy::type_repetition_in_bounds)]
514#[inline]
515pub fn primitive_float_sum<T: PrimitiveFloat>(xs: &[T]) -> T
516where
517    Float: From<T> + PartialOrd<T>,
518    for<'a> T: ExactFrom<&'a Float>,
519{
520    emulate_float_slice_to_float_fn(Float::sum_prec, xs)
521}
522
523impl Sum<Self> for Float {
524    /// Adds up all the [`Float`]s in an iterator.
525    ///
526    /// The result has the maximum of the precisions of the inputs (or 1 if there are no inputs),
527    /// and the sum is rounded to nearest. Only a single rounding is performed, no matter how many
528    /// inputs there are: the result is the correctly-rounded exact sum, with no intermediate
529    /// rounding, overflow, or underflow.
530    ///
531    /// $$
532    /// f((x_i)_ {i=0}^{n-1}) = \sum_ {i=0}^{n-1} x_i + \varepsilon.
533    /// $$
534    /// - If $\sum_ {i=0}^{n-1} x_i$ is infinite, zero, or `NaN`, $\varepsilon$ may be ignored or
535    ///   assumed to be 0.
536    /// - If $\sum_ {i=0}^{n-1} x_i$ is finite and nonzero, then $|\varepsilon| \leq
537    ///   2^{\lfloor\log_2 |\sum_ {i=0}^{n-1} x_i|\rfloor-p}$, where $p$ is the maximum precision of
538    ///   the inputs.
539    ///
540    /// See [`Float::sum_prec_round`] for a description of the special cases and of overflow and
541    /// underflow behavior.
542    ///
543    /// # Worst-case complexity
544    /// $T(n, m) = O(m (n + m))$
545    ///
546    /// $M(n, m) = O(n + m)$
547    ///
548    /// where $T$ is time, $M$ is additional memory, $n$ is `xs.count()`, and $m$ is
549    /// `u64::sum(xs.map(Float::significant_bits))`.
550    ///
551    /// # Examples
552    /// ```
553    /// use core::iter::Sum;
554    /// use malachite_base::num::basic::traits::{One, Two};
555    /// use malachite_float::Float;
556    ///
557    /// let sum = Float::sum([Float::ONE, Float::TWO, Float::from(3)].into_iter());
558    /// assert_eq!(sum.to_string(), "6.0");
559    ///
560    /// // All twenty inputs have precision 1, so the result has precision 1, but the sum is still
561    /// // exact: it is computed as if in infinite precision and rounded only once.
562    /// let sum = Float::sum(vec![Float::ONE; 20].into_iter());
563    /// assert_eq!(sum.to_string(), "16.0");
564    /// ```
565    fn sum<I>(xs: I) -> Self
566    where
567        I: Iterator<Item = Self>,
568    {
569        let xs: Vec<Self> = xs.collect();
570        let refs: Vec<&Self> = xs.iter().collect();
571        sum_prec_round_helper(&refs, max_prec(xs.iter()), Nearest).0
572    }
573}
574
575impl<'a> Sum<&'a Self> for Float {
576    /// Adds up all the [`Float`]s in an iterator of [`Float`] references.
577    ///
578    /// The result has the maximum of the precisions of the inputs (or 1 if there are no inputs),
579    /// and the sum is rounded to nearest. Only a single rounding is performed, no matter how many
580    /// inputs there are: the result is the correctly-rounded exact sum, with no intermediate
581    /// rounding, overflow, or underflow.
582    ///
583    /// $$
584    /// f((x_i)_ {i=0}^{n-1}) = \sum_ {i=0}^{n-1} x_i + \varepsilon.
585    /// $$
586    /// - If $\sum_ {i=0}^{n-1} x_i$ is infinite, zero, or `NaN`, $\varepsilon$ may be ignored or
587    ///   assumed to be 0.
588    /// - If $\sum_ {i=0}^{n-1} x_i$ is finite and nonzero, then $|\varepsilon| \leq
589    ///   2^{\lfloor\log_2 |\sum_ {i=0}^{n-1} x_i|\rfloor-p}$, where $p$ is the maximum precision of
590    ///   the inputs.
591    ///
592    /// See [`Float::sum_prec_round`] for a description of the special cases and of overflow and
593    /// underflow behavior.
594    ///
595    /// # Worst-case complexity
596    /// $T(n, m) = O(m (n + m))$
597    ///
598    /// $M(n, m) = O(n + m)$
599    ///
600    /// where $T$ is time, $M$ is additional memory, $n$ is `xs.count()`, and $m$ is
601    /// `u64::sum(xs.map(Float::significant_bits))`.
602    ///
603    /// # Examples
604    /// ```
605    /// use core::iter::Sum;
606    /// use malachite_base::num::basic::traits::{One, Two};
607    /// use malachite_float::Float;
608    ///
609    /// let xs = vec![Float::ONE, Float::TWO, Float::from(3)];
610    /// assert_eq!(Float::sum(xs.iter()).to_string(), "6.0");
611    /// ```
612    fn sum<I>(xs: I) -> Self
613    where
614        I: Iterator<Item = &'a Self>,
615    {
616        let xs: Vec<&Self> = xs.collect();
617        sum_prec_round_helper(&xs, max_prec(xs.iter().copied()), Nearest).0
618    }
619}