pub fn solve_linear_system(a_flat: &[f64], b: &[f64], n: usize) -> Vec<f64>
Solve A*x = b via Gaussian elimination with partial pivoting.