pub fn solve_linear(a: &Matrix, b: &[f64]) -> Option<Vec<f64>>
Solve Ax = b via LU decomposition. Returns None if singular.