# What's wrong with this C code for Matrix Multiplication using Strassen's Algorithm?

``````// Probably some mistake in "conquer" part of the recursive function "multiplySquareMatrices".

#include<stdio.h>

void showSquareMatrix(int n, int[][n]);
void copyMatrix_Complete2Part(int n, int[][n], int, int, int, int, int m,
int[][m]);
void subSquareMatrices(int n, int[][n], int[][n], int[][n]);
void addSquareMatrices(int n, int[][n], int[][n], int[][n]);
void copyMatrix_Part2Complete(int n, int[][n], int m, int[][m], int, int, int, int);
void multiplySquareMatrices(int n, int[][n], int[][n], int[][n]);
void acceptSquareMatrix(int n, int[][n]);

int main()
{
int n = 2;
int matrixA[n][n];
int matrixB[n][n];
int matrixAB[n][n];

acceptSquareMatrix(n, matrixA);
acceptSquareMatrix(n, matrixB);

multiplySquareMatrices(n, matrixAB, matrixA, matrixB);
showSquareMatrix(n, matrixAB);
}

void acceptSquareMatrix(int n, int matrix[][n])     // matrix[][columns] should not be changed.
{
int i;
int j;
for(i = 0; i < n; i++)
{
for(j = 0; j < n; j++)
{
printf("\n[%d][%d]: ", i, j);
scanf("%d", &matrix[i][j]);
}
}
}

void multiplySquareMatrices(int n, int matrixAB[][n], int matrixA[][n], int matrixB[][n])
{
if(n == 1)
{
matrixAB[0][0] = matrixA[0][0] * matrixB[0][0];
}
else
{
//divide
int matrix_a[n / 2][n / 2];
int matrix_b[n / 2][n / 2];
int matrix_c[n / 2][n / 2];
int matrix_d[n / 2][n / 2];

int matrix_e[n / 2][n / 2];
int matrix_f[n / 2][n / 2];
int matrix_g[n / 2][n / 2];
int matrix_h[n / 2][n / 2];

int rowStart_a = 0;     int rowEnd_a = (n/2) - 1;     int columnStart_a = 0;        int columnEnd_a = (n/2) - 1;
int rowStart_b = 0;     int rowEnd_b = (n/2) - 1;     int columnStart_b = n/2;      int columnEnd_b = n - 1;
int rowStart_c = n/2;   int rowEnd_c = n - 1;         int columnStart_c = 0;        int columnEnd_c = (n/2) - 1;
int rowStart_d = n/2;   int rowEnd_d = n - 1;         int columnStart_d = n/2;      int columnEnd_d = n - 1;

int rowStart_e = 0;     int rowEnd_e = (n/2) - 1;     int columnStart_e = 0;        int columnEnd_e = (n/2) - 1;
int rowStart_f = 0;     int rowEnd_f = (n/2) - 1;     int columnStart_f = n/2;      int columnEnd_f = n - 1;
int rowStart_g = n/2;   int rowEnd_g = n - 1;         int columnStart_g = 0;        int columnEnd_g = (n/2) - 1;
int rowStart_h = n/2;   int rowEnd_h = n - 1;         int columnStart_h = n/2;      int columnEnd_h = n - 1;

copyMatrix_Part2Complete(n / 2, matrix_a, n, matrixA, rowStart_a, rowEnd_a, columnStart_a, columnEnd_a);
copyMatrix_Part2Complete(n / 2, matrix_b, n, matrixA, rowStart_b, rowEnd_b, columnStart_b, columnEnd_b);
copyMatrix_Part2Complete(n / 2, matrix_c, n, matrixA, rowStart_c, rowEnd_c, columnStart_c, columnEnd_c);
copyMatrix_Part2Complete(n / 2, matrix_d, n, matrixA, rowStart_d, rowEnd_d, columnStart_d, columnEnd_d);

copyMatrix_Part2Complete(n / 2, matrix_e, n, matrixB, rowStart_e, rowEnd_e, columnStart_e, columnEnd_e);
copyMatrix_Part2Complete(n / 2, matrix_f, n, matrixB, rowStart_f, rowEnd_f, columnStart_f, columnEnd_f);
copyMatrix_Part2Complete(n / 2, matrix_g, n, matrixB, rowStart_g, rowEnd_g, columnStart_g, columnEnd_g);
copyMatrix_Part2Complete(n / 2, matrix_h, n, matrixB, rowStart_h, rowEnd_h, columnStart_h, columnEnd_h);

//conquer
/*
//conquer for divide and conquer algorithm
int matrix_r[n / 2][n / 2];
int matrix_s[n / 2][n / 2];
int matrix_t[n / 2][n / 2];
int matrix_u[n / 2][n / 2];

int matrix_ae[n / 2][n / 2];
int matrix_bg[n / 2][n / 2];
int matrix_af[n / 2][n / 2];
int matrix_bh[n / 2][n / 2];
int matrix_ce[n / 2][n / 2];
int matrix_dg[n / 2][n / 2];
int matrix_cf[n / 2][n / 2];
int matrix_dh[n / 2][n / 2];

multiplySquareMatrices(n / 2, matrix_ae, matrix_a, matrix_e);
multiplySquareMatrices(n / 2, matrix_bg, matrix_b, matrix_g);
multiplySquareMatrices(n / 2, matrix_af, matrix_a, matrix_f);
multiplySquareMatrices(n / 2, matrix_bh, matrix_b, matrix_h);
multiplySquareMatrices(n / 2, matrix_ce, matrix_c, matrix_e);
multiplySquareMatrices(n / 2, matrix_dg, matrix_d, matrix_g);
multiplySquareMatrices(n / 2, matrix_cf, matrix_c, matrix_f);
multiplySquareMatrices(n / 2, matrix_dh, matrix_d, matrix_h);

addSquareMatrices(n / 2, matrix_r, matrix_ae, matrix_bg);
addSquareMatrices(n / 2, matrix_s, matrix_af, matrix_bh);
addSquareMatrices(n / 2, matrix_t, matrix_ce, matrix_dg);
addSquareMatrices(n / 2, matrix_u, matrix_cf, matrix_dh);
*/

//conquer for Strassen's algorithm
int matrix_fSubH[n / 2][n / 2];
int matrix_aAddB[n / 2][n / 2];
int matrix_cAddD[n / 2][n / 2];
int matrix_gSubE[n / 2][n / 2];
int matrix_aAddD[n / 2][n / 2];
int matrix_eAddH[n / 2][n / 2];
int matrix_bSubD[n / 2][n / 2];
int matrix_gAddH[n / 2][n / 2];
int matrix_aSubC[n / 2][n / 2];
int matrix_eAddF[n / 2][n / 2];

subSquareMatrices(n / 2, matrix_fSubH, matrix_f, matrix_h); //S1 = B12 - B22
subSquareMatrices(n / 2, matrix_gSubE, matrix_g, matrix_e); //S4 = B21 - B11
subSquareMatrices(n / 2, matrix_bSubD, matrix_b, matrix_d); //S7 = A12 - A22
subSquareMatrices(n / 2, matrix_aSubC, matrix_a, matrix_c); //S9 = A11 - A21

int matrix_p1[n / 2][n / 2];
int matrix_p2[n / 2][n / 2];
int matrix_p3[n / 2][n / 2];
int matrix_p4[n / 2][n / 2];
int matrix_p5[n / 2][n / 2];
int matrix_p6[n / 2][n / 2];
int matrix_p7[n / 2][n / 2];

multiplySquareMatrices(n / 2, matrix_p1, matrix_a,     matrix_fSubH); //P1 = A11*S1
multiplySquareMatrices(n / 2, matrix_p2, matrix_aAddB, matrix_h);     //P2 = S2*B22
multiplySquareMatrices(n / 2, matrix_p3, matrix_cAddD, matrix_e);     //P3 = S3*B11
multiplySquareMatrices(n / 2, matrix_p4, matrix_d,     matrix_gSubE); //P4 = A22*S4
multiplySquareMatrices(n / 2, matrix_p6, matrix_bSubD, matrix_gAddH); //P6 = S7*S8
multiplySquareMatrices(n / 2, matrix_p7, matrix_aSubC, matrix_eAddF); //P7 = S9*S10

int matrix_r[n / 2][n / 2];
int matrix_s[n / 2][n / 2];
int matrix_t[n / 2][n / 2];
int matrix_u[n / 2][n / 2];

int matrix_rTemp1[n / 2][n / 2];
int matrix_rTemp2[n / 2][n / 2];
int matrix_uTemp1[n / 2][n / 2];
int matrix_uTemp2[n / 2][n / 2];

addSquareMatrices(n / 2, matrix_rTemp1,        matrix_p5,            matrix_p4); //P5 + P4
subSquareMatrices(n / 2, matrix_rTemp2,        matrix_rTemp1,        matrix_p2); //(P5 + P4) - P2
addSquareMatrices(n / 2, matrix_r,             matrix_rTemp2,        matrix_p6); //C11 = (P5 + P4 - P2) + P6

addSquareMatrices(n / 2, matrix_s,             matrix_p1,            matrix_p2); //C12 = P1 + P2

addSquareMatrices(n / 2, matrix_t,             matrix_p3,            matrix_p4); //C21 = P3 + P4

addSquareMatrices(n / 2, matrix_uTemp1,        matrix_p5,            matrix_p1); //P5 + P1
subSquareMatrices(n / 2, matrix_uTemp2,        matrix_uTemp1,        matrix_p3); //(P5 + P1) - P3
subSquareMatrices(n / 2, matrix_u,             matrix_uTemp2,        matrix_p7); //C22 = (P5 + P1 - P3) - P7

//combine
int rowStart_r = 0;     int rowEnd_r = (n/2) - 1;     int columnStart_r = 0;        int columnEnd_r = (n/2) - 1;
int rowStart_s = 0;     int rowEnd_s = (n/2) - 1;     int columnStart_s = n/2;      int columnEnd_s = n - 1;
int rowStart_t = n/2;   int rowEnd_t = n - 1;         int columnStart_t = 0;        int columnEnd_t = (n/2) - 1;
int rowStart_u = n/2;   int rowEnd_u = n - 1;         int columnStart_u = n/2;      int columnEnd_u = n - 1;

copyMatrix_Complete2Part(n, matrixAB, rowStart_r, rowEnd_r, columnStart_r, columnEnd_r, (n / 2), matrix_r);
copyMatrix_Complete2Part(n, matrixAB, rowStart_s, rowEnd_s, columnStart_s, columnEnd_s, (n / 2), matrix_s);
copyMatrix_Complete2Part(n, matrixAB, rowStart_t, rowEnd_t, columnStart_t, columnEnd_t, (n / 2), matrix_t);
copyMatrix_Complete2Part(n, matrixAB, rowStart_u, rowEnd_u, columnStart_u, columnEnd_u, (n / 2), matrix_u);

}
}

void copyMatrix_Part2Complete(int destN, int destMatrix[][destN], int srcN,
int srcMatrix[][srcN], int rowS, int rowE, int colS, int colE)
{
int i; int j;
int k; int l;

for(i = rowS, k = 0; k < destN; i++, k++)
{
for(j = colS, l = 0; l < destN; j++, l++)
{
destMatrix[k][l] = srcMatrix[i][j];
}
}
}

{
int i; int j;
for(i = 0; i < n; i++)
{
for(j = 0; j < n; j++)
{
}
}
}

void subSquareMatrices(int n, int subtractedMatrix[][n], int matrixA[][n], int matrixB[][n])
{
int i; int j;
for(i = 0; j < n; i++)
{
for(j = 0; j < n; j++)
{
subtractedMatrix[i][j] = matrixA[i][j] - matrixB[i][j];
}
}
}

void copyMatrix_Complete2Part(int destN, int destMatrix[][destN], int rowStart, int rowEnd, int columnStart, int columnEnd, int srcN, int srcMatrix[][srcN])
{
int i; int j;
int k; int l;

for(i = 0, k = rowStart; k <= rowEnd; i++, k++)
{
for(j = 0, l = columnStart; l <= columnEnd; j++, l++)
{
destMatrix[k][l] = srcMatrix[i][j];
}
}
}

void showSquareMatrix(int n, int matrix[][n])
{
int i;
int j;
for(i = 0; i < n; i++)
{
printf("\n");
for(j = 0; j < n; j++)
{
printf("\t%d", matrix[i][j]);
}
}
}
``````

In the conquer part of the recurseve fuction "multiplySquareMatrices", there may be some problem.

If I replace the "conquer for Strassen's algorithm" part of the code with "conquer for the Divide & Conquer algorithm" by decommenting the "conquer for the divide & conquer algorithm" part & commenting the "conquer for Strassen's algorithm" part, the code works. I'm unable to understand why?

calgorithmmatrix-multiplication

answered 3 months ago rafix07 #1

In conquer for the divide & conquer algorithm algorithm you don't use `subSquareMatrices` function. Now look at definition of this function

``````int i; int j;
for(i = 0; j < n; i++)  // <-------
{
for(j = 0; j < n; j++)
{
subtractedMatrix[i][j] = matrixA[i][j] - matrixB[i][j];
}
}
``````

`j` should be replaced by `i` in outer for loop.