inconsistent results from mkl_dcscmv and mkl_dcsrmv

Kuo-ling_Huang · ‎06-30-2014

Hello,

I tried to calculate the (sparse) matrix-vector product using mkl_dcscmv and mkl_dcsrmv. However, sometimes they gave me different result. Here is the example I was using:

Let

M : = [ 1 1 0 0 ], x := [ 1 ], and sols := [ 1 ]
[ 1 0 1 0 ] [ 1 ] [ 1 ]
[ 1 0 -1 0 ] [ 1 ] [ 1 ]
[ 0 1 0 0 ] [ 1 ] [ 1 ]
[ 0 0 1 0 ] [ 1 ].

I tried to calculate

test 1: sols := M x - sols,
test 2: sols := Mx, and
test 3: sols := Mx + sols

by using either mkl_dcscmv or mkl_dcsrmv. Here are the generated result:

=== Test 1: sols := M x - sols ===
(Mx - sols) using mkl_dcscmv (expected: [1 1 -1 0 0]^T)
sols[0] = 1.000000E+00.
sols[1] = 1.000000E+00.
sols[2] = -1.000000E+00.
sols[3] = 0.000000E+00.
sols[4] = 2.000000E+00.
(Mx - sols) using mkl_dcsrmv (expected: [1 1 -1 0 0]^T)
sols[0] = 1.000000E+00.
sols[1] = 1.000000E+00.
sols[2] = -1.000000E+00.
sols[3] = 0.000000E+00.
sols[4] = 0.000000E+00.

=== Test 2: sols := M x + 0 * sols ===
result 2.1: (Mx) using mkl_dcscmv (expected: [2 2 0 1 1]^T)
sols[0] = 2.000000E+00.
sols[1] = 2.000000E+00.
sols[2] = 0.000000E+00.
sols[3] = 1.000000E+00.
sols[4] = 2.000000E+00.
result 2.2: (Mx) using mkl_dcsrmv (expected: [2 2 0 1 1]^T)
sols[0] = 2.000000E+00.
sols[1] = 2.000000E+00.
sols[2] = 0.000000E+00.
sols[3] = 1.000000E+00.
sols[4] = 1.000000E+00.

=== Test 3: sols := M x + 1 * sols ===
result 3.1: (Mx + sol) using mkl_dcscmv (expected: [3 3 1 2 2]^T)
sols[0] = 3.000000E+00.
sols[1] = 3.000000E+00.
sols[2] = 1.000000E+00.
sols[3] = 2.000000E+00.
sols[4] = 2.000000E+00.
result: 3.2: (Mx + sol) using mkl_dcsrmv (expected: [3 3 1 2 2]^T)
sols[0] = 3.000000E+00.
sols[1] = 3.000000E+00.
sols[2] = 1.000000E+00.
sols[3] = 2.000000E+00.
sols[4] = 2.000000E+00.

===============================

My computational environment:

OS: Red Hat Enterprise Linux Server release 6.4 (Santiago)

Compiler and linker: Intel Composer-XE version 2013.2.146

GNU libc version: 2.12

===============================

Following is the complete code I was using for the test:

//===================

#include <stdlib.h>

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

int main (void)
{
#define NUM_ROWS 5
#define NUM_COLS 4
#define NUM_ENTS 8
    
    MKL_INT num_rows = NUM_ROWS;
    MKL_INT num_cols = NUM_COLS;
    MKL_INT num_ents = NUM_ENTS;
    MKL_INT i;
    
    double minus_one = -1.0;
    double one = +1.0;
    double zero = 0.0;
    char   notran = 'N';
    char   matdescra[4] = {'G', 'L', 'N', 'C'}; 
    
    /* M : = [ 1  1  0  0 ] x := [ 1 ] sols := [ 1 ]  
     *       [ 1  0  1  0 ]      [ 1 ]         [ 1 ]
     *       [ 1  0 -1  0 ]      [ 1 ]         [ 1 ]
     *       [ 0  1  0  0 ]      [ 1 ]         [ 1 ]
     *       [ 0  0  1  0 ]                    [ 1 ]
     * 
     */
    MKL_INT    M_bgn[NUM_COLS + 1] = {0, 3, 5, 8, 8};
    MKL_INT    M_idx[NUM_ENTS]     = {0, 1, 2, 0, 3, 1,  2, 4};
    double     M_val[NUM_ENTS]     = {1, 1, 1, 1, 1, 1, -1, 1};
   
    MKL_INT    MT_bgn[NUM_ROWS + 1] = {0, 2, 4, 6, 7, 8};
    MKL_INT    MT_idx[NUM_ENTS]     = {0, 1, 0, 0, 0, 2,  1, 2};
    double     MT_val[NUM_ENTS]     = {1, 1, 1, 1, 1, -1, 1, 1};
    
    double     sols[NUM_ROWS];
    double     x[NUM_COLS];
    
    /* initialize the solution and x */
#define INIT_SOL                                \
    for (i = 0; i < num_rows; ++i)              \
    {                                           \
        sols = 1.0;                          \
    }                                           \
    for (i = 0; i < num_cols; ++i)              \
    {                                           \
        x = 1.0;                             \
    }
    
#define PRINT_SOL(MSG)                          \
    printf("%s\n", MSG);                        \
    for (i = 0; i < num_rows; ++i)              \
    {                                           \
        printf("sols[%d] = %E.\n", i, sols); \
    }    
    
    /* test 1: compute sols := M x - sols
     *                       = [2 2 0 1 1]^T - [1 1 1 1 1]^T
     *                       = [1 1 -1 0 0]^T.
     */
    printf("\n=== Test 1: sols := M x - sols ===\n");
    INIT_SOL;
    /* test 1.1: using mkl_dcscmv */
    mkl_dcscmv(&notran, &num_rows, &num_cols, &one, matdescra,
        M_val, M_idx, M_bgn, M_bgn + 1, x, &minus_one, sols);     
    PRINT_SOL("(Mx - sols) using mkl_dcscmv (expected: [1 1 -1 0 0]^T)");
 
    INIT_SOL;    
    /* test 1.2: using mkl_dcsrmv */
    mkl_dcsrmv(&notran, &num_rows, &num_cols, &one, matdescra,
        MT_val, MT_idx, MT_bgn, MT_bgn + 1, x, &minus_one, sols);  
    PRINT_SOL("(Mx - sols) using mkl_dcsrmv (expected: [1 1 -1 0 0]^T)");
     
    /* test 2: compute sols := M x + 0 * sols
     *                       = [2 2 0 1 1]^T - 0 * [1 1 1 1 1]^T
     *                       = [2 2 0 1 1]^T.
     */
    printf("\n=== Test 2: sols := M x + 0 * sols ===\n");
    INIT_SOL;    
    /* test 2.1: using mkl_dcscmv */
    mkl_dcscmv(&notran, &num_rows, &num_cols, &one, matdescra, 
        M_val, M_idx, M_bgn, M_bgn + 1, x, &zero, sols);     
    PRINT_SOL("result 2.1: (Mx) using mkl_dcscmv (expected: [2 2 0 1 1]^T)");
    
    INIT_SOL;  
    /* test 2.2: using mkl_dcsrmv */
    mkl_dcsrmv(&notran, &num_rows, &num_cols, &one, matdescra,
        MT_val, MT_idx, MT_bgn, MT_bgn + 1, x, &zero, sols);  
    PRINT_SOL("result 2.2: (Mx) using mkl_dcsrmv (expected: [2 2 0 1 1]^T)");    
         
    /* test 3: compute sols := M x + 1 * sols
     *                       = [2 2 0 1 1]^T + 1 * [1 1 1 1 1]^T
     *                       = [3 3 1 2 2]^T.
     */ 
    printf("\n=== Test 3: sols := M x + 1 * sols ===\n");
    INIT_SOL;    
    /* test 3.1: using mkl_dcscmv */
    mkl_dcscmv(&notran, &num_rows, &num_cols, &one, matdescra,
        M_val, M_idx, M_bgn, M_bgn + 1, x, &one, sols);     
    PRINT_SOL("result 3.1: (Mx + sol) using mkl_dcscmv (expected: [3 3 1 2 2]^T)");
  
    INIT_SOL;  
    /* test 3.2: using mkl_dcsrmv */    
    mkl_dcsrmv(&notran, &num_rows, &num_cols, &one, matdescra,
        MT_val, MT_idx, MT_bgn, MT_bgn + 1, x, &one, sols);  
    PRINT_SOL("result: 3.2: (Mx + sol) using mkl_dcsrmv (expected: [3 3 1 2 2]^T)");    
 
    return 0;
}

//===================

Ying_H_Intel · ‎06-30-2014

Hello,

I have only Composer XE 2013.5 (MKL 11.0.5) installed, just try the code, it seems run fine. It may be a bug in composer XE 2013.2 ( MKL 11.0.2) . See https://software.intel.com/en-us/articles/intel-mkl-110-bug-fixes

Please upgrade the MKL version if possible ( the latest one : MKL 11.1. 3 ).

Best Regards,

Ying

$source /opt/intel/composer_xe_2013.5.192/mkl/bin/mklvars.sh intel64
$ gcc dcsrmv.c -lmkl_intel_lp64 -lmkl_intel_thread -lmkl_core -liomp5 -lm -lpthread

$ ./a.out

=== Test 1: sols := M x - sols ===
(Mx - sols) using mkl_dcscmv (expected: [1 1 -1 0 0]^T)
sols[0] = 1.000000E+00.
sols[1] = 1.000000E+00.
sols[2] = -1.000000E+00.
sols[3] = 0.000000E+00.
sols[4] = 0.000000E+00.
(Mx - sols) using mkl_dcsrmv (expected: [1 1 -1 0 0]^T)
sols[0] = 1.000000E+00.
sols[1] = 1.000000E+00.
sols[2] = -1.000000E+00.
sols[3] = 0.000000E+00.
sols[4] = 0.000000E+00.

=== Test 2: sols := M x + 0 * sols ===
result 2.1: (Mx) using mkl_dcscmv (expected: [2 2 0 1 1]^T)
sols[0] = 2.000000E+00.
sols[1] = 2.000000E+00.
sols[2] = 0.000000E+00.
sols[3] = 1.000000E+00.
sols[4] = 1.000000E+00.
result 2.2: (Mx) using mkl_dcsrmv (expected: [2 2 0 1 1]^T)
sols[0] = 2.000000E+00.
sols[1] = 2.000000E+00.
sols[2] = 0.000000E+00.
sols[3] = 1.000000E+00.
sols[4] = 1.000000E+00.

=== Test 3: sols := M x + 1 * sols ===
result 3.1: (Mx + sol) using mkl_dcscmv (expected: [3 3 1 2 2]^T)
sols[0] = 3.000000E+00.
sols[1] = 3.000000E+00.
sols[2] = 1.000000E+00.
sols[3] = 2.000000E+00.
sols[4] = 2.000000E+00.
result: 3.2: (Mx + sol) using mkl_dcsrmv (expected: [3 3 1 2 2]^T)
sols[0] = 3.000000E+00.
sols[1] = 3.000000E+00.
sols[2] = 1.000000E+00.
sols[3] = 2.000000E+00.
sols[4] = 2.000000E+00.