Hi Tianxiong,

Tianxiong_Lu · ‎09-14-2015

Hi all,

I have two questions about memory management of sparse_matrix_t :

1. When I create a matrix handle using routine mkl_sparse_?_create_csr, the handle will copy all data to new memory blocks? Or it only save the pointers in handle? If it copy all data to new block, and the routine mkl_sparse_destroy will free these memory blocks or not ?

2. As the result of routine mkl_sparse_spmm , is it created in mkl routines？I will call destroy routine for it ? And then the internal data like column indexes array and matrix values array should be free expicit or not?

By the way, I found routines declared in header file : mkl_sparse_?_set_value, but they are not found in manual. And there are comment about these routines:

/* update existing value in the matrix ( for internal storage only, should not work with user-allocated matrices) */

It means I can not using them in my code?

Thanks,

Tianxiong Lu

Alexander_K_Intel2 · ‎09-14-2015

Hi,

Please see my comments below

Thanks,

Alex

> When I create a matrix handle using routine mkl_sparse_?_create_csr, the handle will copy all data to new memory blocks? Or it only save the pointers in handle?

It save pointer in handle only.

> 2. As the result of routine mkl_sparse_spmm , is it created in mkl routines？I will call destroy routine for it ? And then the internal data like column indexes array and matrix values array should be free expicit or not?

Didn't got it. This function return sparse_matrix_t *C which contain result of multiplication and could be used in further routines and/or export routine.

> By the way, I found routines declared in header file : mkl_sparse_?_set_value, but they are not found in manual. And there are comment about these routines:

/* update existing value in the matrix ( for internal storage only, should not work with user-allocated matrices) */

It means I can not using them in my code?

This function allow you to change value of matrix in sparse_matrix structure so i don't see any barrier to use it in your code

Thanks,

Tianxiong Lu

Tianxiong_Lu · ‎09-14-2015

Hi Alex,

1. What does the means "should not work with user-allocated matrices“ ？Maybe I can ignore it.

2. -- Didn't got it. This function return sparse_matrix_t *C which contain result of multiplication and could be used in further routines and/or export routine.

Did you mean that MKL will free struct C ? Or I should destroy it while program finished.

Thanks,

Tianxiong Lu

Alexander_K_Intel2 · ‎09-24-2015

Hi,

1. What does the means "should not work with user-allocated matrices“ ？Maybe I can ignore it. - yes, just ignore i

2. Did you mean that MKL will free struct C.

In the end of the program you need to call mkl_sparse_destroy ( C ) routine to free struct C

Thanks,

Alex

Ying_H_Intel · ‎05-22-2016

Hi Tianxiong,

Just update the thread.

There is new example, sparse_sppm.c added in MKL 11.3 update 3, in which we show how to use mkl_sparse_spmm and mkl_sparse_export_csr.

There also is shown how to work with attendant arrays and memory for them.

You are welcomed to try it.

Cheers,

Ying

/*******************************************************************************
* Copyright 2013-2016 Intel Corporation All Rights Reserved.
*
* The source code,  information  and material  ("Material") contained  herein is
* owned by Intel Corporation or its  suppliers or licensors,  and  title to such
* Material remains with Intel  Corporation or its  suppliers or  licensors.  The
* Material  contains  proprietary  information  of  Intel or  its suppliers  and
* licensors.  The Material is protected by  worldwide copyright  laws and treaty
* provisions.  No part  of  the  Material   may  be  used,  copied,  reproduced,
* modified, published,  uploaded, posted, transmitted,  distributed or disclosed
* in any way without Intel's prior express written permission.  No license under
* any patent,  copyright or other  intellectual property rights  in the Material
* is granted to  or  conferred  upon  you,  either   expressly,  by implication,
* inducement,  estoppel  or  otherwise.  Any  license   under such  intellectual
* property rights must be express and approved by Intel in writing.
*
* Unless otherwise agreed by Intel in writing,  you may not remove or alter this
* notice or  any  other  notice   embedded  in  Materials  by  Intel  or Intel's
* suppliers or licensors in any way.
*******************************************************************************/

/******************************************************************************
*
* Consider the matrix A
*
*                 |  10     11      0     0     0   |
*                 |   0      0     12    13     0   |
*   A    =        |  15      0      0     0    14   |,
*                 |   0     16     17     0     0   |
*                 |   0      0      0    18    19   |
*
* and diagonal matrix B
*
*                 |   5      0      0     0     0   |
*                 |   0      6      0     0     0   |
*   B    =        |   0      0      7     0     0   |.
*                 |   0      0      0     8     0   |
*                 |   0      0      0     0     9   |
*
*  Both matrices A and B are stored in a zero-based compressed sparse row (CSR) storage
*  scheme with three arrays (see 'Sparse Matrix Storage Schemes' in the
*  Intel Math Kernel Library Developer Reference) as follows:
*
*           values_A = ( 10  11  12  13  15  14  16  17  18  19 )
*          columns_A = (  0   1   2   3   0   4   1   2   3   4 )
*         rowIndex_A = (  0       2       4       6       8      10 )
*
*           values_B = ( 5  6  7  8  9  )
*          columns_B = ( 0  1  2  3  4  )
*         rowIndex_B = ( 0  1  2  3  4  5 )
*
*  The example computes two scalar products :
*
*         < (A*B)*x ,       y > = left,   using MKL_SPARSE_D_SPMM and CBLAS_DDOT.
*         <     B*x , (A^t)*y > = right,  using MKL_SPARSE_D_MV and CBLAS_DDOT.
*
*         These products should result in the same value. To obtain matrix C,
*         use MKL_SPARSE_D_EXPORT_CSR and print the result.
*
******************************************************************************/

#include <stdio.h>
#include <assert.h>
#include <math.h>
#include "mkl.h"

int main() {

#define M 5
#define NNZ 10
#define ALIGN 128

/* To avoid constantly repeating the part of code that checks inbound SparseBLAS functions' status,
   use macro CALL_AND_CHECK_STATUS */
#define CALL_AND_CHECK_STATUS(function, error_message) do { \
          if(function != SPARSE_STATUS_SUCCESS)             \
          {                                                 \
          printf(error_message); fflush(0);                 \
          status = 1;                                       \
          goto memory_free;                                 \
          }                                                 \
} while(0)

/* Declaration of values */
    double  *values_A = NULL, *values_B = NULL, *values_C = NULL;
    MKL_INT *columns_A = NULL, *columns_B = NULL, *columns_C = NULL;
    MKL_INT *rowIndex_A = NULL, *rowIndex_B = NULL, *pointerB_C = NULL, *pointerE_C = NULL;

    double  *rslt_mv = NULL, *rslt_mv_trans = NULL, *x = NULL, *y = NULL;

    double   left, right, residual;
    MKL_INT  rows, cols, i, j, ii, status;

    sparse_index_base_t    indexing;
    struct matrix_descr    descr_type_gen;
    sparse_matrix_t        csrA = NULL, csrB = NULL, csrC = NULL;

/* Allocation of memory */
    values_A = (double *)mkl_malloc(sizeof(double) * NNZ, ALIGN);
    columns_A = (MKL_INT *)mkl_malloc(sizeof(MKL_INT) * NNZ, ALIGN);
    rowIndex_A = (MKL_INT *)mkl_malloc(sizeof(MKL_INT) * (M + 1), ALIGN);

    values_B = (double *)mkl_malloc(sizeof(double) * M, ALIGN);
    columns_B = (MKL_INT *)mkl_malloc(sizeof(MKL_INT) * M, ALIGN);
    rowIndex_B = (MKL_INT *)mkl_malloc(sizeof(MKL_INT) * (M + 1), ALIGN);

    x = (double *)mkl_malloc(sizeof(double) * M, ALIGN);
    y = (double *)mkl_malloc(sizeof(double) * M, ALIGN);
    rslt_mv = (double *)mkl_malloc(sizeof(double) * M, ALIGN);
    rslt_mv_trans = (double *)mkl_malloc(sizeof(double) * M, ALIGN);

/* Set values of the variables*/
    descr_type_gen.type = SPARSE_MATRIX_TYPE_GENERAL;
    status = 0, ii = 0;
 //Matrix A 
    for( i = 0; i < NNZ; i++ )
          values_A = i + 10;
    for( i = 0; i < NNZ; i++ )
          columns_A = i % 5;
    rowIndex_A[0] = 0;
    for( i = 1; i < M + 1; i++ )
          rowIndex_A = rowIndex_A[i - 1] + 2;

 //Matrix B
    for( i = 0; i < M; i++ )
          values_B = i + 5;
    for( i = 0; i < M; i++ )
          columns_B = i % 5;
    for( i = 0; i < M + 1; i++ )
          rowIndex_B = i;

 //Vectors x and y
    for( i = 0; i < M; i++ )
    {
          x = 1.0; y = 1.0;
    }
/* Printing usable data */
    printf( "\n\n_______________Example program for MKL_SPARSE_D_SPMM_________________\n\n" );
    printf( " COMPUTE  A * B = C, where matrices are stored in CSR format\n" );
    printf( "\n MATRIX A:\nrow# : (value, column) (value, column)\n" );
    for( i = 0; i < M; i++ )
    {
        printf("row#%d:", i + 1); fflush(0);
        for( j = rowIndex_A; j < rowIndex_A[i+1]; j++ )
        {
            printf(" (%5.0f, %6d)", values_A[ii], columns_A[ii] ); fflush(0);
            ii++;
        }
        printf( "\n" );
    }
    ii = 0;
    printf( "\n MATRIX B:\nrow# : (value, column)\n" );
    for( i = 0; i < M; i++ )
    {
        printf("row#%d:", i + 1); fflush(0);
        for( j = rowIndex_B; j < rowIndex_B[i+1]; j++ )
        {
            printf(" (%5.0f, %6d)", values_B[ii], columns_B[ii] ); fflush(0);
            ii++;
        }
        printf( "\n" );
    }
    printf( "\n Check the resultant matrix C, using two scalar products\n" );
    printf( " (values of these scalar products must match).\n" );

/* Prepare arrays, which are related to matrices.
   Create handles for matrices A and B stored in CSR format */
    CALL_AND_CHECK_STATUS(mkl_sparse_d_create_csr( &csrA, SPARSE_INDEX_BASE_ZERO, M, M, rowIndex_A, rowIndex_A+1, columns_A, values_A ),
                          "Error after MKL_SPARSE_D_CREATE_CSR, csrA \n");
    CALL_AND_CHECK_STATUS(mkl_sparse_d_create_csr( &csrB, SPARSE_INDEX_BASE_ZERO, M, M, rowIndex_B, rowIndex_B+1, columns_B, values_B ),
                          "Error after MKL_SPARSE_D_CREATE_CSR, csrB \n");

/* Compute C = A * B  */
    CALL_AND_CHECK_STATUS(mkl_sparse_spmm( SPARSE_OPERATION_NON_TRANSPOSE, csrA, csrB, &csrC ),
                          "Error after MKL_SPARSE_SPMM \n");

/* Analytic Routines for MKL_SPARSE_D_MV.
   HINTS: provides estimate of number and type of upcoming matrix-vector operations
   OPTIMIZE: analyze sparse matrix; choose proper kernels and workload balancing strategy */
    CALL_AND_CHECK_STATUS(mkl_sparse_set_mv_hint( csrA, SPARSE_OPERATION_TRANSPOSE,     descr_type_gen, 1 ),
                          "Error after MKL_SPARSE_SET_MV_HINT, csrA \n");
    CALL_AND_CHECK_STATUS(mkl_sparse_set_mv_hint( csrB, SPARSE_OPERATION_NON_TRANSPOSE, descr_type_gen, 1 ),
                          "Error after MKL_SPARSE_SET_MV_HINT, csrB \n");
    CALL_AND_CHECK_STATUS(mkl_sparse_set_mv_hint( csrC, SPARSE_OPERATION_NON_TRANSPOSE, descr_type_gen, 1 ),
                          "Error after MKL_SPARSE_SET_MV_HINT, csrC \n");

    CALL_AND_CHECK_STATUS(mkl_sparse_optimize( csrA ),
                          "Error after MKL_SPARSE_OPTIMIZE, csrA \n");
    CALL_AND_CHECK_STATUS(mkl_sparse_optimize( csrB ),
                          "Error after MKL_SPARSE_OPTIMIZE, csrB \n");
    CALL_AND_CHECK_STATUS(mkl_sparse_optimize( csrC ),
                          "Error after MKL_SPARSE_OPTIMIZE, csrC \n");

/* Execution Routines */
/* Step 1:
          Need to compute the following variables:
                 rslt_mv = C * x
                    left = <rslt_mv, y>              */
    CALL_AND_CHECK_STATUS(mkl_sparse_d_mv( SPARSE_OPERATION_NON_TRANSPOSE, 1.0, csrC, descr_type_gen, x, 0.0, rslt_mv ),
                          "Error after MKL_SPARSE_D_MV, csrC*x  \n");
    left = cblas_ddot( M, rslt_mv, 1, y, 1 );

/* Step 2:
          Need to compute the following variables:
           rslt_mv       =     B * x
           rslt_mv_trans = (A)^t * y
                   right = <rslt_mv, rslt_mv_trans>  */

    CALL_AND_CHECK_STATUS(mkl_sparse_d_mv( SPARSE_OPERATION_NON_TRANSPOSE, 1.0, csrB, descr_type_gen, x, 0.0, rslt_mv ),
                          "Error after MKL_SPARSE_D_MV, csrB*x  \n");
    CALL_AND_CHECK_STATUS(mkl_sparse_d_mv( SPARSE_OPERATION_TRANSPOSE,     1.0, csrA, descr_type_gen, y, 0.0, rslt_mv_trans),
                          "Error after MKL_SPARSE_D_MV, csrA*y  \n");
    right = cblas_ddot( M, rslt_mv, 1, rslt_mv_trans, 1);

/* Step 3:
          Compare values obtained for left and right  */
    residual = fabs(left - right)/(fabs(left)+1);

    printf( "\n The difference between < C*x , y > and < B*x , (A^t)*y > = %g,\n", residual );
    printf( " which means that MKL_SPARSE_D_SPMM arrived correct at a solution.\n" );
/* Printing OUTPUT DATA */
    CALL_AND_CHECK_STATUS(mkl_sparse_d_export_csr( csrC, &indexing, &rows, &cols, &pointerB_C, &pointerE_C, &columns_C, &values_C ),
                          "Error after MKL_SPARSE_D_EXPORT_CSR  \n");

    printf( "\n RESULTANT MATRIX C:\nrow# : (value, column) (value, column)\n" );
    ii = 0;
    for( i = 0; i < M; i++ )
    {
        printf("row#%d:", i + 1); fflush(0);
        for( j = pointerB_C; j < pointerE_C; j++ )
        {
            printf(" (%5.0f, %6d)", values_C[ii], columns_C[ii] ); fflush(0);
            ii++;
        }
        printf( "\n" );
    }
    printf( "_____________________________________________________________________  \n" );

/* Deallocate memory */
memory_free:
 //Release matrix handle. Not necessary to deallocate arrays for which we don't allocate memory: values_C, columns_C, pointerB_C, and pointerE_C.
 //These arrays will be deallocated together with csrC structure.
    if( mkl_sparse_destroy( csrC ) != SPARSE_STATUS_SUCCESS)
    { printf(" Error after MKL_SPARSE_DESTROY, csrC \n");fflush(0); status = 1; }

 //Deallocate arrays for which we allocate memory ourselves.
    mkl_free(rslt_mv_trans); mkl_free(rslt_mv); mkl_free(x); mkl_free(y);

 //Release matrix handle and deallocate arrays for which we allocate memory ourselves.
    if( mkl_sparse_destroy( csrA ) != SPARSE_STATUS_SUCCESS)
    { printf(" Error after MKL_SPARSE_DESTROY, csrA \n");fflush(0); status = 1; }
    mkl_free(values_A); mkl_free(columns_A); mkl_free(rowIndex_A);

    if( mkl_sparse_destroy( csrB ) != SPARSE_STATUS_SUCCESS)
    { printf(" Error after MKL_SPARSE_DESTROY, csrB \n");fflush(0); status = 1; }
    mkl_free(values_B); mkl_free(columns_B); mkl_free(rowIndex_B);


    return status;
}

What is the memory management strategy in sparse_matrix_t ?