Intel® oneAPI Math Kernel Library
Ask questions and share information with other developers who use Intel® Math Kernel Library.

Matrix Multiplication problem: Wrong results in a 1st row of a resulting matrix C ( C[8][8] = A[8][8] x B[8][8] )

SergeyKostrov
Valued Contributor II
403 Views

I detected some problem with SGEMM and DGEMM functions. In essence, there are wrong numbers in a 1st row of a resulting matrix C ( C[8][8] = A[8][8] x B[8][8] ).

It is verified with 32-bit and 64-bit versions of MKL:
...
#define __INTEL_MKL__ 10
#define __INTEL_MKL_MINOR__ 3
#define __INTEL_MKL_UPDATE__ 12
...
and
...
#define __INTEL_MKL__ 11
#define __INTEL_MKL_MINOR__ 0
#define __INTEL_MKL_UPDATE__ 2
...
using Intel, Microsoft, MinGW and Borland C++ compilers.

 

0 Kudos
13 Replies
SergeyKostrov
Valued Contributor II
403 Views
Here are matricies: A[8][8]= 5 2 3 4 5 6 7 8 1 2 3 4 5 6 7 8 1 2 3 4 5 6 7 8 1 2 3 4 5 6 7 8 1 2 3 4 5 6 7 8 1 2 3 4 5 6 7 8 1 2 3 4 5 6 7 8 1 2 3 4 5 6 7 8 x B[8][8]= 7 2 3 4 5 6 7 8 1 2 3 4 5 6 7 8 1 2 3 4 5 6 7 8 1 2 3 4 5 6 7 8 1 2 3 4 5 6 7 8 1 2 3 4 5 6 7 8 1 2 3 4 5 6 7 8 1 2 3 4 5 6 7 8 Here is a result calculated usimg an online matrix-calculator: [ www.bluebit.gr/matrix-calculator/matrix_multiplication.aspx ] C[8][8]= 70.0 80.0 120.0 160.0 200.0 240.0 280.0 320.0 42.0 72.0 108.0 144.0 180.0 216.0 252.0 288.0 42.0 72.0 108.0 144.0 180.0 216.0 252.0 288.0 42.0 72.0 108.0 144.0 180.0 216.0 252.0 288.0 42.0 72.0 108.0 144.0 180.0 216.0 252.0 288.0 42.0 72.0 108.0 144.0 180.0 216.0 252.0 288.0 42.0 72.0 108.0 144.0 180.0 216.0 252.0 288.0 42.0 72.0 108.0 144.0 180.0 216.0 252.0 288.0 Here is a result calculated by: [ MKL SGEMM and DGEMM functions ] C[8][8]= 70.0 84.0 126.0 168.0 210.0 252.0 294.0 336.0 40.0 72.0 108.0 144.0 180.0 216.0 252.0 288.0 40.0 72.0 108.0 144.0 180.0 216.0 252.0 288.0 40.0 72.0 108.0 144.0 180.0 216.0 252.0 288.0 40.0 72.0 108.0 144.0 180.0 216.0 252.0 288.0 40.0 72.0 108.0 144.0 180.0 216.0 252.0 288.0 40.0 72.0 108.0 144.0 180.0 216.0 252.0 288.0 40.0 72.0 108.0 144.0 180.0 216.0 252.0 288.0 The C/C++ test case is based on a well known MKL example for these functions ( nothing is special / alpha=1 and beta=0 ). Please do your own verifications. Thanks in advance.
0 Kudos
SergeyKostrov
Valued Contributor II
403 Views
Also, if I use 'T' and 'T' for first two input parameters in both functions then a resulting transposed matrix C is correct. In overall, wrong numbers look like some rounding error ( I could be wrong here ) and it needs to be investigated.
0 Kudos
mecej4
Honored Contributor III
403 Views

You have, apparently, overlooked the fact that C stores matrices row-by-row, whereas Fortran stores them column-by-column. You did not show any code, but your results are consistent with passing C-native 2-D arrays to routines (such as SGEMM/DGEMM) that expect Fortran-native 2-D arrays.

There are routines in MKL with the prefix cblas_  that you may consider for using in combination with C-native 2-D arrays.

0 Kudos
SergeyKostrov
Valued Contributor II
403 Views
Mecej4, Just run the test with values for matricies I gave. If I don't put values 5 and 7 repectively ( 1st values ) it calculates the C properly.
0 Kudos
SergeyKostrov
Valued Contributor II
403 Views
>>...You did not show any code... Just run the test with values for matricies I gave. If I don't put values 5 and 7 repectively ( 1st values ) it calculates the C properly.
0 Kudos
SergeyKostrov
Valued Contributor II
403 Views
A message to a Moderator: Many of my posts are delayed and please take a look. Thank you in advance.
0 Kudos
SergeyKostrov
Valued Contributor II
403 Views
Note: Corrected a 1st value in the 2nd B matrix. Just run the test with: A[8][8]= 5 2 3 4 5 6 7 8 1 2 3 4 5 6 7 8 1 2 3 4 5 6 7 8 1 2 3 4 5 6 7 8 1 2 3 4 5 6 7 8 1 2 3 4 5 6 7 8 1 2 3 4 5 6 7 8 1 2 3 4 5 6 7 8 x B[8][8]= 7 2 3 4 5 6 7 8 1 2 3 4 5 6 7 8 1 2 3 4 5 6 7 8 1 2 3 4 5 6 7 8 1 2 3 4 5 6 7 8 1 2 3 4 5 6 7 8 1 2 3 4 5 6 7 8 1 2 3 4 5 6 7 8 and after that with: A[8][8]= 1 2 3 4 5 6 7 8 1 2 3 4 5 6 7 8 1 2 3 4 5 6 7 8 1 2 3 4 5 6 7 8 1 2 3 4 5 6 7 8 1 2 3 4 5 6 7 8 1 2 3 4 5 6 7 8 1 2 3 4 5 6 7 8 x B[8][8]= 1 2 3 4 5 6 7 8 1 2 3 4 5 6 7 8 1 2 3 4 5 6 7 8 1 2 3 4 5 6 7 8 1 2 3 4 5 6 7 8 1 2 3 4 5 6 7 8 1 2 3 4 5 6 7 8 1 2 3 4 5 6 7 8
0 Kudos
SergeyKostrov
Valued Contributor II
403 Views
>>...I detected some problem with SGEMM and DGEMM functions. In essence, there are wrong numbers in a 1st row of >>a resulting matrix C ( C[8][8] = A[8][8] x B[8][8] )... A manual verification was also completed for all values of the 1st row ( here are calculations for 1st and 8th values ): ( 5*7 ) + ( 2*1 ) + ( 3*1 ) + ( 4*1 ) + ( 5*1 ) + ( 6*1 ) + ( 7*1 ) + ( 8*1 ) = 70 ... ( 5*8 ) + ( 2*8 ) + ( 3*8 ) + ( 4*8 ) + ( 5*8 ) + ( 6*8 ) + ( 7*8 ) + ( 8*8 ) = 320 Note: Not 336 and I expect that Intel software developer needs to verify that simple test case with MKL functions for matrix multiplication in its own test application ( that is, independent ).
0 Kudos
SergeyKostrov
Valued Contributor II
403 Views
>>...I detected some problem with SGEMM and DGEMM functions. In essence, there are wrong numbers in >>a 1st row of a resulting matrix C ( C[8][8] = A[8][8] x B[8][8] )... I wonder if Intel software engineers had a chance to verify ( independently ) these two functions with 8x8 matricies and data I've provided? Thanks in advance. PS: I'm very sorry but my test case can not be provided.
0 Kudos
Murat_G_Intel
Employee
403 Views

Hi Sergey,

As mecej4 already hinted, C programming language stores the 2D arrays in row-major fashion. If you would like to use 2D C arrays for dgemm, you need to call cblas_dgemm in row-major ordering. Namely, cblas_dgemm(order=CblasRowMajor, transa, transb, m, n, k ...). If you just call dgemm, it assumes column-major matrix ordering.

I created a simple program which calls cblas_dgemm and dgemm for your input values. As you can see from the output, cblas_dgemm gives expected values, whereas, dgemm call gives the "wrong" results mentioned in your post. However, the results are not actually wrong, it's just the dgemm uses what you think as rows of your matrix as columns. I'm attaching the simple program.

Here is the output from it:

Results from CBLAS row-major ordering
70, 80, 120, 160, 200, 240, 280, 320,
42, 72, 108, 144, 180, 216, 252, 288,
42, 72, 108, 144, 180, 216, 252, 288,
42, 72, 108, 144, 180, 216, 252, 288,
42, 72, 108, 144, 180, 216, 252, 288,
42, 72, 108, 144, 180, 216, 252, 288,
42, 72, 108, 144, 180, 216, 252, 288,
42, 72, 108, 144, 180, 216, 252, 288,

Results from Fortran-style col-major ordering
70, 84, 126, 168, 210, 252, 294, 336,
40, 72, 108, 144, 180, 216, 252, 288,
40, 72, 108, 144, 180, 216, 252, 288,
40, 72, 108, 144, 180, 216, 252, 288,
40, 72, 108, 144, 180, 216, 252, 288,
40, 72, 108, 144, 180, 216, 252, 288,
40, 72, 108, 144, 180, 216, 252, 288,
40, 72, 108, 144, 180, 216, 252, 288,

Best Wishes,

Efe

0 Kudos
SergeyKostrov
Valued Contributor II
403 Views
>>...However, the results are not actually wrong, it's just the dgemm uses what you think as rows of your matrix as columns... Hi Efe, Thank you very much and this is exactly what I expected to confirm. Best regards, Sergey
0 Kudos
Abdul__J_
Beginner
403 Views

Hi,

can you tell me that what would be the values for leading dimensions lda, ldb and ldc in case of non-square matrix-multiplication? If A is m x n B is n x k then C will be m x k so what would be lda, ldb and ldc for cblas_gemm?

0 Kudos
mecej4
Honored Contributor III
403 Views

Abdul J wrote:
can you tell me that what would be the values for leading dimensions lda, ldb and ldc 

There is no single answer to this question, as these dimensions depend on whether the matrices are stored by rows or by columns. Please see the example .../mkl/examples/cblas/source/cblas_sgemmx.c.

If the matrix is stored by rows, as is more common in C, lda should be equal to the max. possible row length = max. number of columns of A. If you allocate the matrices to be the exact size needed and not larger, lda = the the number of columns of A.

And, it goes similarly for the other matrices.

0 Kudos
Reply