LU Decomposition: difference between between hand calculation and solver? I have a $3 \times 3$ matrix $A$ and have to perform the $LU$ Factorization (1)
$$ A = \begin{bmatrix}
1 & 2 & 3 \\ 
1 & -1 & 3 \\
-2 & -10 & -2
\end{bmatrix}$$
Using row reduction, I would first substract the second row by (1/1) the first row which gives (2): 
$$ \begin{bmatrix}
1 & 2 & 3 \\ 
0 & -3 & 0 \\
-2 & -10 & -2
\end{bmatrix}$$
I would then substract the third row by (-2/1) the first row which gives (3)
$$ \begin{bmatrix}
1 & 2 & 3 \\ 
0 & -3 & 0 \\
0 & -6 & 4
\end{bmatrix}$$
Finally I would substract the third row by (-6/2) the second row which would give (4): 
$$ \begin{bmatrix}
1 & 2 & 3 \\ 
0 & -3 & 0 \\
0 & 0 & 4
\end{bmatrix} = U$$
with the factor found in (2) (3) and (4) I can get the matrix $U$
$$L =  \begin{bmatrix}
1 & 0 & 0 \\ 
1 & 1 & 0 \\
-2 & 2 & 1
\end{bmatrix}$$
I can verify that $L \times U = A$
However when I compute [L, U] = lu(A) with Octave, I get 
$$ U = 
\begin{bmatrix}
-2 & -10 & -2 \\ 
0 & -6 & 2 \\
0 & 0 & 1
\end{bmatrix}$$
and 
$$ L= 
\begin{bmatrix}
-0.5 & 0.5 & 1 \\ 
-0.5 & 1 & 0 \\
1 & 0 & 0
\end{bmatrix}$$
Here is the octave / matlab code: 
A = [1 2 3; 1 -1 3; -2 -10 -2];
U_hand = [1 2 3; 0 -3 0; 0 0 4];
L_hand = [1 0 0; 1 1 0; -2 2 1];

L_hand * U_hand

[L, U] = lu(A)

How can how explain the differences? I'm probably wrong somewhere but where?
 A: You shouldn't have got that for your LU decomp. I used python which uses the same LAPACK
import scipy.linalg
import 

A = scipy.array([[1 ,2,3],[1, -1, 3 ] ,[-2,-10,-2]])
P,L,U = scipy.linalg.lu(A)


L
Out[3]: 
array([[ 1. ,  0. ,  0. ],
       [-0.5,  1. ,  0. ],
       [-0.5,  0.5,  1. ]])

U
Out[4]: 
array([[ -2., -10.,  -2.],
       [  0.,  -6.,   2.],
       [  0.,   0.,   1.]])

Upon further inspection. The reference says
It's the product of the permutation matrix and the L matrix 

With two output arguments, returns the permuted forms of the upper and
  lower triangular matrices, such that A = L * U. With one output
  argument y, then the matrix returned by the LAPACK routines is
  returned. If the input matrix is sparse then the matrix L is embedded
  into U to give a return value similar to the full case. For both full
  and sparse matrices, lu loses the permutation information.

Which leads to testing it in Python.
Atest = np.dot(P,L)

Out[2]: 
array([[-0.5,  0.5,  1. ],
       [-0.5,  1. ,  0. ],
       [ 1. ,  0. ,  0. ]])

