Intel® oneAPI Math Kernel Library
Ask questions and share information with other developers who use Intel® Math Kernel Library.
Announcements
Welcome to the Intel Community. If you get an answer you like, please mark it as an Accepted Solution to help others. Thank you!
6429 Discussions

question about MKL matrix-matrix multiplier.

SUKLA
Beginner
726 Views

I am using Intel MKL (compilers_and_libraries_2020.4.311).  "mkl_sparse_d_mm" produces wrong result while "mkl_dcsrmv" and "mkl_sparse_d_mv" produces the correct result. Can anybody tell me what I am doing wrong? I want to use "SPARSE_LAYOUT_COLUMN_MAJOR" for dense matrix. Basically MKL is accessing B matrix differently that I intended to use.

Sukla

 

 


#include <sstream>
#include <iostream>

#include <mkl.h>
#include <stdio.h>
#include <string.h>

#include <sstream>

int main(int argc, char ** argv)
{
MKL_INT mkl_nrow = 2;
MKL_INT mkl_ncol = 3;
MKL_INT mkl_numVec = 2;
MKL_INT amklIaPtr[4] = {0,2,4,6};
MKL_INT amklJaPtr[6] = { 0,1,0,1,0,1 };


double aValuePtr[6] = { 1.,2.0 ,0.0 ,5.,-1,-2.};
double bValuePtr[] = { 1.,2.0,5.,2.,3.,-5 ,-.1,-0.2,0.3};
double resultProd[10] ;

double alpha = 1.0, beta = 0.;
char matdescra[6] = { 'S','U','N','C',' ',' ' };
matdescra[0] = 'G';
char transa='T';
for (int icol = 0; icol < mkl_numVec; icol++) {
mkl_dcsrmv(&transa, &mkl_ncol, &mkl_nrow, &alpha, matdescra, aValuePtr, &amklJaPtr[0],
&amklIaPtr[0], &amklIaPtr[1], &bValuePtr[icol*mkl_ncol],
&beta, &resultProd[icol*mkl_nrow]);
}
// ************ CORRECT ANSWER *********************
// resultProd[4] = { -4.0,2.0,7.0,29.0 };
// ************ CORRECT ANSWER *********************

sparse_matrix_t aMatrix;
sparse_index_base_t indexing = SPARSE_INDEX_BASE_ZERO;
sparse_status_t status = mkl_sparse_d_create_csc(
&aMatrix, indexing, /* indexing: C-style or Fortran-style */
mkl_nrow, mkl_ncol,
amklIaPtr, &amklIaPtr[1], amklJaPtr, aValuePtr);
sparse_operation_t spOp = SPARSE_OPERATION_NON_TRANSPOSE;
matrix_descr descr = { SPARSE_MATRIX_TYPE_GENERAL, SPARSE_FILL_MODE_LOWER, SPARSE_DIAG_NON_UNIT };
for(int icol=0; icol< mkl_numVec; icol++){
status = mkl_sparse_d_mv(
spOp,
alpha,
aMatrix,
descr, /* sparse_matrix_type_t + sparse_fill_mode_t + sparse_diag_type_t */
&bValuePtr[icol*mkl_ncol],
beta,
&resultProd[icol*mkl_nrow]);
}

// ************ CORRECT ANSWER *********************
// resultProd[4] = { -4.0,2.0,7.0,29.0 };
// ************ CORRECT ANSWER *********************

status = mkl_sparse_d_mm(
spOp,
alpha,
aMatrix,
descr, /* sparse_matrix_type_t + sparse_fill_mode_t + sparse_diag_type_t */
SPARSE_LAYOUT_COLUMN_MAJOR, /* storage scheme for the dense matrix: C-style or Fortran-style */
bValuePtr,
mkl_numVec,
mkl_ncol, //---> Something Wrong here
beta,
resultProd,
mkl_nrow);
// ************ WRONG ANSWER *********************
// resultProd[4] = { 1.1,2.2,12.2,19.4 };
// ************ WRONG ANSWER *********************

return 0;

0 Kudos
9 Replies
mecej4
Black Belt
710 Views

I have not inspected every line of your code, but it strikes me that the problem you have is that the sparse matrix represented is the transpose of A, rather than A itself. If mkl_nrow = 2, the arrays (using notation in MKL Ref. Man.) pntre[] and pntrb[] should each have 2 entries. The representation is required to be the CSR representation of A, not of its transpose.

SUKLA
Beginner
701 Views

Mecej,

    However, I am using mkl_sparse_d_create_csc instead of mkl_sparse_d_create_csr. I do not think I need to do transpose here. Please help.

 

SUKLA
Beginner
699 Views

Furthermore, MatrixVector multiplication works without transposing.

 

mecej4
Black Belt
686 Views

Sorry, I overlooked the part where you used  mkl_sparse_d_create_csc rather than the _csr version.

GouthamK_Intel
Moderator
646 Views

Hi Sukla,

We are able to reproduce the issue which you are facing on our end. We are escalating this thread to a Subject Matter Expert(SME) who will guide you further.

Have a Good day!


Thanks & Regards

Goutham


SUKLA
Beginner
613 Views

Goutam,

   Let me know when I can get new MKL library to test it. I have another question about EFAST. Does the B matrix have to be positive definite? Can it not be semi-positive definite? I can not to get eigenvalue from FEAST if the B matrix is semi-positive definite (I have zero row and columns for particulation degree-of-freedom in B matrix.).

Sukla

Gennady_F_Intel
Moderator
558 Views

B has to be positive devined


Gennady_F_Intel
Moderator
252 Views

The cause of the issue has been investigated:

The status returned by the mkl_sparse_d_mm() call is ==  SPARSE_STATUS_SUCCESS even though it should return SPARSE_STATUS_NOT_SUPPORTED.


Please check the documentation - https://software.intel.com/content/www/us/en/develop/documentation/onemkl-developer-reference-c/top/...


You called SPARSE_INDEX_BASE_ZERO and CSC format but this configuration doesn’t support.


Gennady_F_Intel
Moderator
236 Views

The issue is closing and we will no longer respond to this thread. If you require additional assistance from Intel, please start a new thread. Any further interaction in this thread will be considered community only.



Reply