/* 
   matrix3.c -  Matrix Multiplication
                matrix mult expressed by transpose and dot product

   compile: gcc -Wall -O -o matrix3 matrix3.c
   run:     ./matrix3 file
*/

#include <stdio.h>
#include <stdlib.h>
#include <time.h>

double startTime, stopTime;

/* Timing ------------------------------------------------------- */
void showElapsed(int id, char * m)
{  printf("%d: %s %f secs\n",id,m,(clock()-startTime)/CLOCKS_PER_SEC);
}

/* Aux fcts ------------------------------------------------------- */
int * allocVector(int n)
{
  return (int *)malloc(n*sizeof(int));
}

int ** allocMatrix(int m, int n)
{
  int ** newM = (int **)malloc(m*sizeof(int *));
  int i;
  for (i = 0; i < m; i++)
    newM[i] = allocVector(n);
  return newM;
}

void readMatrix(FILE * fin, int ** M, int m, int n)
{
  int i, j;
  for (i = 0; i < m; i++)
    for (j = 0; j < n; j++)
      fscanf(fin, "%d", &(M[i][j]));
}

void writeMatrix(FILE * fout, int ** M, int m, int n)
{
  int i, j;
  for (i = 0; i < m; i++) {
    for (j = 0; j < n; j++)
      fprintf(fout, "%d ", M[i][j]);
    putc('\n', fout);
  }
}


/* Aux fct: transpose m*n matrix into n*m matrix */
int ** transpose(int ** M, int m, int n)
{
  int ** MT;
  int i, j;
  MT = allocMatrix(n, m);
  for (i = 0; i < m; i++)
    for (j = 0; j < n; j++)
      MT[j][i] = M[i][j];
  return MT;
}


/* Aux fct: dot product of two vectors of length n */
int dotProd(int * V1, int * V2, int n)
{
  int dp = 0;
  int i;
  for (i = 0; i < n; i++)
    dp = dp + V1[i] * V2[i];
  return dp;
}


/* Matrix multiplication; second argument assumed to be transposed */
int ** matrixProd(int ** M1, int ** M2T, int m, int n)
{
  int i, j;
  int ** M3 = allocMatrix(m, m);
  for (i = 0; i < m; i++)
    for (j = 0; j < m; j++)
      M3[i][j] = dotProd(M1[i], M2T[j], n);
  return M3;
}

/* Main fct ------------------------------------------------------- */

int main(int argc, char ** argv)
{
  FILE * fin;
  FILE * fout;
  int ** M1;
  int ** M2;
  int ** M2T;
  int ** M3;
  int m, n;

  fin = fopen(argv[1], "r");
  fscanf(fin, "%d %d", &m, &n);
  M1 = allocMatrix(m, n);
  readMatrix(fin, M1, m, n);
  M2 = allocMatrix(n, m);
  readMatrix(fin, M2, n, m);
  fclose(fin);

  startTime = clock();
  M2T = transpose(M2, n, m);
  // writeMatrix(stdout, M2T, m, n);

  M3 = matrixProd(M1, M2T, m, n);
  stopTime = clock();
  printf("%d * %d; SEQUENTIAL; %f secs\n",
	  m,n, (stopTime-startTime)/CLOCKS_PER_SEC);
  // Print output:
  fout = fopen("/dev/null", "w");
  //writeMatrix(stdout,M3,m,m);
  /*
  fout = fopen("MRESULT_SEQ3","w");
  writeMatrix(fout,M3,m,m);
  fclose(fout);
  */

  return 0;
}
