arithmetic_analysis/lu_decomposition is wrong

[spamegg1]

The provided example matrix in the file is

[[2, -2, 1], 
[0, 1, 2], 
[5, 3, 1]]

The algorithm returns

L = [[1.  0.  0. ]
 [0.  1.  0. ]
 [2.5 0.  1. ]]
U = [[ 2.  -2.   1. ]
 [ 0.   1.   2. ]
 [ 0.   8.  -1.5]]

Multiplying L and U does give us back the original matrix.

However as you can see the matrix U is NOT upper-triangular.

Wikipedia says that a matrix admits LU-decomposition if and only if all its principal minors (in this case, the determinants of all the 2x2 matrices obtained by removing a row and a column) are nonzero. Upon inspection we can see that this holds true of the original matrix, so it should definitely admit LU-decomposition.

An online calculator gives

L = [[1.  0.  0. ]
       [0.  1.  0. ]
       [2.5 8.  1. ]]
U = [[ 2.  -2.   1. ]
       [ 0.   1.   2. ]
       [ 0.   0.  -17.5]]

Just to be sure I tried it on another matrix (it satisfies the principal minors condition):

[[1, 2, 3], 
[4, 5, 6], 
[7, 8, 9]]

The algorithm gives

L = [[1.  0.  0. ]
       [0.  1.  0. ]
       [7. 0.  1. ]]
U = [[ 1.  2.   3. ]
       [ 4.   5.   6. ]
       [ 0.   -6.  -12.]]

Once again U is not upper-triangular.

The online calculator gives

L = [[1.  0.  0. ]
       [4.  1.  0. ]
       [7. 2.  1. ]]
U = [[ 1.  2.   3. ]
       [ 0.   -3.   -6. ]
       [ 0.   0.  0.]]

I'm not sure how to fix this algorithm, so help is definitely wanted.