how to implement matrix multiplication in C?-CodePudding

I was trying to create an algorithm to perform matrix multiplication.

I've designed matrix as follows:

// matrix.h  
 #pragma once  
    #include <stdlib.h>
    #include <string.h>
    struct matrix {
        size_t rows, cols;
        double* data;
    };
    extern struct matrix* mat_mul(const struct matrix* m1, const struct matrix* m2);

// matrix.c 
#include "matrix.h"

void mat_constr(struct matrix* m, size_t rows, size_t cols) {
    m->rows = rows; m->cols = cols; 
    m->data = calloc(rows * cols, sizeof(double)); 
    if (!m->data) {
        return; 
    }
}

void mat_destr(struct matrix* m) {
    free(m->data); 
}

mat_constr is matrix constructor, and mat_destr is mat_destructor. To test the algorithm I've used this main

 // main
 int main(void) {
        struct matrix A; 
        mat_constr(&A, 2, 3); 
        memcpy(A.data, (double[6]) { 1, 2, 3, 4, 5, 6 }, 6 * sizeof(double)); 
        struct matrix B; 
        mat_constr(&B, 3, 2); 
        memcpy(B.data, (double[6]) { 7, 8, 9, 10, 11, 12 }, 6 * sizeof(double)); 
        struct matrix* C = mat_mul(&A, &B); 
        mat_destr(&A); mat_destr(&B); 
        mat_destr(C); 
        return 0; 
    }

and this is the mat_mul function

struct matrix* mat_mul(const struct matrix* m1, const struct matrix* m2) {
    if ((m1 == NULL) || (m2 == NULL)) {
        return NULL; 
    }
    if (m1->cols != m2->rows) {
        return NULL; 
    }
    struct matrix* result = malloc(sizeof(struct matrix)); 
    if (!result) {
        return NULL; 
    }
    mat_constr(result, m1->rows, m2->cols); 

    
    size_t k = 1; 
    for (size_t r = 0; r < m1->rows; r  ) {
        for (size_t c = 0; c < m1->cols; c  ) {
            result->data[r * result->cols   c] = m1->data[r * m1->cols   k] * m2->data[k * m2->cols   c]; 
        }
        k  ; 
        }
    return result; 
}

in order to perform matrix multiplication, I have to use this sum: sum from k = 1 to m1->cols of a_i k-th column * a_j k-th row (in this forum I don't know how to write using mathjax, because symbols like this $$ $$ doesn't work here).

this is the minimal reproducible example:

// matrix.h
    #pragma once  
    #include <stdlib.h>
    #include <string.h>
    struct matrix {
        size_t rows, cols;
        double* data;
    };
    extern struct matrix* mat_mul(const struct matrix* m1, const struct matrix* m2);


// matrix.c  
    #include "matrix.h"
    
    void mat_constr(struct matrix* m, size_t rows, size_t cols) {
        m->rows = rows; m->cols = cols; 
        m->data = calloc(rows * cols, sizeof(double)); 
        if (!m->data) {
            return; 
        }
    }
    
    void mat_destr(struct matrix* m) {
        free(m->data); 
    }
    
    struct matrix* mat_mul(const struct matrix* m1, const struct matrix* m2) {
        if ((m1 == NULL) || (m2 == NULL)) {
            return NULL; 
        }
        if (m1->cols != m2->rows) {
            return NULL; 
        }
        struct matrix* result = malloc(sizeof(struct matrix)); 
        if (!result) {
            return NULL; 
        }
        mat_constr(result, m1->rows, m2->cols); 
    
        
        size_t k = 1; 
        for (size_t r = 0; r < m1->rows; r  ) {
            for (size_t c = 0; c < m1->cols; c  ) {
                result->data[r * result->cols   c] = m1->data[r * m1->cols   k] * m2->data[k * m2->cols   c]; 
            }
            k  ; 
            }
            
    
    
    
        return result; 
    
    }
    
    
    
    
    int main(void) {
        struct matrix A; 
        mat_constr(&A, 2, 3); 
        memcpy(A.data, (double[6]) { 1, 2, 3, 4, 5, 6 }, 6 * sizeof(double)); 
        struct matrix B; 
        mat_constr(&B, 3, 2); 
        memcpy(B.data, (double[6]) { 7, 8, 9, 10, 11, 12 }, 6 * sizeof(double)); 
        struct matrix* C = mat_mul(&A, &B); 
        mat_destr(&A); mat_destr(&B); 
        mat_destr(C); 
        return 0; 
    }

this solution allocates enough memory, and return the pointer of the new allocated matrix correctly. But the problem is in the last for-loops. According to my linear algebra knowledge, I have to scroll columns by columns the first matrix, and scroll rows by rows the second matrix. But these for-loops are not correct.

I have only one question: "why is this method of computing matrix multiplication wrong? how can I solve it?"

note that r * cols c gives exactly the index of the i-th entry of the matrix.

CodePudding user response：

You are simply missing a nested for loop and a =; The correct code for multiplying your 2 matrices would be something like:

for (size_t r = 0; r < m1->rows;   r)
    for (size_t c = 0; c < m2->cols;   c)
        for (size_t k = 0; k < m2->rows;   k)
            result->data[r * result->cols   c]  = m1->data[r * m1->cols   k] * m2->data[k * m2->cols   c]