Sparse Matrix mkl_?csrmultd problem

[forrest_R_]

Hello,

I want to use the mkl_?csrmultd to do 2 matrix product and the output is a dense matrix.

But i am confused when i read the manual, The ldc (leading dimension of dense matrix C) is a output  parameters( not a input parameter as usual ?), and the length of ib is m+1. The definition of ia is also very different with other sparse matrix routines, because the length is not m+1. 

Is there any problem about this part?

Any help and comment will be appreciated.

Jiajun

[mecej4]

I agree that 'ldc' is incorrectly described as an output argument to mkl_?csr_multd(). Rather, the matrix C is the output argument, and ldC is its leading dimension.

The description is complicated by the presence of an option to transpose A before multiplying into B, which requires that A or A^T should conform with B for multiplication, depending on the value provided for trans.

Among the MKL examples you will find (you may have to unzip the examples zip file if you have not done so already) mkl_?csrmult.f, where ? stands for s, d, c, z. Look at the example and see if that helps.

[forrest_R_]

mecej4 wrote:

I agree that 'ldc' is incorrectly described as an output argument to mkl_?csr_multd(). Rather, the matrix C is the output argument, and ldC is its leading dimension.

The description is complicated by the presence of an option to transpose A before multiplying into B, which requires that A or A^T should conform with B for multiplication, depending on the value provided for trans.

Among the MKL examples you will find (you may have to unzip the examples zip file if you have not done so already) mkl_?csrmult.f, where ? stands for s, d, c, z. Look at the example and see if that helps.

Thank you for your help. I have seen the example. Since in the example the matrix is a square matrix, the m,n,k are the same number,  there will be no mismatch. 

I think the ia is always an array of length m+1. And ib is k+1, no matter when the trans is 'N' or 'T'. Is it right?

[mecej4]

Here is an example with rectangular matrices. There are two cases, one with trans = 'N' and the other with 'T'. You and I know that At = A^T, but it is not easy to detect that in the CSR representation. Therefore, seeing the same output for C with the different values for trans is a good test that the calls are being made properly.

program xscsrmultd
use mkl_mod
!
! illustrate usage of mkl_scsrmultd with trans = 'N' and 'Y'
! A is 4 X 3, B is 3 X 2, so C is 4 X 2.
! At is the transpose of A, stored in CSR format
!
implicit none
integer :: i,j
integer :: mA=4, nA=3, mB=3, nB=2, mC=4, nC=2
integer :: iA(5)  = [1, 3, 4, 6, 8],         jA(7) = [1,2, 2, 1,3, 1,2]
integer :: iB(4)  = [1, 2, 4, 5],            jB(4) = [2, 1,2, 2]
real    :: A(7)   = [4, -2, 4, 2, 3, -1, -1], B(4) = [2, 2, -1, 3], C(4,2)
integer :: iAt(4) = [1,4,7,8], jAt(7) = [1,3,4, 1,2,4, 3]
real    :: At(7)  = [4,2,-1,-2,4,-1,3]
! C = A.B
call mkl_scsrmultd('N',mA,nA,nB, A,   ja, ia,  B, jB, iB, C, mC)
write(*,10)((C(i,j),j=1,nC),i=1,mC)
10 format(/(2F10.4))
! C = (At)^T.B
call mkl_scsrmultd('T',nA,mA,nB, At, jAt, iAt, B, jB, iB, C, mC)
write(*,10)((C(i,j),j=1,nC),i=1,mC)
end program

Here is the source code of mkl_mod, which I use to avoid recompiling mkl.fi everytime it is used:

      module mkl_mod
      include 'mkl.fi'
      end module

(I do have six leading blanks on each of these three lines, but the forum software removes them from the display).