Community
cancel
Showing results for
Did you mean: Beginner
144 Views

## Problem with sparse BLAS Level 2 function mkl_scsrgemv()

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.
1 Solution New Contributor I
144 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

8 Replies New Contributor I
144 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 New Contributor I
145 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 Black Belt
144 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. Beginner
144 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. Employee
144 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 New Contributor I
144 Views

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

Victor Beginner
144 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.  