# LUDecomposition for non-square matrices

How do we go about finding the LUDecomposition for non-square matrices. When i try to input the standard LUDecomposition command, it tells me that the matrix is not-invertible (understandable). Is there still a way to somehow find the LUDecomposition?

This is the matrix:

a = {{2, 4, -1, 5, -2}, {-4, -5, 3, -8, 1}, {2, -5, -4, 1, 8}, {-6, 0, 7, -3, 1}}

This is what i did:

{lu, p, c} = LUDecomposition[a]

Result:

LUDecomposition::sing: "Matrix {{2,4,-1,5,-2},{-4,-5,3,-8,1},{2,-5,-4,1,8},{-6,0,7,-3,1}} is singular"

{{2, 4, -1, 5, -2}, {-2, 3, 1, 2, -3}, {1, -3, 0, 2, 1}, {-3, 4, 0, 2, 5}}, {1, 2, 3, 4}, 0}

The output data is different from the input data so i figured i may be able to still go ahead and use the upper and lower triangularize commands, but i get an error saying that Mathematica can't combine objects of unequal length.

(This answer is incomplete, but I don't know if I'll get time to finish it soon.)

LUDecomposition nominally requires its input matrix to be square, as per the documentation. As far as I can see, the output in this case is correct except in the last row:

{lu, a, b} =
LUDecomposition[{{2, 4, -1, 5, -2}, {-4, -5, 3, -8, 1},
{2, -5, -4, 1, 8}, {-6, 0, 7, -3, 1}}];

ReplacePart[
LowerTriangularize[lu[[All, 1 ;; -2]]],
{i_, i_} -> 1] .
UpperTriangularize[lu] // MatrixForm

Outputs $\left( \begin{array}{ccccc} 2 & 4 & -1 & 5 & -2 \\ -4 & -5 & 3 & -8 & 1 \\ 2 & -5 & -4 & 1 & 8 \\ -6 & 0 & 7 & -5 & -1 \\ \end{array} \right)$.

I'm not sure what the fix is.

A wasteful method:

mat = {{2, 4, -1, 5, -2}, {-4, -5, 3, -8, 1},
{2, -5, -4, 1, 8}, {-6, 0, 7, -3, 1}};
{m, n} = Dimensions[mat];
{lu, piv} = Most[LUDecomposition[mat]] // Quiet;
lm = Drop[LowerTriangularize[lu, -1], None, -1] + IdentityMatrix[m];
um = LinearSolve[lm, mat];
lm.um == IdentityMatrix[m][[piv]].mat
True