pub const fn column_basis<const R: usize, const C: usize>(
rows: &[[f64; C]; R],
tolerance: Tolerance,
) -> Result<ColumnBasis<R, C>, LaError>Expand description
Select columns of a rectangular array by partial-pivoted elimination.
Columns are visited from left to right. At each step the largest absolute
entry among remaining rows is the candidate pivot; ties select the first
row in the current row order. A pivot is accepted only when its magnitude
is strictly greater than tolerance.get(). Rejected columns are skipped;
selection ends after all columns or R accepted pivots. Empty and zero
matrices return rank zero; tall, square, and wide shapes are supported.
The tolerance is absolute, in input coefficient units. There is no automatic
scaling or normalization. To use a relative policy, compute a finite absolute
threshold from the original input scale and pass it through
Tolerance::try_new. Rescaling rows or columns can change the selection.
Elimination uses partial row pivoting and fused multiply-add updates
(Gaussian elimination, REFERENCES.md [11, 12]). No absolute rounding-error
bound is provided. Rounding, underflow, growth, and ill-conditioning can alter
the selected columns and rank even at zero tolerance. The result is provisional:
exact certification of a square selected basis remains the caller’s job.
Failure to find a pivot never produces an exact singularity claim.
The input is borrowed and preserved. A local R × C copy and R indices
use stack storage; no allocation or optional feature is required. Dimensions
have no dispatch cap. Work is O(R (1 + C (1 + min(R, C)))), including
validation and initialization for empty shapes.
§Errors
All inputs are checked in row-major order before elimination, including
columns beyond a full basis. NaN/infinity returns LaError::NonFinite
with input origin and the original matrix cell. A non-finite update returns
the same error with ArithmeticOperation::ColumnBasis and its cell in the
row-permuted work array. Insufficient numerical rank is a successful partial
selection, not an error. No partially computed result is returned on error.
§Examples
use la_stack::prelude::*;
let rows = [[0.0, 1.0, 2.0, 0.0], [0.0, 2.0, 4.0, 1.0]];
let basis: ColumnBasis<2, 4> = column_basis(&rows, Tolerance::try_new(1e-12)?)?;
assert_eq!(basis.columns(), &[1, 3]);
assert_eq!(basis.rank(), 2);