Conversion from coordinate format to CSR (MKL oneAPI in C)

[Dan78]

Hi all,

I would like to convert my matrices which are stored in coordinate format to CSR to be able to call Pardiso. I have investigated a lot and the only solution that I got was to create a matrix with internal format using mkl_sparse_z_create_coo and then export it by calling mkl_sparse_z_export_csr. I have unsymmetric and complex_16 matrix. However, calling the first routine works fine, but calling the second one always returns invalid input parameters. I wanted to know whether there is any better way to directly convert coordinate format to CSR in oneAPI MKL=

Thanks,

Dan

[Gajanan_Choudhary]

Hi @Dan78,

Thanks for reaching out to us about your problem. You are on the right track, but maybe you are missing a call to mkl_sparse_convert_csr() between your calls to mkl_sparse_?_create_coo() and mkl_sparse_?_export_csr(). It would look something like the following assuming some of the variables have been declared with correct types:

sparse_matrix_t cooA, csrA;
status = mkl_sparse_z_create_coo(&cooA, cooIndexBase, cooNrows, cooNcols, cooNnz, cooRowIndex, cooColIndex, cooValues);
/* Do not forget to check "status" and handle any errors here */
status = mkl_sparse_convert_csr(cooA, SPARSE_OPERATION_NON_TRANSPOSE, &csrA);
/* Do not forget to check "status" and handle any errors here again */
status = mkl_sparse_z_export_csr(csrA, &csrIndexBase, &csrNrows, &csrNcols, &csrRowStart, &csrRowEnd, &csrColIndex, &csrValues);
/* Do not forget to check "status" and handle any errors here yet again */
/* Also do not forget to call mkl_sparse_destroy() when you no longer need csrA, cooA*/

Please let us know if that works for you.

Regards,

Gajanan Choudhary

[Dan78]

Hi Gajanan,

Thanks for your help. It solved my problem. Now I have another question. Like you said, calling mkl_sparse_destroy() shall be called to release all memory allocated for internal sparse data type used by MKL. Now I am wondering how memory management is done for csrValues for instance. This variable shall be defined as "MKL_Complex16 *csrValues", so does it mean that mkl_sparse_z_export_csr will allocate internally memory for this to store values in CSR format? Then shall I call "delete[] csrValues" when the matrix is not needed anymore? Actually I did and my code crashed, so I am wondering how this memory management is done? Maybe calling "mkl_sparse_destroy()" will release all memory including csrValue and the memory allocated for row and col start/end and indices?

Regards,

Dan

[shb]

Hi @Dan78,

If you only used csrValues pointer as the output of mkl_sparse_z_export_csr(), then you shouldn't call "delete[] csrValues".

mkl_sparse_?_export_csr() gives users to access those memory objects that are allocated by oneMKL library but still they are library owned memory objects, so you shouldn't release them by yourself.  Refer to below link:

https://www.intel.com/content/www/us/en/docs/onemkl/developer-reference-c/2024-0/mkl-sparse-export-csr.html

In the above example, the "csrA" sparse handle owns the csrValues memory object, so the memory will be released when you call "mkl_sparse_destroy(csrA)".  That's why your code crashed when you explicitly call "delete[] csrValues".

So, the rule of thumb is that you need to release only the memory objects that you allocated explicitly. Internally allocated memory objects will be handled by the oneMKL library.

Please let us know if this helps.  Thanks for using oneMKL library!

Best,

Seung-hee

[Dan78]

Hi,

Any update on this? One more thing I would like to add. When I pass the converted vectors from mkl_sparse_z_csr, to Pardiso solver, the solver always crashes. I check that the returned dimensions are correct and the matrix is correctly 1-based after exporting which matches with the original data in coordinate format. Is there anything I have missed during exporting to CSR format? It would be far easier if Pardiso solver would accept coordinate format but it apparently doesn't. I think there is memory management issues that I have missed in my code so the solver cannot correctly access the data that is set using MKL's export routine.

Regards,

Dan

[JohnNichols]

@mecej4 kindly provided a working solution for making PARDISO work with the translation from a standard matrix to the form required by Pardiso, it is contained in the attached solution.  There are examples. 

This has worked for more than a decade. 

[Dan78]

Hi @JohnNichols 

Thanks but the attached code is written in Fortran while my code is written in C. But that would be OK if I knew where this conversion is done. As far as I can see the a routine called "linsolve" takes care of input arrays before they are used by Pardiso but I didn't find any implementation on that. I imply have a sparse matrix in coordinate format and would like to convert it to CSR format as efficient as possible and then send the data to Pardiso. Intel MKL has routines to convert to internal format, then change the format from COO to CSR and finally export to any other format, here CSR. This works, but Pardiso crashes when I send the exported routines to Pardiso even though the conversion seems to be correct.

I think there are some memory issues that I am not following. That is why I asked how should I release the memory allocated for CSR format and when it is allocated? Is is the export routine that allocates memory? If yes then why calling "delete[]" on the pointer fails?

Regards,

Dan

[Dan78]

Thanks all for your replies and helps.

I actually decided to convert COO format to CSR myself using following lines of code:

    for (auto i = 0; i < nnz; i++) {
        csrVals[i] = cooVals[i];
        csrCols[i] = cooCols[i];
        csrRows[cooRows[i] + 1]++;
    }
    for (auto i = 0; i < rows; i++) {
        csrRows[i + 1] += csrRows[i];
    }

Then I know where to allocate and when to release the memory for CSR format and I Pardiso works fins and doesn't crash anymore.

Intel MKL had a routine to directly convert COO to CSR but this is not deprecated. I think this was a better to do the job instead of converting to an internal format, but you know better. Maybe this is more efficient way.

Regards,

Dan

[mecej4]

Most routines that use the CSR representation follow this convention:

• the entries in row-i must occur after the entries in row-(i-1), if i > 1
• the entries in each row must be sorted in order of increasing column number
• no duplicate column numbers within a row -- you have to group by column number and add up the cooVals associated with the same column number

Your simple code will fail if the input data in COO format has not been sorted and grouped beforehand.

[Dan78]

Hi @mecej4,

Thanks for your reply. So you mean the COO format shall be created by (logically) traversing the original matrix row-wise? Meaning all elements in one row are stored in COO format and then next row, and so on? I think this is the way that I compose the spare matrix in COO format.

Another question is, if I convert from COO to CSR using MKL routines, then does sorting the elements in each row mater, or elements in COO format can then be in any arbitrary order?

Regards,

Dan

[mecej4]

The conversion from COO to CSR is not to be undertaken without paying attention to how the CSR data structure is to be used subsequently.

For instance, if the intention is to pass the newly formed CSR data to the sparse linear solver PARDISO, the documentation stipulates as follows for the matrix nonzero values argument :

The matrix must be stored in the three-array variant of the compressed sparse row (CSR3) or in the three-array variant of the block compressed sparse row (BSR3) format, and the matrix must be stored with increasing values of ja for each row.

Note that not only is ordering required but there can be only one entry for each column within the present row.

The COO format is often used with finite element discretization, and so the entries can be in any order, and there can be multiple non-contiguous entries corresponding to a specific row, column pair. For example, if four quadrilaterals have a common corner in the FEA mesh, there will be four contributions from the elements to the stiffness matrix elements related to that corner node.

[VarshaS_Intel]

Hi Dan,

Has the information provided by the Black Belt user mecej4 helped? Could you please provide us with an update on your issue?

Thanks & Regards,

Varsha

[Dan78]

Hi Varsha,

Yes it was very useful. Actually I create COO format by traversing row-wise the original matrix, so I think the CSR format is correctly created. No I observe that the residual of the solution returned by Pardiso is very high. I am going to try on different use cases which are all within EM simulation solving matrices from integral equations, so the matrices are ill-conditioned. I think this thread can be closed. Thanks all for your kind help.

Regards,

Dan

[VarshaS_Intel]

Hi Dans,

>>Yes it was very useful. I think this thread can be closed.

Thanks for the explanation and confirmation. It’s great to know that the issue has been resolved, in case you run into any other issues please feel free to create a new thread.

Have a Good Day!

Thanks & Regards,

Varsha