pub fn dot(a: &[f64], b: &[f64], backend: ComputeBackend) -> f64Expand description
Inner product sum_i a_i * b_i, on the chosen backend.
Krylov methods evaluate two or three of these per iteration (conjugate
gradient: r.r and p.Ap), so on a large mesh they are worth threading even
though the operation is memory-bandwidth bound rather than arithmetic bound.
§Arguments
a,b— equal-length dimensionless vectors in cell order.backend— requested execution backend; seevecop_backend_forfor what will actually run.
§Returns
The dimensionless scalar product. Returns 0.0 for empty inputs.
§Determinism, and how this differs from crate::krylov::vecops::dot
Sums in fixed blocks of REDUCTION_BLOCK elements, combining block partials
in ascending block order. The result is therefore bitwise identical between
backends and at any thread count, and reproducible run to run. It differs
from the flat left-to-right sum in crate::krylov::vecops::dot only by
floating-point non-associativity.
Measured deviation. Methodology: two fixed-seed xorshift64*
pseudorandom vectors with elements uniform on [-1, 1), at lengths 1 024
through 4 194 304 in powers of four; compare against
crate::krylov::vecops::dot on the same inputs. Two measures are reported,
because the raw relative difference is taken against a heavily cancelled sum
(terms of order 1, total of order sqrt(n)) and so overstates the error:
the raw relative difference |blocked - flat| / max(|blocked|, |flat|), and
the conditioning-aware |blocked - flat| / sum |a_i b_i|. Pass criteria:
worst raw <= 1e-12 and worst conditioned <= 1e-15.
Results, measured 2026-08-12 by the test dot_matches_flat_reference in
parallel/tests.rs:
| n | blocked | flat | raw rel diff | vs sum abs terms |
|---|---|---|---|---|
| 1 024 | -1.07439412213983019e1 | -1.07439412213983019e1 | 0 | 0 |
| 4 096 | 2.27305615265860794e0 | 2.27305615265858707e0 | 9.1824e-15 | 2.0525e-17 |
| 16 384 | -4.10927502452607953e1 | -4.10927502452608877e1 | 2.2479e-15 | 2.2574e-17 |
| 65 536 | -1.48877508980721871e2 | -1.48877508980720563e2 | 8.7817e-15 | 7.9476e-17 |
| 262 144 | -2.57469059787491972e2 | -2.57469059787487652e2 | 1.6779e-14 | 6.6111e-17 |
| 1 048 576 | 2.06465247477648091e2 | 2.06465247477646983e2 | 5.3687e-15 | 4.2292e-18 |
| 4 194 304 | -1.50102185250877881e2 | -1.50102185250842240e2 | 2.3744e-13 | 3.3997e-17 |
Worst raw relative difference 2.3744e-13 at n = 4 194 304; worst conditioned difference 7.9476e-17 at n = 65 536.
Interpretation: the difference is summation reassociation at rounding
level, not a defect — relative to the arithmetic actually performed it never
exceeds about a third of one f64 epsilon. The raw figure grows with n
because the exact sum of n random signed terms grows only like sqrt(n)
while the rounding error grows faster, so the quotient inflates; that is a
property of the test vectors, not of either kernel. Blocked summation is the
two-level, more accurate form — a flat sum of n terms has error growing
like n * eps, a blocked sum like (block + n / block) * eps — so the
blocked result is if anything closer to the exact value.
§Panics
Panics if a and b have different lengths.
§Example
use outram_foam_basic_lib::compute::ComputeBackend;
use outram_foam_basic_lib::ldu_matrix::parallel::dot;
let a = [1.0, 2.0, 3.0];
let b = [4.0, -5.0, 6.0];
// 1*4 + 2*(-5) + 3*6 = 12
assert_eq!(dot(&a, &b, ComputeBackend::Serial), 12.0);
// Bitwise identical, not merely close.
assert_eq!(
dot(&a, &b, ComputeBackend::CpuMulti),
dot(&a, &b, ComputeBackend::Serial),
);