I have implement the PLU algorithm for matrix decomposition. The computational complexity is now reduced to n^3 for both the inversion and the computation of the determinant (while it was n! before). I have tested the class for matrices of up to 300 x 300 and the inversion works fine and FAST (though I haven't checked if the solutions were accurate for big matrices).
The PLU decomposition can be accessed by the function $this->lu(). See the example file.
While programming this, I have discovered a couple of bugs on secondary functions such as get_columns, which where corrected. I also discovered more significant bugs that could arise if matrices had a full column/row of zeros. Then, it was likely that at some point in a sequence of computations, the dimensions would go wrong.
I have corrected this by overloading the matrix constructor, which can now can one or three arguments. The one argument version is the original one, while the three arguments version takes the number of rows and the number of columns as aditionnal information.
Finally, I have added a file of the historic of the changes in the package.
As usual, feel free to comment.