matdescra(1) = 'D' in mkl_dcsrmv doesn't work?

[Bob_G_1]

Hi,

I'm trying to write the equivalent of this matlab code in F#, with ILP64 MKL:

S = temp + (I - diag(diag(temp)))

where S, temp are matrices in Rn, I is the identity matrix in Rn, and diag(diag(temp)) is the matrix containing only the main diagonal of temp

In other words, given the matrix

[ a b c
d e f
g h i ]

I want the result

[1 b c
d 1 f
g h 1 ]

for any square sparse matrix.

I think I can perform part of this operation using mkl_dcsrmv, which computes

y := alpha*A*x + beta*y

where x and y are the values array of the identity matrix, alpha is -1.0 and beta is 1.0. According to the manual, I can operate on the diagonal by setting the first element of matdescra to 'D'.

Given the matrix

[1 2 3
4 5 6
7 8 9 ]

my calculation correctly returns the vector [-5 -14 -3] for a General matrix, but when I set matdescra(1) to D, I get [1 1 1] where I expect [0 -4 -8]

It's as though setting matdescra(1) to 'D' returns, not the main diagonal of my matrix, but [0 0 0].

I include some F# code that reproduces the issue, error handling and memory cleanup elided.

let test_dcsrmv (structure:char) =
   let p x = PinnedArray.of_array x

   let mutable transOption = 'N'
   let mutable m = 3
   let mutable k = 3
   let mutable alpha = -1.0
   let mutable beta = 1.0
   let madescra = [| structure; '_'; 'N'; 'F'|] |> PinnedArray.of_array

   let indx = p [|1; 2; 3;
                     1; 2; 3;
                     1; 2; 3; |]

   let vals = p [|
                        1.0; 2.0; 3.0;
                        4.0; 5.0; 6.0;
                        7.0; 8.0; 9.0; |]

   let pntrb = p [|1; 4; 7;|]
   let pntre = p [|4; 7; 10|]

   let x = p [| 1.0; 1.0; 1.0 |]
   let y = [| 1.0; 1.0; 1.0 |]
   let y_handle = p y
   mkl_dcsrmv(&&transOption, &&m, &&k, &&alpha, madescra.Ptr, vals.Ptr, indx.Ptr, pntrb.Ptr, pntre.Ptr, x.Ptr, &&beta, y_handle.Ptr)

   y

let general = test_dcsrmv 'G' // returns [-5 -14 -23]
let diagonal = test_dcsrmv 'D' // returns [1 1 1]

[Noah_C_Intel]

Hi, let me give this some though and I will get back to you.

[Alexander_K_Intel2]

Hi,

It's a really interesting question how using MKL functionality without direct work with arrays change a diagonal of matrix... Nevertheless your approach a bit incorrect, you try to use

y := alpha*A*x + beta*y

where x and y are the values array of the identity matrix, alpha is -1.0 and beta is 1.0

but it is equal to Y=1.*A-1*I = 

[1 2 3 ]     [1 0 0 ]     [0 2 3 ]     [1 2 3 ]
[4 5 6 ]   - [0 1 0 ] =  [4 4 6 ] != [4 1 6 ]
[7 8 9 ]     [0 0 1 ]     [7 8 8 ]     [7 8 1 ]

The parameter 'D' & 'G' in your program specify characteristic of matrix "vals" as diagonal matrix or general - in you example 'G' is more correct. 

Currently i don't see the way how to change one diagonal of matrix using SparseBals functionality only but I will investigate it. Is it possible to provide the reason of such request? Is it possible to provide additional information about it - based on it we could provide some suitable workaround of implementation such algorithm. If it is a kind of secret you can write to me in private post.

With best regards,

Alexander Kalinkin

[Bob_G_1]

Hi Alexander, thanks for the reply. I'm a little confused, because this forum post [http://software.intel.com/en-us/forums/topic/302746] seems to claim that mkl_dcsrmv WILL operate over the main diagonal if the matdescra parameter is set up correctly. Moreover, the reference manual says this will work, too. http://software.intel.com/sites/products/documentation/hpc/mkl/mklman/GUID-34C8DB79-0139-46E0-8B53-99F3BEE7B2D4.htm """ One of the distinctive feature of the Intel MKL Sparse BLAS routines is a possibility to perform operations only on partial matrices composed of certain parts (triangles and the main diagonal) of the input sparse matrix. It can be done by setting properly first three elements of the parameter matdescra. An arbitrary sparse matrix A can be decomposed as A = L + D + U where L is the strict lower triangle of A, U is the strict upper triangle of A, D is the main diagonal. """ The reason I want to do this is that I'm building an iterative implementation of the SimRank algorithm [http://ilpubs.stanford.edu:8090/508/1/2001-41.pdf]. The iterative form requires that the elements of the main diagonal are reset to 1 after each iteration, since every item is exactly similar to itself. I do have a paper describing a non-iterative implementation, which uses the kronecker product, but I'm not sure how I'd calculate that in MKL. I think the operation does work: x and y are vectors, containing a one for each row in the input matrix, we use the values array of an identity matrix because I want the matrix Y which has the vector y as its main diagonal. This technique only requires that matdescra works as described, and can operate on the diagonal of a sparse matrix. Given alpha = -1, beta = 1, x = [1 1 1], y = [1 1 1], A = [1 2 3 4 5 6 7 8 9 ] alpha * Diag(A) = [ -1 0 0 0 -5 0 0 0 -9] * [1 1 1] = [-1 -5 -9] + [1 1 1] = [0 -4 -8] A + [0 0 0 0 -4 0 0 0 -8 ] = [ 1 2 3 4 1 6 7 8 1 ]