Definition of LU descomposition of a matrix depends on the base ring

[tmonteil]

As reported on ​this ask question, there is an inconsistency in the definition of the LU decomposition of a matrix depending on its base ring. If a matrix A belongs to ZZ, QQ, AA, QQbar, A.LU() returns a triple (P,L,U) such that A=PLU. If a matrix A belongs to RDF, A.LU() returns a triple (P,L,U) such that PA=LU. For example:

sage: A = random_matrix(ZZ,4)
sage: A.LU()
(
[0 1 0 0]  [    1     0     0     0]  [   72    -4   -38     0]
[0 0 1 0]  [    0     1     0     0]  [    0    -1    -2    10]
[1 0 0 0]  [ 1/36   8/9     1     0]  [    0     0  17/6 -80/9]
[0 0 0 1], [-1/36   1/9    -1     1], [    0     0     0   -17]
)
sage: B = A.change_ring(RDF)
sage: B.LU()
(
[0.0 0.0 1.0 0.0]
[1.0 0.0 0.0 0.0]
[0.0 1.0 0.0 0.0]
[0.0 0.0 0.0 1.0],

[             1.0              0.0              0.0              0.0]
[             0.0              1.0              0.0              0.0]
[ 0.0277777777778   0.888888888889              1.0              0.0]
[-0.0277777777778   0.111111111111             -1.0              1.0],

[          72.0           -4.0          -38.0            0.0]
[           0.0           -1.0           -2.0           10.0]
[           0.0            0.0  2.83333333333 -8.88888888889]
[           0.0            0.0            0.0          -17.0]
)


sage: B.LU()[0] == A.LU()[0]
False
sage: B.LU()[0] == A.LU()[0].transpose()
True
sage: 

The aim of this ticket is to fix this inconsistency and choose a common definition for all rings.

[tmonteil]

It is worth noticing that the source code for matrices over RDF calls scipy.linalg.lu which returns a triple (P,L,U) such that A=PLU and then transpose P with the documentation:

# Numpy has a different convention than we had with GSL
# So we invert (transpose) the P to match our prior behavior
# TODO: It's an awful waste to store a huge matrix for P, which
# is just a simple permutation, really.

So, i guess this extra transposition should be reverted, so that it becomes consistent with both scipy and the implementation for exact rings.