UNB/ CS/ David Bremner/ teaching/ cs3383/ examples/ row-mult.cc
#include <random>
#include <vector>
#include <iostream>
#include <cilk/cilk.h>
#include <cilk/cilk_api.h>
#include <sys/time.h>
#include <unistd.h>
#include <climits>
#include "mat.hpp"

using namespace std;
double get_wall_time(){
    struct timeval time;
    if (gettimeofday(&time,NULL)){
      abort();
    }
    return (double)time.tv_sec + (double)time.tv_usec * .000001;
}

void
RowMult(const mat &A, const vec &x, vec &y, int i) {
  y[i]=0;
  for (int j=0; j < A[i].size(); j++){
    y[i] += A[i][j]*x[j];
  }
}

void MatVec(const mat &A, const vec &x, vec &y) {
  cilk_for(int i=0; i<A.size(); i++){
    RowMult(A, x, y, i);
  }
}

void usage (int argc, char** argv) {
  cerr << "usage: " << argv[0] << " [-d] [-s seed] <rows> <cols> " << endl;
  exit(1);
}

int main(int argc, char** argv){

  int rows,cols;

  int opt;
  int debug = 0, seed_found = 0;
  long seed;

  while ((opt = getopt(argc, argv, "ds:")) != -1) {
    switch (opt) {
    case 'd':
      debug = true;
      break;
    case 's':
      seed = atol (optarg);
      seed_found = 1;
      break;
    default:
      usage(argc, argv);
    }
  }

  if (optind + 1 >= argc)
    usage(argc, argv);

  rows = atol(argv[optind]);
  cols = atol(argv[optind+1]);
  if (argv[optind+2])
    __cilkrts_set_param("nworkers",argv[optind+2]);

  if (! seed_found) {
    random_device rd;
    seed = rd();
  }

  double start = get_wall_time();

  mat A = random_mat(seed, rows, cols);
  vec x,y;

  if (debug)
    cout << "x = " << endl;

  mt19937 gen(seed+rows+cols);
  for(int j=0; j < cols; j++) {
    scalar val = dist(gen());
    x.push_back(val);
    if (debug) {
      cout << " " << val;
    }
  }
  if (debug)
    cout << endl;

  y.resize(rows);

  double init = get_wall_time();

  MatVec(A,x,y);

  double end = get_wall_time();

  cout << "init = " << (double)(init - start) << endl << "mult = " << (double)(end - init)  << endl;

  if (debug) {
    cout << "y =" << endl;
    for(int j=0; j < rows; j++) {
      cout << " "  << y[j];
    }
    cout << endl;
  }

}