Cblas_sgemm produces wrong results

Question

I have a problem using the sgemm function of cblas.

Here is the code:

#include <math.h>
#include <stdio.h>
#include <stdlib.h>
#include <string.h>
#include <cblas.h>

#define MATRIX_DIM 5

int main(){

    float *a_mat = calloc(MATRIX_DIM*MATRIX_DIM, sizeof(float));
    float *b_mat = calloc(MATRIX_DIM, sizeof(float));
    float *c_mat = calloc(MATRIX_DIM, sizeof(float));
    int i,j;

    for(i=0; i<MATRIX_DIM*MATRIX_DIM; i++) {
        a_mat[i] = 1.0f;
        b_mat[i] = 1.0f;
        c_mat[i] = 0.0f;
    }

    cblas_sgemm(CblasRowMajor, CblasNoTrans,
                 CblasNoTrans, MATRIX_DIM, MATRIX_DIM,
                 MATRIX_DIM, 1.0, a_mat,
                 MATRIX_DIM, b_mat, MATRIX_DIM,
                 1.0, c_mat, MATRIX_DIM);

    //RESULT PRINTING
    printf("Printing A MATRIX: \n");
    for(i=0; i<MATRIX_DIM; i++) {
        for(j=0; j<MATRIX_DIM; j++){
                printf("%0.1f ", a_mat[i*MATRIX_DIM+j]);
        }
        printf("\n");
    }
    printf("Printing B MATRIX: \n");
    for(i=0; i<MATRIX_DIM; i++) {
            for(j=0; j<MATRIX_DIM; j++){
                printf("%0.1f ", b_mat[i*MATRIX_DIM+j]);
        }
        printf("\n");
    }

    printf("\nPrinting the Results: \n");
    for(i=0; i<MATRIX_DIM;i++){
        for(j=0; j<MATRIX_DIM; j++){
                printf("%0.1f ", c_mat[i*MATRIX_DIM+j]);
        }
        printf("\n");
    }

    free(a_mat);
    free(b_mat);
    free(c_mat);

    return 0;
}

I am pretty some of the arguments I put in are wrong, but I really don't know which. The results should be a 5x5 matrix filled with 5.0. Instead the program responds with this:

6.0 6.0 6.0 16.0 86.0 
6.0 6.0 6.0 16.0 86.0 
16.0 36.0 6.0 46.0 86.0 
16.0 36.0 5.0 45.0 85.0 
20.0 80.0 5.0 45.0 85.0

I know the rowmajor order or the transpose arguments might be wrong and I will figure those out later, but in this particular multiplication the answer should be 5.0 either way.

Hmm... Shouldn't you `calloc` for `sizeof(float)`? How is `sizeof(*a_mat)` even valid? And shouldn't your result be a *vector* of length 5, all elements equal to 5? And shouldn't that result be written in *the vector* `c_mat`? How are you not getting any segfaults? — Andras Deak -- Слава Україні, Jan 29 '16 at 23:48
@AndrasDeak Well, you might be right about the `sizeof(float)` part but for some reason it works like this just fine. Nevertheless I will change it. Since it is Matrix-Matrix Multiplication the result should be a matrix of MxN (M being the rows of the first matrix and N being the columns of the second Matrix) which would be a 5x5 matrix in this case. I think you are referring to Matrix-Vector Multiplication. — kmentis, Jan 29 '16 at 23:56
Where can I see that the second matrix has 5 columns? `float *b_mat = calloc(MATRIX_DIM, sizeof(*b_mat));` is only a columnful of values allocated. Hence my surprise that you don't get a segfault a bit later. One possibility is that `sizeof(*a_mat)` (whatever that is) has more room than a `float`. And your problem is exactly that it's *not* working fine:) — Andras Deak -- Слава Україні, Jan 30 '16 at 00:00
Ok you are right! The problem was I forgot to add one more `*MATRIX_DIM` on the b_mat and c_mat. I copy-pasted the code from my previous implementation of sgemv and totally forgot about it. The `sizeof(*a_mat)` still worked btw, I am not sure if it is a problem though, so I changed it to what you said aswell. Thank you for pointing it out! :) — kmentis, Jan 30 '16 at 00:07

kmentis · Accepted Answer · 2016-01-30T00:37:59.127

0

Thanks to @AndrasDeak in the comments, all I needed was to allocate more space on the two matrices, which I previously overlooked.

So basically change:

float *b_mat = calloc(MATRIX_DIM, sizeof(float));
float *c_mat = calloc(MATRIX_DIM, sizeof(float));

To:

float *b_mat = calloc(MATRIX_DIM*MATRIX_DIM, sizeof(float));
float *c_mat = calloc(MATRIX_DIM*MATRIX_DIM, sizeof(float));

Since those are supposed to be 2-dimensional matrices and not vectors.

edited Jan 30 '16 at 00:37

answered Jan 30 '16 at 00:09

kmentis

23
5

I'd still change those `sizeof()`s to float:) And since this is essentially a typo in the code, you should consider whether future users can find this Q&A useful. If you *don't* think so, than it might be worth deleting the question. If you decide to keep the Q&A up, then accept your answer after a while when you'll be able to (24 or 48 hours, I think) – Andras Deak -- Слава Україні Jan 30 '16 at 00:13
I changed the `calloc` to float like you said, should I change the `*MATRIX_DIM` part as well to have it fully functioning? I am new to stackoverflow, so I am not familiar with the procedure. – kmentis Jan 30 '16 at 00:32
Ooooh no no no, sorry, you misunderstood:) I meant just changing your actual code, but not your question above (your answer wasn't clear about whether you changed the `sizeof` as well, that's why I noted that too in my comment). The procedure is "leave the question as-is, unless adding clarification", and any fixes should go into answers. Then everything is present exactly once: one clear question statement, and one (or more) clear answers. – Andras Deak -- Слава Україні Jan 30 '16 at 00:33
Alright I will leave it with the changed `sizeof(float)` then. Thanks! :) – kmentis Jan 30 '16 at 00:34

Cblas_sgemm produces wrong results

1 Answers1