# How do I extract the output of Aasen's algorithm into a usable form?

I tried implementing the algorithm in Aasen's 1971 paper on factorizing symmetric indefinite matrices. I've translated the code verbatim from Algol into Python, and I used the test example given in the paper. Their print-out from the test example matches mine. I've tried to extract the output into a usable form: A tridiagonal matrix T, a lower unit triangular matrix L, and a permutation P. The problem is that when I execute pprint(L * T * L.T) I get nothing like the input matrix M. How do I extract the output into a usable form?

Anyway, the input is:

⎡ -1   -0.2  0.5   0.2   -0.7⎤
⎢                            ⎥
⎢-0.2   1    -0.8   0    0.8 ⎥
⎢                            ⎥
⎢0.5   -0.8  0.4   -0.9  0.9 ⎥
⎢                            ⎥
⎢0.2    0    -0.9   1    0.5 ⎥
⎢                            ⎥
⎣-0.7  0.8   0.9   0.5   0.8 ⎦


The final state of M is the same as in the paper:

⎡        -1                -0.2                0.5                0.2
⎢
⎢       -0.7                0.8               -0.8                 0
⎢
⎢-0.714285714285714  1.47142857142857   2.09387755102041          -0.9
⎢
⎢0.285714285714286   0.388349514563107  -1.4621359223301    1.42801395041946
⎢
⎣-0.285714285714286  0.495145631067961  0.792828685258964  -0.340239043824701

-0.7      ⎤
⎥
0.8       ⎥
⎥
0.9       ⎥
⎥
0.5       ⎥
⎥
1.62752019809209⎦


and then finally, I attempt to convert the final state of M into a usable form: L becomes

⎡        1                   0                  0          0  0⎤
⎢                                                              ⎥
⎢        0                   1                  0          0  0⎥
⎢                                                              ⎥
⎢-0.714285714285714          0                  1          0  0⎥
⎢                                                              ⎥
⎢0.285714285714286   0.388349514563107          0          1  0⎥
⎢                                                              ⎥
⎣-0.285714285714286  0.495145631067961  0.792828685258964  0  1⎦


and T becomes

⎡ -1         -0.7               0                  0                   0
⎢
⎢-0.7        0.8         1.47142857142857          0                   0
⎢
⎢ 0    1.47142857142857  2.09387755102041   -1.4621359223301           0
⎢
⎢ 0           0          -1.4621359223301   1.42801395041946   -0.340239043824
⎢
⎣ 0           0                 0          -0.340239043824701   1.627520198092

⎤
⎥
⎥
⎥
⎥
⎥
701⎥
⎥
09 ⎦


[I've removed the code because I thought it would confuse things.]

The elements of L need to be shifted rightward. That is, this:

⎡        1                   0                  0          0  0⎤
⎢                                                              ⎥
⎢        0                   1                  0          0  0⎥
⎢                                                              ⎥
⎢-0.714285714285714          0                  1          0  0⎥
⎢                                                              ⎥
⎢0.285714285714286   0.388349514563107          0          1  0⎥
⎢                                                              ⎥
⎣-0.285714285714286  0.495145631067961  0.792828685258964  0  1⎦


needs to become this

⎡1          0                   0                  0          0⎤
⎢                                                              ⎥
⎢0          1                   0                  0          0⎥
⎢                                                              ⎥
⎢0  -0.714285714285714          1                  0          0⎥
⎢                                                              ⎥
⎢0  0.285714285714286   0.388349514563107          1          0⎥
⎢                                                              ⎥
⎣0  -0.285714285714286  0.495145631067961  0.792828685258964  1⎦


This wasn't clear from the original paper.