#include "sopnamsp.h"
#include "machdefs.h"
#include "JTC/JTC.h"

#include <math.h>
#include <iostream>

#include "tarrinit.h"
#include "array.h"
#include "timing.h"

// Programme de test des threads JThreadsC++ 
// avec les arrays de Sophya

// A global monitor for print synchronisation 
JTCMonitor prtmon;  // <JThreadC++> 

// A thread class for doing matrix computation 
class MtxComputer : public JTCThread   // <JThreadC++> 
{
public:
  MtxComputer(int id, int sz, int nloop);
  virtual void run();
protected:
  int id_;
  int sz_;
  int nloop_;
};

MtxComputer::MtxComputer(int id, int sz, int nloop)
{
  id_ = id;  sz_ = sz;  nloop_ = nloop;  
}

void MtxComputer::run()
{
  int n = sz_;
  double seuil = 1.e-6;
  Matrix id;
  id = IdentityMatrix(1.,n);
  double gmax = -1.e99;
  double gmin = 1.e99;
  int npb = 0;
  // Loop creating a random matrix, inverting it
  // and checking the result 
  for(int k=0; k<nloop_; k++) {
    try {
      Matrix a(n,n);
      a = Sequence(RandomSequence(RandomSequence::Gaussian, 0., 2.5));
      Matrix inva = Inverse(a);
      Matrix diff = id-a*inva;
      double max = -1.e99;
      double min = 1.e99;
      double x;
      int nerr = 0;
      for(int i=0; i<n; i++) 
	for(int j=0; j<n; j++) {
	  x = diff(i,j);
	  if (x < min)  min = diff(i,j); 
	  if (x > max)  max = diff(i,j); 
	  if (fabs(x) > seuil)  nerr++;
	}
      if (min < gmin) gmin = min;
      if (max > gmax) gmax = max;
      if (nerr > 0) npb++;
      {  // Synchronized writing to cout stream  
	JTCSynchronized sync(prtmon);  // <JThreadC++> 
	cout << " ------- Thread[" << id_ << "] K= " 
	     << k << "  NErr = " << nerr << endl;
	cout << "  Min(Diff) = " << min << " Max(Diff) = " << max << endl;
	if (k == nloop_-1) {
	  double frac = (double)npb*100./(double)nloop_;
	  cout << " ...... Thread[" << id_ << "] End NPb= " << npb 
	       << " / NTot= " << nloop_ << " ( = " << frac << " %) " << endl;
	cout << "  GMin(Diff) = " << gmin << " GMax(Diff) = " << gmax << endl;
	cout << " ..................................................... " << endl;
	}
      } 
    }
    catch (PException const & e) {
	cerr << " Catched PException in  Thread(" << (string)getName() 
	     << ") Msg= " << e.Msg() << endl;
    }
    catch (JTCException const & e) {  // <JThreadC++>
	cerr << " Catched JTCException in  Thread(" << (string)getName() 
	     << ") Msg= " << e.getMessage() << endl;
    }
  }
}


// ------- The main program ---------
int main(int narg, char* arg[])
{

#define MAXNTHR 10    
  if (narg < 4) {
    cout << " jtctarr - JThreadsC++/Sophya::TArray<T> Test \n " 
	 << " ... Usage jtctarr NThreads NLoop MtxSize [-JTCss ... \n" 
	 << " NThreads=1...4  NLoop=10...10^4  MtxSize=10...10^3 \n"
	 << endl;
    exit(0);
  }

  InitTim();   // Initializing the CPU timer

  int nthr = atoi(arg[1]);
  if (nthr > MAXNTHR) nthr = MAXNTHR;
  int nloop = atoi(arg[2]);
  int size = atoi(arg[3]);
  cout << " ::::: jtctarr - JThreadsC++/Sophya::TArray<T> Test ::::: \n"
       << " NThreads= " << nthr << " NLoop= " << nloop << " MtxSize= " 
       << size << endl;   

  try {
    SophyaInit();  // Sophya Initialization
    JTCInitialize iniJTC(narg, arg);  // <JThreadC++>  Initialize library

    JTCThreadHandle t[MAXNTHR];  // <JThreadC++>  
    int k;
    for(k=0; k<nthr; k++)  t[k] = new MtxComputer(k, size, nloop);
    cout << " Starting threads ... " << endl;
    for(k=0; k<nthr; k++) t[k]->start();  // <JThreadC++> 
    // Waiting for threads to end
    for(k=0; k<nthr; k++) t[k]->join();   // <JThreadC++>
  }
  catch (PThrowable const & e) {
    cerr << " Catched PThrowable in main Msg= " << e.Msg() << endl;
  }
  catch (JTCException const & e) {  // <JThreadC++>
    cerr << " Catched JTCException in main Msg= " << e.getMessage() << endl;
  }
  catch(...) {
    cerr << " Catched ... exception in main " << endl;
  }

  PrtTim("End of jtctarr");
}