topic CBLAS level 3 routine problem in Intel® oneAPI Math Kernel Library

CBLAS level 3 routine problem

stanleyimko — Sun, 20 Nov 2011 07:33:00 GMT

Hi, all. I would like to use the BLAS routine dgemm to do matrix multiplication. I made a testing program and compared the results of CBLAS and BLAS. I am confused that CBLAS produced a wrong result. I check from time to time but still don't know where goes wrong. Can anyone give me some suggestion? The code and results are pasted below. Thanks in advance.

int main()
{
const int n = 3, m = 2;
double A, B;
double d[n*m], e[n*m];

/*************************
two dimension storage:
A = 1, 2 B = 1, 1
2, 3 1, 1
3, 4 1, 1
*************************/

for( int i = 0; i < n; i++ )
for( int j = 0; j < m; j++ ){
A = i+j+1;
B = 1;
}

/*************************
one dimension storage
with column-major:
d = {1, 2, 3, 2, 3, 4}
e = {1, 1, 1, 1, 1, 1}
*************************/

for( int j = 0; j < m; j++ )
for( int i = 0; i < n; i++ ){
d[j*n+i] = A;
e[j*n+i] = B;
}

/*********************************
Call dgemm:
C := alpha*op(A)*op(B) + beta*C
********************************/

/************CBLAS***************/

CBLAS_ORDER order;
order = CblasRowMajor; //specify row major
CBLAS_TRANSPOSE transa;
transa = CblasTrans; // transpose A
CBLAS_TRANSPOSE transb;
transb = CblasNoTrans;

int lda = n; //the first dimension of A
int ldb = n; //the first dimension of B
int ldc = m; //the first dimension of C

double alpha = 1.0, beta = 0.0;
double C;

cblas_dgemm( order, transa, transb, m, m, n,
alpha, *A, lda, *B, ldb, beta, *C, ldc);

for( int i = 0; i < m; i++ ){
for( int j = 0; j < m; j++ ){
cout << C << ", ";
}
cout << endl;
}

/*********direct Fotran call*********/

double f[m*m];
char transaa = 'T';
char transbb = 'N';
dgemm( &transaa, &transbb, &m, &m, &n,
α, d, &lda, e, &ldb, β, f, &ldc);

cout << endl;
for( int i = 0; i < m*m; i++ )
cout << f << ", ";
cout << endl;

return 0;
}

Result:

5, 6,
6, 7,

6, 9, 6, 9,

CBLAS level 3 routine problem

Murat_G_Intel — Mon, 21 Nov 2011 05:45:16 GMT

Hello,

The leading dimensions for the CBLAS call are set incorrectly. For the row-major matrix format, the leading dimensions should stride rows of the matrix. Since your 2D arrays are A[m] and B[m], lda and ldb should be m. Please see themodified code below:

[cpp]    CBLAS_ORDER order;
    order = CblasRowMajor; //specify row major
    CBLAS_TRANSPOSE transa;
    transa = CblasTrans; // transpose A
    CBLAS_TRANSPOSE transb;
    transb = CblasNoTrans;

    int lda = m; //modified leading dimension for A
    int ldb = m; //modified leading dimension for B
    int ldc = m; //the first dimension of C
[/cpp]

Then, before the FORTRAN dgemm call, we need to reset lda and ldb to n since the FORTRAN call is in column-major format.

[cpp]    /*********direct Fotran call*********/

    lda = n; //reset the leading dimension of A
    ldb = n; //reset the leading dimension of B
    double f[m*m];
    char transaa = 'T';
    char transbb = 'N';
    dgemm( &transaa, &transbb, &m, &m, &n,
            α, d, &lda, e, &ldb, β, f, &ldc);[/cpp]

After these modifications, I get same results from CBLAS and FORTRAN dgemm:

Result:

6, 6,

9, 9,

6, 9, 6, 9,

Best wishes,

Efe

CBLAS level 3 routine problem

stanleyimko — Mon, 21 Nov 2011 12:53:36 GMT

Efe,

Thanks so much. I am kind of misunderstanding the dimension stuff. I should check carefully before posting the thread next time.

Regards,

Stanley