Incomplete LU decomposition of sparse matrix

Question

I have a sparse matrix stored in CSR format. For this matrix, I would like to get the incomplete LU decomposition. I tried to find algorithms which can utilize the CSR format but I could not find anything. All work, papers, software seem to access the matrix by giving row index and column index like for COO format (A[i,j]).

Yes, this can be done, but I think that doing it this way costs a lot of time for looking for the correct matrix entry in CSR format as one has to iterate through a complete row to find the correct column. Any ideas or papers or git repositories which I haven't found?

If there is no other possibility, what would be the fastest way to access elements in a CSR matrix? Simply iterate through the given row until the defined column is reached?

Answer 1

In terms of implementation, you definitely should take a look at Intel MKL's one. Say, dcsrilu0 and dcsrilut (with threshold) accept the matrix in a DSS format (I linked to a nonsymmetric one), which is a wrapper on top of a three-array CSR storage. (both pages would contain the same matrix as an example).

So, implementation-wise, you should be able to get the ILU from CSR without you reordering yourself.

It also certainly is worth taking a look at SuperLU. However, it works with CCS, natively. Excerpt from the documentation:

The factorization and solve routines in SuperLU are designed to handle column-wise storage only. If the input matrix $A$ is in row-oriented storage, i.e., in SLU_NR format, then the driver routines (dgssv() and dgssvx()) actually perform the LU decomposition on $A^T$, which is column-wise, and solve the system using the $L^T$ and $U^T$ factors. The data structures holding $L$ and $U$ on output are different (swapped) from the data structures you get from column-wise input.

...

Alternatively, the users may call a utility routine dCompRow_to_CompCol() to convert the input matrix in SLU_NR format to another matrix in SLU_NC format, before calling SuperLU.

While the excerpt talks about the full factorize and solve, the same should be true for the incomplete "preconditioner"-related counterparts.