Description | MA is decomposed into a product MA = L * U, where L is lower-triangular with the diagonal elements equal to 1, and U is upper-triangular. As the combined number of non-trivial elements of L and U just fit into a matrix of the same dimensions as MA, the result is stored in a single matrix MLU rather than in two distinct matrices L and U. Actually, it is not the "true" matrices L and U which are combined into MLU, but rather a row-wise permutation, whose indices are output in the vector Ind.
MA may or may not be overwritten by MLU.
There are applications where it makes sense to make near-singular matrices decomposable by "pivot editing": If, in the partial pivoting process employed in the decomposition and inversion, no diagonal element larger than the threshold is available, the scaling of the matrix is done with the threshold value rather than with the tiny diagonal element. That way, divisions by (near-)zero are avoided.
By default, this pivot editing is switched off. If you wish to switch it on for all calls to MF_LUdecompose and to the functions based upon it, namely MF_inv and MF_solve, you can do so by calling MF_LUDsetEdit. However, as this method is not thread-safe, you cannot use it in order to set different thresholds for different calls to the functions mentioned. Instead, use their "wEdit" variants, i.e. MF_LUdecomposewEdit, MF_invwEdit and MF_solvewEdit. They take the desired threshold as the additional argument thresh. Note that thresh is always real, also in the complex versions.
The return value indicates if the factorization was successful and if the number of permuations was even or uneven:
Return value | Meaning |
0 | Matrix is singular; result cannot be used |
+1 | Factorization successful; even number of permutations |
-1 | Factorization successful; uneven number of permutations |
+2 | Matrix (nearly) singular; could be factorized by pivot editing; result partially usable; even number of permutations |
-2 | Matrix (nearly) singular; could be factorized by pivot editing; result partially usable; uneven number of permutations |
To check if MF_LUdecompose was successful, in single-thread programs, you may also call MF_LUDresult, whose return value will be FALSE (0), if the MA could be decomposed, and TRUE (1) for singular MA. In multi-thread programs, on the other hand, it would not be clear wich instance of MF_LUdecompose the call to MF_LUDresult would refer to. So, here, inspection of the return value of MF_LUdecompose is the only option.
The LU decomposition functions need buffer memory. The "normal" versions (prefixes MF_, MCF_ etc.) allocate it themselves, whereas the version with the prefixes MFb_, MCFb_ etc. take a vector Buf as an additional argument. The required size of Buf can be obtained by calling MFb_LUdecompose_sizeBuf() etc., whereby the size is given as the number of elements of Buf in the respective data type (rather than in bytes). |