Skip to main content

malachite_float/float/arithmetic/
dot.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::{Finite, Infinity, NaN, Zero};
10use crate::emulate_float_slice_float_slice_to_float_fn;
11use crate::float::arithmetic::sum::{complete_sum_result, max_prec, update_zero_sign};
12use crate::{
13    Float, float_infinity, float_nan, float_negative_infinity, float_negative_zero, float_zero,
14};
15use alloc::vec::Vec;
16use core::cmp::Ordering::{self, *};
17use malachite_base::num::basic::floats::PrimitiveFloat;
18use malachite_base::num::conversion::traits::ExactFrom;
19use malachite_base::num::logic::traits::SignificantBits;
20use malachite_base::rounding_modes::RoundingMode::{self, *};
21use malachite_nz::natural::Natural;
22use malachite_nz::natural::arithmetic::float::sum::{FloatSumInput, sum_float_significands};
23
24/// Computes the dot product of two equal-length slices of primitive floats, with a single rounding.
25///
26/// The result is correctly rounded to the nearest value: the products are exact, the sum is
27/// computed as if in infinite precision, and only a single rounding is performed, at the end. This
28/// includes gradual underflow: results in the subnormal range are correctly rounded to their
29/// reduced precisions. Intermediate overflow and underflow cannot occur.
30///
31/// $$
32/// f((x_i)_ {i=0}^{n-1}, (y_i)_ {i=0}^{n-1}) = \sum_ {i=0}^{n-1} x_i y_i + \varepsilon.
33/// $$
34/// - If $\sum_ {i=0}^{n-1} x_i y_i$ is infinite, zero, or `NaN`, $\varepsilon$ may be ignored or
35///   assumed to be 0.
36/// - If $\sum_ {i=0}^{n-1} x_i y_i$ is finite and nonzero, then $|\varepsilon| \leq
37///   2^{\lfloor\log_2 |\sum_ {i=0}^{n-1} x_i y_i|\rfloor-p}$, where $p$ is the precision of the
38///   output (typically 24 if `T` is a [`f32`] and 53 if `T` is a [`f64`], but less if the output is
39///   subnormal).
40///
41/// See [`Float::dot_prec_round`] for a description of the special cases, which follow the rules of
42/// multiplication for each term and the rules of addition for their combination.
43///
44/// If the result overflows, $\pm\infty$ is returned, and if it underflows, $\pm0.0$ is returned.
45///
46/// # Worst-case complexity
47/// $T(n) = O(n)$
48///
49/// $M(n) = O(n)$
50///
51/// where $T$ is time, $M$ is additional memory, and $n$ is `xs.len()`: the products are
52/// constant-size, and a primitive float's exponent range is bounded, so the summation window is
53/// repositioned only a constant number of times.
54///
55/// # Panics
56/// Panics if `xs` and `ys` have different lengths.
57///
58/// # Examples
59/// ```
60/// use malachite_base::num::float::NiceFloat;
61/// use malachite_float::float::arithmetic::dot::primitive_float_dot;
62///
63/// // A naive fold overflows on the first product; the correctly-rounded dot product does not.
64/// let xs = [1.0e300f64, 1.0e300];
65/// let ys = [1.0e300f64, -1.0e300];
66/// assert_eq!(NiceFloat(primitive_float_dot(&xs, &ys)), NiceFloat(0.0));
67/// ```
68#[allow(clippy::type_repetition_in_bounds)]
69#[inline]
70pub fn primitive_float_dot<T: PrimitiveFloat>(xs: &[T], ys: &[T]) -> T
71where
72    Float: From<T> + PartialOrd<T>,
73    for<'a> T: ExactFrom<&'a Float>,
74{
75    emulate_float_slice_float_slice_to_float_fn(Float::dot_prec, xs, ys)
76}
77
78// A Float's finite fields: the sign, the exponent, the precision, and a reference to the
79// significand.
80fn parts(f: &Float) -> (bool, i32, u64, &Natural) {
81    let Float(Finite {
82        sign,
83        exponent,
84        precision,
85        significand,
86    }) = f
87    else {
88        unreachable!()
89    };
90    (*sign, *exponent, *precision, significand)
91}
92
93impl Float {
94    /// Computes the dot product of two equal-length slices of [`Float`]s, rounding the result to
95    /// the specified precision and with the specified rounding mode. An [`Ordering`] is also
96    /// returned, indicating whether the rounded dot product is less than, equal to, or greater than
97    /// the exact dot product. Although `NaN`s are not comparable to any [`Float`], whenever this
98    /// function returns a `NaN` it also returns `Equal`.
99    ///
100    /// The products are never rounded, and only a single rounding is performed, at the end: the
101    /// result is the correctly-rounded exact dot product. Intermediate overflow and underflow
102    /// cannot occur: each product is computed exactly at the significand level, with its exponent
103    /// tracked over a range twice as wide as a [`Float`]'s, and only the final sum is subject to
104    /// the exponent range check. (MPFR's `mpfr_dot`, which is documented as experimental, computes
105    /// each product at full precision and requires those multiplications to be exact, so it does
106    /// not handle inputs whose products leave the exponent range.)
107    ///
108    /// See [`RoundingMode`] for a description of the possible rounding modes.
109    ///
110    /// $$
111    /// f((x_i)_ {i=0}^{n-1}, (y_i)_ {i=0}^{n-1}, p, m) = \sum_ {i=0}^{n-1} x_i y_i +
112    /// \varepsilon.
113    /// $$
114    /// - If $\sum_ {i=0}^{n-1} x_i y_i$ is infinite, zero, or `NaN`, $\varepsilon$ may be ignored
115    ///   or assumed to be 0.
116    /// - If $\sum_ {i=0}^{n-1} x_i y_i$ is finite and nonzero, and $m$ is not `Nearest`, then
117    ///   $|\varepsilon| < 2^{\lfloor\log_2 |\sum_ {i=0}^{n-1} x_i y_i|\rfloor-p+1}$.
118    /// - If $\sum_ {i=0}^{n-1} x_i y_i$ is finite and nonzero, and $m$ is `Nearest`, then
119    ///   $|\varepsilon| \leq 2^{\lfloor\log_2 |\sum_ {i=0}^{n-1} x_i y_i|\rfloor-p}$.
120    ///
121    /// If the output has a precision, it is `prec`.
122    ///
123    /// Each term $x_iy_i$ follows the rules of [`Float`] multiplication, and the terms are then
124    /// combined following the rules of [`Float`] addition:
125    /// - The dot product of empty slices is $0.0$.
126    /// - If any term is a `NaN` — because an input is `NaN`, or because a zero is paired with an
127    ///   infinity — the dot product is `NaN`.
128    /// - If two infinite terms have different signs, the dot product is `NaN`. Otherwise, if any
129    ///   term is infinite, the dot product is an infinity of that sign.
130    /// - If every term is a zero and all the terms have the same sign, the dot product is a zero of
131    ///   that sign. If they do not all have the same sign, the dot product is $0.0$, unless $m$ is
132    ///   `Floor`, in which case it is $-0.0$.
133    /// - If some terms are nonzero but the exact dot product is zero, the dot product is $0.0$,
134    ///   unless $m$ is `Floor`, in which case it is $-0.0$.
135    ///
136    /// Overflow and underflow:
137    /// - If $f((x_i)_ {i=0}^{n-1},(y_i)_ {i=0}^{n-1},p,m)\geq 2^{2^{30}-1}$ and $m$ is `Ceiling`,
138    ///   `Up`, or `Nearest`, $\infty$ is returned instead.
139    /// - If $f((x_i)_ {i=0}^{n-1},(y_i)_ {i=0}^{n-1},p,m)\geq 2^{2^{30}-1}$ and $m$ is `Floor` or
140    ///   `Down`, $(1-(1/2)^p)2^{2^{30}-1}$ is returned instead.
141    /// - If $f((x_i)_ {i=0}^{n-1},(y_i)_ {i=0}^{n-1},p,m)\leq -2^{2^{30}-1}$ and $m$ is `Floor`,
142    ///   `Up`, or `Nearest`, $-\infty$ is returned instead.
143    /// - If $f((x_i)_ {i=0}^{n-1},(y_i)_ {i=0}^{n-1},p,m)\leq -2^{2^{30}-1}$ and $m$ is `Ceiling`
144    ///   or `Down`, $-(1-(1/2)^p)2^{2^{30}-1}$ is returned instead.
145    /// - If $0<f((x_i)_ {i=0}^{n-1},(y_i)_ {i=0}^{n-1},p,m)<2^{-2^{30}}$, and $m$ is `Floor` or
146    ///   `Down`, $0.0$ is returned instead.
147    /// - If $0<f((x_i)_ {i=0}^{n-1},(y_i)_ {i=0}^{n-1},p,m)<2^{-2^{30}}$, and $m$ is `Ceiling` or
148    ///   `Up`, $2^{-2^{30}}$ is returned instead.
149    /// - If $0<f((x_i)_ {i=0}^{n-1},(y_i)_ {i=0}^{n-1},p,m)\leq2^{-2^{30}-1}$, and $m$ is
150    ///   `Nearest`, $0.0$ is returned instead.
151    /// - If $2^{-2^{30}-1}<f((x_i)_ {i=0}^{n-1},(y_i)_ {i=0}^{n-1},p,m)<2^{-2^{30}}$, and $m$ is
152    ///   `Nearest`, $2^{-2^{30}}$ is returned instead.
153    /// - If $-2^{-2^{30}}<f((x_i)_ {i=0}^{n-1},(y_i)_ {i=0}^{n-1},p,m)<0$, and $m$ is `Ceiling` or
154    ///   `Down`, $-0.0$ is returned instead.
155    /// - If $-2^{-2^{30}}<f((x_i)_ {i=0}^{n-1},(y_i)_ {i=0}^{n-1},p,m)<0$, and $m$ is `Floor` or
156    ///   `Up`, $-2^{-2^{30}}$ is returned instead.
157    /// - If $-2^{-2^{30}-1}\leq f((x_i)_ {i=0}^{n-1},(y_i)_ {i=0}^{n-1},p,m)<0$, and $m$ is
158    ///   `Nearest`, $-0.0$ is returned instead.
159    /// - If $-2^{-2^{30}}<f((x_i)_ {i=0}^{n-1},(y_i)_ {i=0}^{n-1},p,m)<-2^{-2^{30}-1}$, and $m$ is
160    ///   `Nearest`, $-2^{-2^{30}}$ is returned instead.
161    ///
162    /// If you know you'll be using `Nearest`, consider using [`Float::dot_prec`] instead. If you
163    /// know that your target precision is the maximum of the precisions of the inputs, consider
164    /// using [`Float::dot_round`] instead. If both of these things are true, consider using
165    /// [`Float::dot`] instead.
166    ///
167    /// # Worst-case complexity
168    /// $T(n, m, p) = O(n + m (n + p) + m \log m \log\log m)$
169    ///
170    /// $M(n, m, p) = O(n + p + m \log m)$
171    ///
172    /// where $T$ is time, $M$ is additional memory, $n$ is `xs.len()`, $m$ is the sum of the
173    /// significant bits of the elements of `xs` and `ys`, and $p$ is `prec`: each term is an exact
174    /// significand product (mul-class in the pair's bits), and the terms then feed the summation
175    /// kernel, which inherits the summation bound.
176    ///
177    /// # Panics
178    /// Panics if `prec` is zero, if `xs` and `ys` have different lengths, or if `rm` is `Exact` and
179    /// the exact dot product is not exactly representable with `prec` bits.
180    ///
181    /// # Examples
182    /// ```
183    /// use malachite_base::num::arithmetic::traits::PowerOf2;
184    /// use malachite_base::num::basic::traits::{One, Two};
185    /// use malachite_base::rounding_modes::RoundingMode::*;
186    /// use malachite_float::Float;
187    /// use std::cmp::Ordering::*;
188    ///
189    /// let xs = [Float::ONE, Float::TWO, Float::from(3)];
190    /// let ys = [Float::from(4), Float::from(5), Float::from(6)];
191    /// let (dot, o) = Float::dot_prec_round(&xs, &ys, 10, Floor);
192    /// assert_eq!(dot.to_string(), "32.000");
193    /// assert_eq!(o, Equal);
194    ///
195    /// // 0.25 * 0.25 + 2 * 5 = 10.0625
196    /// let xs = [Float::power_of_2(-2i64), Float::TWO];
197    /// let ys = [Float::power_of_2(-2i64), Float::from(5)];
198    /// let (dot, o) = Float::dot_prec_round(&xs, &ys, 3, Floor);
199    /// assert_eq!(dot.to_string(), "10.0");
200    /// assert_eq!(o, Less);
201    ///
202    /// let (dot, o) = Float::dot_prec_round(&xs, &ys, 3, Ceiling);
203    /// assert_eq!(dot.to_string(), "12.0");
204    /// assert_eq!(o, Greater);
205    /// ```
206    pub fn dot_prec_round(
207        xs: &[Self],
208        ys: &[Self],
209        prec: u64,
210        rm: RoundingMode,
211    ) -> (Self, Ordering) {
212        // The dot product is correctly rounded: the products are never rounded, and only a single
213        // rounding is performed, at the end. Unlike MPFR's experimental `mpfr_dot`, intermediate
214        // overflow and underflow cannot occur: each product is computed exactly at the significand
215        // level, with its exponent tracked in an `i64` (twice the `Float` exponent range fits
216        // comfortably), and only the final sum is subject to the exponent range check.
217        //
218        // The `Exact` rounding mode is handled by computing with `Nearest` and panicking if the
219        // result is inexact.
220        assert_ne!(prec, 0);
221        assert_eq!(
222            xs.len(),
223            ys.len(),
224            "dot product requires slices of equal length"
225        );
226        let n = xs.len();
227        if n == 0 {
228            return (float_zero!(), Equal);
229        } else if n == 1 {
230            return xs[0].mul_prec_round_ref_ref(&ys[0], prec, rm);
231        }
232        // Classify each term x * y according to the multiplication rules, then combine the terms
233        // according to the addition rules, determining the sign of an infinite result, the sign of
234        // an all-zero result, and the regular terms.
235        let mut sign_inf = 0i8;
236        let mut sign_zero = 0i8;
237        let mut regulars: Vec<(&Self, &Self)> = Vec::new();
238        for (x, y) in xs.iter().zip(ys.iter()) {
239            if x.is_nan() || y.is_nan() {
240                return (float_nan!(), Equal);
241            }
242            let term_sign = if x.is_sign_negative() == y.is_sign_negative() {
243                1
244            } else {
245                -1
246            };
247            if x.is_infinite() || y.is_infinite() {
248                if *x == 0u32 || *y == 0u32 {
249                    // A zero times an infinity is NaN.
250                    return (float_nan!(), Equal);
251                }
252                if sign_inf == 0 {
253                    sign_inf = term_sign;
254                } else if sign_inf != term_sign {
255                    // Infinite terms of opposite signs add to NaN.
256                    return (float_nan!(), Equal);
257                }
258            } else if *x == 0u32 || *y == 0u32 {
259                if regulars.is_empty() {
260                    // This choice is sticky when new zeros are considered.
261                    update_zero_sign(&mut sign_zero, term_sign, rm);
262                }
263            } else {
264                regulars.push((x, y));
265            }
266        }
267        // At this point the result cannot be NaN.
268        if sign_inf != 0 {
269            return if sign_inf > 0 {
270                (float_infinity!(), Equal)
271            } else {
272                (float_negative_infinity!(), Equal)
273            };
274        }
275        // At this point every term is finite.
276        if regulars.is_empty() {
277            // All the terms were zeros (and there is at least one). The dot product is zero with
278            // sign sign_zero.
279            assert_ne!(sign_zero, 0);
280            return if sign_zero > 0 {
281                (float_zero!(), Equal)
282            } else {
283                (float_negative_zero!(), Equal)
284            };
285        }
286        // Optimize the case where there are only one or two regular terms, delegating to the
287        // correctly-rounded multiplication and fused multiply-add-multiply.
288        if regulars.len() == 1 {
289            return regulars[0]
290                .0
291                .mul_prec_round_ref_ref(regulars[0].1, prec, rm);
292        } else if regulars.len() == 2 {
293            return regulars[0].0.mul_add_mul_prec_round_ref_ref_ref_ref(
294                regulars[0].1,
295                regulars[1].0,
296                regulars[1].1,
297                prec,
298                rm,
299            );
300        }
301        let (kernel_rm, exact) = if rm == Exact {
302            (Nearest, true)
303        } else {
304            (rm, false)
305        };
306        // Compute each product exactly at the significand level. A Float's significand is stored
307        // limb-aligned with its top bit set, and its value is significand * 2^(exponent - 64 *
308        // len); the product of two such significands has its top bit either exactly at the combined
309        // width (in which case the product is already aligned) or one position below it (in which
310        // case a shift by 1 restores the alignment and the exponent decreases by 1).
311        let terms: Vec<(bool, i64, u64, Natural)> = regulars
312            .iter()
313            .map(|&(x, y)| {
314                let (sx, ex, px, sig_x) = parts(x);
315                let (sy, ey, py, sig_y) = parts(y);
316                let mut s = sig_x * sig_y;
317                let mut exp = i64::from(ex) + i64::from(ey);
318                let full = sig_x.significant_bits() + sig_y.significant_bits();
319                let mut term_prec = px + py;
320                if s.significant_bits() < full {
321                    s <<= 1;
322                    exp -= 1;
323                    term_prec -= 1;
324                }
325                (sx == sy, exp, term_prec, s)
326            })
327            .collect();
328        let inputs: Vec<FloatSumInput> = terms
329            .iter()
330            .map(|(sign, exp, term_prec, s)| FloatSumInput {
331                sign: *sign,
332                exp: *exp,
333                prec: *term_prec,
334                significand: s,
335            })
336            .collect();
337        complete_sum_result(
338            sum_float_significands(&inputs, prec, kernel_rm),
339            prec,
340            rm,
341            exact,
342            "Inexact Float dot product",
343        )
344    }
345
346    /// Computes the dot product of two equal-length slices of [`Float`]s, rounding the result to
347    /// the nearest value of the specified precision. An [`Ordering`] is also returned, indicating
348    /// whether the rounded dot product is less than, equal to, or greater than the exact dot
349    /// product. Although `NaN`s are not comparable to any [`Float`], whenever this function returns
350    /// a `NaN` it also returns `Equal`.
351    ///
352    /// The products are never rounded, and only a single rounding is performed, at the end: the
353    /// result is the correctly-rounded exact dot product, and intermediate overflow and underflow
354    /// cannot occur.
355    ///
356    /// If the dot product is equidistant from two [`Float`]s with the specified precision, the
357    /// [`Float`] with fewer 1s in its binary expansion is chosen. See [`RoundingMode`] for a
358    /// description of the `Nearest` rounding mode.
359    ///
360    /// $$
361    /// f((x_i)_ {i=0}^{n-1}, (y_i)_ {i=0}^{n-1}, p) = \sum_ {i=0}^{n-1} x_i y_i + \varepsilon.
362    /// $$
363    /// - If $\sum_ {i=0}^{n-1} x_i y_i$ is infinite, zero, or `NaN`, $\varepsilon$ may be ignored
364    ///   or assumed to be 0.
365    /// - If $\sum_ {i=0}^{n-1} x_i y_i$ is finite and nonzero, then $|\varepsilon| \leq
366    ///   2^{\lfloor\log_2 |\sum_ {i=0}^{n-1} x_i y_i|\rfloor-p}$.
367    ///
368    /// If the output has a precision, it is `prec`.
369    ///
370    /// See [`Float::dot_prec_round`] for a description of the special cases and of overflow and
371    /// underflow behavior.
372    ///
373    /// If you know that your target precision is the maximum of the precisions of the inputs,
374    /// consider using [`Float::dot`] instead.
375    ///
376    /// # Worst-case complexity
377    /// $T(n, m, p) = O(n + m (n + p) + m \log m \log\log m)$
378    ///
379    /// $M(n, m, p) = O(n + p + m \log m)$
380    ///
381    /// where $T$ is time, $M$ is additional memory, $n$ is `xs.len()`, $m$ is the sum of the
382    /// significant bits of the elements of `xs` and `ys`, and $p$ is `prec`.
383    ///
384    /// # Panics
385    /// Panics if `prec` is zero or if `xs` and `ys` have different lengths.
386    ///
387    /// # Examples
388    /// ```
389    /// use malachite_base::num::basic::traits::{One, Two};
390    /// use malachite_float::Float;
391    /// use std::cmp::Ordering::*;
392    ///
393    /// let xs = [Float::ONE, Float::TWO, Float::from(3)];
394    /// let ys = [Float::from(4), Float::from(5), Float::from(6)];
395    /// let (dot, o) = Float::dot_prec(&xs, &ys, 10);
396    /// assert_eq!(dot.to_string(), "32.000");
397    /// assert_eq!(o, Equal);
398    ///
399    /// let (dot, o) = Float::dot_prec(&xs, &ys, 3);
400    /// assert_eq!(dot.to_string(), "32.0");
401    /// assert_eq!(o, Equal);
402    /// ```
403    #[inline]
404    pub fn dot_prec(xs: &[Self], ys: &[Self], prec: u64) -> (Self, Ordering) {
405        Self::dot_prec_round(xs, ys, prec, Nearest)
406    }
407
408    /// Computes the dot product of two equal-length slices of [`Float`]s, rounding the result with
409    /// the specified rounding mode. The precision of the result is the maximum of the precisions of
410    /// the inputs (or 1 if there are no inputs). An [`Ordering`] is also returned, indicating
411    /// whether the rounded dot product is less than, equal to, or greater than the exact dot
412    /// product. Although `NaN`s are not comparable to any [`Float`], whenever this function returns
413    /// a `NaN` it also returns `Equal`.
414    ///
415    /// The products are never rounded, and only a single rounding is performed, at the end: the
416    /// result is the correctly-rounded exact dot product, and intermediate overflow and underflow
417    /// cannot occur.
418    ///
419    /// See [`RoundingMode`] for a description of the possible rounding modes.
420    ///
421    /// $$
422    /// f((x_i)_ {i=0}^{n-1}, (y_i)_ {i=0}^{n-1}, m) = \sum_ {i=0}^{n-1} x_i y_i + \varepsilon.
423    /// $$
424    /// - If $\sum_ {i=0}^{n-1} x_i y_i$ is infinite, zero, or `NaN`, $\varepsilon$ may be ignored
425    ///   or assumed to be 0.
426    /// - If $\sum_ {i=0}^{n-1} x_i y_i$ is finite and nonzero, and $m$ is not `Nearest`, then
427    ///   $|\varepsilon| < 2^{\lfloor\log_2 |\sum_ {i=0}^{n-1} x_i y_i|\rfloor-p+1}$, where $p$ is
428    ///   the maximum precision of the inputs.
429    /// - If $\sum_ {i=0}^{n-1} x_i y_i$ is finite and nonzero, and $m$ is `Nearest`, then
430    ///   $|\varepsilon| \leq 2^{\lfloor\log_2 |\sum_ {i=0}^{n-1} x_i y_i|\rfloor-p}$, where $p$ is
431    ///   the maximum precision of the inputs.
432    ///
433    /// See [`Float::dot_prec_round`] for a description of the special cases and of overflow and
434    /// underflow behavior.
435    ///
436    /// If you know you'll be using `Nearest`, consider using [`Float::dot`] instead.
437    ///
438    /// # Worst-case complexity
439    /// $T(n, m) = O(m (n + m))$
440    ///
441    /// $M(n, m) = O(n + m \log m)$
442    ///
443    /// where $T$ is time, $M$ is additional memory, $n$ is `xs.len()`, and $m$ is the sum of the
444    /// significant bits of the elements of `xs` and `ys`.
445    ///
446    /// # Panics
447    /// Panics if `xs` and `ys` have different lengths, or if `rm` is `Exact` and the exact dot
448    /// product is not exactly representable with the maximum of the precisions of the inputs.
449    ///
450    /// # Examples
451    /// ```
452    /// use malachite_base::num::arithmetic::traits::PowerOf2;
453    /// use malachite_base::num::basic::traits::Two;
454    /// use malachite_base::rounding_modes::RoundingMode::*;
455    /// use malachite_float::Float;
456    /// use std::cmp::Ordering::*;
457    ///
458    /// // 0.25 * 0.25 + 2 * 5 = 10.0625, whose inputs have maximum precision 3
459    /// let xs = [Float::power_of_2(-2i64), Float::TWO];
460    /// let ys = [Float::power_of_2(-2i64), Float::from(5)];
461    /// let (dot, o) = Float::dot_round(&xs, &ys, Floor);
462    /// assert_eq!(dot.to_string(), "10.0");
463    /// assert_eq!(o, Less);
464    ///
465    /// let (dot, o) = Float::dot_round(&xs, &ys, Ceiling);
466    /// assert_eq!(dot.to_string(), "12.0");
467    /// assert_eq!(o, Greater);
468    /// ```
469    #[inline]
470    pub fn dot_round(xs: &[Self], ys: &[Self], rm: RoundingMode) -> (Self, Ordering) {
471        Self::dot_prec_round(xs, ys, max_prec(xs.iter().chain(ys.iter())), rm)
472    }
473
474    /// Computes the dot product of two equal-length slices of [`Float`]s. The precision of the
475    /// result is the maximum of the precisions of the inputs (or 1 if there are no inputs), and the
476    /// dot product is rounded to nearest.
477    ///
478    /// The products are never rounded, and only a single rounding is performed, at the end: the
479    /// result is the correctly-rounded exact dot product, and intermediate overflow and underflow
480    /// cannot occur.
481    ///
482    /// $$
483    /// f((x_i)_ {i=0}^{n-1}, (y_i)_ {i=0}^{n-1}) = \sum_ {i=0}^{n-1} x_i y_i + \varepsilon.
484    /// $$
485    /// - If $\sum_ {i=0}^{n-1} x_i y_i$ is infinite, zero, or `NaN`, $\varepsilon$ may be ignored
486    ///   or assumed to be 0.
487    /// - If $\sum_ {i=0}^{n-1} x_i y_i$ is finite and nonzero, then $|\varepsilon| \leq
488    ///   2^{\lfloor\log_2 |\sum_ {i=0}^{n-1} x_i y_i|\rfloor-p}$, where $p$ is the maximum
489    ///   precision of the inputs.
490    ///
491    /// See [`Float::dot_prec_round`] for a description of the special cases and of overflow and
492    /// underflow behavior.
493    ///
494    /// # Worst-case complexity
495    /// $T(n, m) = O(m (n + m))$
496    ///
497    /// $M(n, m) = O(n + m \log m)$
498    ///
499    /// where $T$ is time, $M$ is additional memory, $n$ is `xs.len()`, and $m$ is the sum of the
500    /// significant bits of the elements of `xs` and `ys`.
501    ///
502    /// # Panics
503    /// Panics if `xs` and `ys` have different lengths.
504    ///
505    /// # Examples
506    /// ```
507    /// use malachite_base::num::basic::traits::{One, Two};
508    /// use malachite_float::Float;
509    ///
510    /// let xs = [Float::ONE, Float::TWO, Float::from(3)];
511    /// let ys = [Float::from(4), Float::from(5), Float::from(6)];
512    /// assert_eq!(Float::dot(&xs, &ys).to_string(), "32.0");
513    /// ```
514    #[inline]
515    pub fn dot(xs: &[Self], ys: &[Self]) -> Self {
516        Self::dot_round(xs, ys, Nearest).0
517    }
518}