Implement an in-place PAQ = LU decomposition with partial pivoting and simulated row exchanges.
Requirements and notes:
- The algorithm must perform both row and column exchanges so that pivot columns are moved before non-pivot columns. Column exchanges are virtual: only record the column ordering (a permutation vector
Q), do not physically permute the columns ofA. - Row exchanges must be handled by a permutation vector
P(simulated row exchanges). Do not physically swap rows inA; instead maintain and use the permutationPwhen accessing rows. The matrixAitself will be overwritten in-place and will ultimately contain theLandUfactors (except the unit diagonal ofL, which is implicit). - The implementation must work for rectangular matrices
A(m × n). To derive the algorithm, imagine padding the matrix with zero rows or columns until it is square and then applying the usual square PAQ = LU algorithm. These added zeros are only conceptual: your code must never access entries outside the original matrix bounds — treat the padding as implicit. - The final outputs should include:
- The overwritten matrix
Acontaining theL(strict lower triangle, with implicit ones on the diagonal) andU(upper triangle). - A row permutation vector
Pthat encodes the simulated row exchanges. - A column permutation vector
Q(or equivalent structure) that records the virtual column exchanges / pivot column ordering. - Any additional information needed to identify which columns are pivots (useful for Problem 2).
- The overwritten matrix
Implementation pointers:
- Use partial pivoting: at each step, choose the pivot row (maximum magnitude entry in the current column among available rows), update
Paccordingly (simulated swap), and apply Gaussian elimination on the active submatrix while recording column pivots inQ. - When the matrix is rectangular, only iterate over the valid rows and columns; do not index beyond
morn. - After completion,
Uwill appear in the upper triangular part of the storedA, andLin the strictly lower triangular part (with diagonal ones implied).
Write a solver for linear systems A x = b that uses the PAQ = LU decomposition from Problem 1.
Key ideas:
- After PAQ = LU is computed, pivot columns of
Ucorrespond to basic variables; non-pivot columns correspond to free variables. - The general solution can be written as
x = N x_free + c
where
x_freecontains the free-variable values,Nis a matrix mapping free variables to the full solution vector, andcis a particular solution vector. The separate LaTeX filederivation.texin this folder contains the derivation ofNandc. - Implementation outline:
- Compute PAQ = LU for
A, obtainingAoverwritten withLandU, and permutation vectorsPandQ. - Apply the row permutation
Ptob(simulated via indexing) and solve Ly = Pb by forward substitution, taking into account the implicit unit diagonal ofL. - Identify pivot columns in
U(basic variables). Partition unknowns into basic and free sets according toQ. - Solve U_basic x_basic = y_basic by back substitution for the basic variables; express the dependent variables in terms of free variables using the structure of
U(this yieldsNandc). - Reconstruct the full solution
xin the original column order usingQ(reverse the virtual column permutation).
- Compute PAQ = LU for
- Your solver should support returning:
- One particular solution
c(e.g., with free variables set to zero). - The matrix
N(so all solutions can be generated for arbitrary free-variable choices). - Optionally: a parameterization function that takes
x_freeand returnsx.
- One particular solution
See derivation.tex for the detailed algebra deriving N and c.