pub fn solve_banded<const C: usize>(
a: Vec<Vec<Scalar>>,
first: &[usize],
last: &[usize],
rhs: Vec<Vector<C>>,
) -> GeopResult<Vec<Vector<C>>>Expand description
Solves a x = rhs for the banded m x m matrix a and C
right-hand sides, by Gaussian elimination without pivoting: row k is
nonzero only in the columns first[k] ..= last[k], and the bands never
move left down the rows. Safe without pivoting only where the matrix
says so — a B-spline collocation matrix is totally positive (Piegl &
Tiller §9.2.1), the reason the curve interpolation that builds it does
not pivot either.
Elimination without pivoting keeps the band, so it runs within it:
below a pivot only the rows whose band reaches its column, and in them
only the columns the pivot row’s band reaches. Outside the band every
factor and every entry is an exact zero, and subtracting zero times
zero changes nothing but the enclosure — outward rounding widened each
entry by an ulp per pivot above it, and the dense elimination took
m³ / 3 steps for what the band holds in m (degree + 1)².