How to write a C Program to Sequential Matrix Multiplication in C Programming Language ?
Solution:
/*sequential matrix multiplication*/
#include <stdio.h>
#include <stdlib.h>
#include <omp.h>
#define NRA 100 // number of rows in matrix A
#define NCA 100 // number of columns in matrix A
#define NCB 100 // number of columns in matrix B
// timing routine
double read_timer();
////////////////////////////////////////////////////
int main (int argc, char *argv[]) {
int tid, nthreads, i, j, k;
double **a, **b, **c;
double *a_block, *b_block, *c_block;
double **res;
double *res_block;
double starttime, stoptime;
double x,y;
a = (double **) malloc(NRA*sizeof(double *)); /* matrix a to be multiplied */
b = (double **) malloc(NCA*sizeof(double *)); /* matrix b to be multiplied */
c = (double **) malloc(NRA*sizeof(double *)); /* result matrix c */
a_block = (double *) malloc(NRA*NCA*sizeof(double)); /* Storage for matrices */
b_block = (double *) malloc(NCA*NCB*sizeof(double));
c_block = (double *) malloc(NRA*NCB*sizeof(double));
/* Result matrix for the sequential algorithm */
res = (double **) malloc(NRA*sizeof(double *));
res_block = (double *) malloc(NRA*NCB*sizeof(double));
for (i=0; i<NRA; i++) /* Initialize pointers to a */
a[i] = a_block+i*NRA;
for (i=0; i<NCA; i++) /* Initialize pointers to b */
b[i] = b_block+i*NCA;
for (i=0; i<NRA; i++) /* Initialize pointers to c */
c[i] = c_block+i*NRA;
for (i=0; i<NRA; i++) /* Initialize pointers to res */
res[i] = res_block+i*NRA;
printf("Starting sequential matrix multiplication\n");
printf("Initializing matrices...\n");
/*** Initialize matrices ***/
for (i=0; i<NRA; i++) /* last matrix has been initialized */
for (j=0; j<NCA; j++)
a[i][j]= (double) (i+j);
for (i=0; i<NCA; i++)
for (j=0; j<NCB; j++)
b[i][j]= (double) (i*j);
for (i=0; i<NRA; i++)
for (j=0; j<NCB; j++)
c[i][j]= 0.0;
starttime = read_timer();
// #pragma omp parallel private(i)
//{
//#pragma omp parallel for private(i,x,y)
for(i=0; i<NRA; i++) {
//#pragma omp parallel for
#pragma omp parallel for private(x,y)
for(k=0; k<NCB; k++) {
for(j=0; j<NCA; j++) {
//x = a[i][k];
//y = b[k][j];
//c[i][j] = x * y;
//c[i][j] += a[i][k] * b[k][j];
x += a[i][k] * b[k][j];
}
c[i][j] = x;
}
}
//}
stoptime = read_timer();
printf("Time for sequential matrix multiplication: %3.8f s\n",
stoptime-starttime);
printf ("Done.\n");
exit(0);
}
// timing routine implementation
double read_timer( )
{
static int initialized = 0;
static struct timeval start;
struct timeval end;
if( !initialized )
{
gettimeofday( &start, NULL );
initialized = 1;
}
gettimeofday( &end, NULL );
return (end.tv_sec - start.tv_sec) + 1.0e-6 * (end.tv_usec - start.tv_usec);
}