Skip to main content

primitive_float_dot

Function primitive_float_dot 

Source
pub fn primitive_float_dot<T>(xs: &[T], ys: &[T]) -> T
where Float: From<T> + PartialOrd<T>, for<'a> T: ExactFrom<&'a Float> + PrimitiveFloat,
Expand description

Computes the dot product of two equal-length slices of primitive floats, with a single rounding.

The result is correctly rounded to the nearest value: the products are exact, the sum is computed as if in infinite precision, and only a single rounding is performed, at the end. This includes gradual underflow: results in the subnormal range are correctly rounded to their reduced precisions. Intermediate overflow and underflow cannot occur.

$$ f((x_i)_ {i=0}^{n-1}, (y_i)_ {i=0}^{n-1}) = \sum_ {i=0}^{n-1} x_i y_i + \varepsilon. $$

  • If $\sum_ {i=0}^{n-1} x_i y_i$ is infinite, zero, or NaN, $\varepsilon$ may be ignored or assumed to be 0.
  • If $\sum_ {i=0}^{n-1} x_i y_i$ is finite and nonzero, then $|\varepsilon| \leq 2^{\lfloor\log_2 |\sum_ {i=0}^{n-1} x_i y_i|\rfloor-p}$, where $p$ is the precision of the output (typically 24 if T is a f32 and 53 if T is a f64, but less if the output is subnormal).

See Float::dot_prec_round for a description of the special cases, which follow the rules of multiplication for each term and the rules of addition for their combination.

If the result overflows, $\pm\infty$ is returned, and if it underflows, $\pm0.0$ is returned.

§Worst-case complexity

$T(n) = O(n)$

$M(n) = O(n)$

where $T$ is time, $M$ is additional memory, and $n$ is xs.len(): the products are constant-size, and a primitive float’s exponent range is bounded, so the summation window is repositioned only a constant number of times.

§Panics

Panics if xs and ys have different lengths.

§Examples

use malachite_base::num::float::NiceFloat;
use malachite_float::float::arithmetic::dot::primitive_float_dot;

// A naive fold overflows on the first product; the correctly-rounded dot product does not.
let xs = [1.0e300f64, 1.0e300];
let ys = [1.0e300f64, -1.0e300];
assert_eq!(NiceFloat(primitive_float_dot(&xs, &ys)), NiceFloat(0.0));