Intel® oneAPI Math Kernel Library
Ask questions and share information with other developers who use Intel® Math Kernel Library.
6956 Discussions

Problem with sparse BLAS Level 2 function mkl_scsrgemv()

Wu__Wentao
Beginner
610 Views
Hi!
When programming with mkl sparse BLAS level 2 fucntion to implement a matrix-vector mulitplication, I met a trouble in function mlk_?csrgemv(). A simplified and clear code is given below, in which the same problem happened.
!-----------------------------------------------
program main
implicit none
real, dimension(12) :: val = (/-1, 1,-1, &
1,-1, 1, &
1,-1, 1, &
-1, 1,-1 /)
integer, dimension(12) :: column = (/1,2,3, 4,5,6, 7,8,9, 1,4,7/)
integer, dimension(5) :: rowIndex= (/1,4,7, 10,13/)
real, dimension(9) :: v = 1.0
real, dimension(9) :: res = 0
character :: transa = 'n'
integer :: row_num = 9
call mkl_scsrgemv(transa, row_num, val, rowIndex, column, v, res)
write(*,*) "result = "
write(*,*) res
stop
end program
!------------------------------------------
Error info is: Unhandled exception at 0x0117d008 in SparseBLAS.exe: 0xC0000005: Access violation reading location 0xff25e07c.
It shoud be pointed out that the matrix is not a square matrix of the dimension m*m.
But another non square sparse matrix,
val = (/-1,1,1,-1,1,-1,-1,1/)
column = (/1,2,3,4,5,6,2,4/)
rowIndex = (/1,3,5,7,9/)
worked successfully.
My question is, what is wrong with the code? Can the nonsquare matrix vector multiplication be implemented with CSR(3 array variation) data structure? And How?
(Platform: Intel Pentium T4200, VS2008, MKL 10.3)
Expecting for your help and Thanks a lot.
0 Kudos
1 Solution
Victor_Gladkikh
New Contributor I
610 Views

Or you can keep row_num = 9 and extend RowIndex: integer, dimension(10) ::rowindex=(/1,4,7,10,13,13,13,13,13,13/)

Victor

View solution in original post

0 Kudos
8 Replies
Victor_Gladkikh
New Contributor I
610 Views
Hi David

As I understand your first matrix has 4 rows and 9 columns with next scheme:
(x x x - - - - - -)
(- - - x x x- - -)
(- - -- - - x x x)
(x- - x - - x - -)


In this cases youshould use row_num =4 instead of row_num = 9.

Victor
0 Kudos
Victor_Gladkikh
New Contributor I
611 Views

Or you can keep row_num = 9 and extend RowIndex: integer, dimension(10) ::rowindex=(/1,4,7,10,13,13,13,13,13,13/)

Victor

0 Kudos
mecej4
Honored Contributor III
610 Views
What is the size of your matrix? Could you show the full matrix in standard mathematical notation? Then we could check whether you are describing it in compact row-major storage format.

See also the examples that came with MKL or the ones in the MKL documentation.
0 Kudos
Wu__Wentao
Beginner
610 Views
Victor, Thank you so much!

This is the right answer I am expecting. I have tried and it worked successfully.
But setting the row number to be 4 seems not a perfect idea, though it returns right answer if transa = 'N'. The result is incorrect without any warning and error reports from compiler when transa = 'T'. However, the extended rowIndex works no matter transa = 'N' or 'T'.
Best regards,
David.
0 Kudos
Sergey_P_Intel2
Employee
610 Views

Hi David,

In MKL manual this functionality is described as working with square matrices:

The mkl_?csrgemv routine performs a matrix-vector operation defined as

y := A*x

or

y := A'*x,

where: x and y are vectors, A is an m-by-m sparse square matrix in the CSR format (3-array variation), A' is the transpose of A.

So for rectangular matrices it works improperly in some cases.

Regards,

Sergey

0 Kudos
Victor_Gladkikh
New Contributor I
610 Views

Anyway the second solution have to work because we extend rectangular matrixto square matrix.

Victor

0 Kudos
Wu__Wentao
Beginner
610 Views

Anyway the second solution have to work because we extend rectangular matrixto square matrix.

Victor

Victor, Thank you very much again.
While, When you said we extend the rectanglar matrix to a square one, I get a new idea:
If we have a exact matrix,
(/ x x x x
x 1 2 x
x x x x
3 x 4 x
x x x x
x 5 6 x /)
then it can be presented in a sparse form with CSR (3 array variation)
value = (/1, 2, 3, 4, 5, 6/)
column = (/2, 3, 1, 3, 2, 3/)
rowIndex = (/1, 1, 3, 3, 5, 5, 7/).
That is to say, if there is a trival row(all zero elements row) in the matrix, a duplicated row index can be used to represent this row. I have tried some examples, and it does work!
But according to the MKL reference manual, "The 3-array variation of the CSR format has a restriction: all non-zero elements are stored continuously, that is the set of non-zero elements in the rowJgoes just after the set of non-zero elements in the rowJ-1."
What do you think of it, does my duplicated row index method always work? Or it just fits for the lucky simple cases.
Regrads,
David.
0 Kudos
mecej4
Honored Contributor III
610 Views
Your proposal, it appears to me, goes against the principal reason for the existence of sparse storage schemes: to avoid calculations with zero elements.

Once you allow more than the minimal number of zeros (this number is non-zero for banded matrices, zero in packed storage schemes), the efficiency of the sparse matrix method deteriorates.
0 Kudos
Reply