/* #include "machdefs.h" */

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


/*  ---------------------------------------------------  */
/*  --- Petit programme d'estimation de la puissance CPU */
/*  Compilation: csh> cc -O3 -o cpupower cpupower.c -lm  */
/* OU csh> cc -O -DT_Type=int cpupowerI  cpupower.c -lm  */
/*  R. Ansari   -   LAL   Mai 2004                       */
/*  ---------------------------------------------------  */

/*  Choix de type d'operations  float double int ... */
#ifndef T_Type
  #define  T_Type  double 
#endif

static int SZ=20000;                    /* Taille de tableau */
static double N_OP=0;               /* Nb operations */
static  int OPE=0;                  /* Choix Operation fop_1/2/3/4    */
static int ckprt=0;                 /* > 0 : Print check apres calcul */

/*  Tableaux X,Y,Z   */
static T_Type * x;              
static T_Type * y;
static T_Type * z;

void fop_0()
{
  T_Type s=0;
  int k;
  /*  printf("--- fop_0: Simple Loop: z[k] = x[k]*y[i] --- \n"); */
  for(k=0; k<SZ; k++) z[k] = x[k]*y[k];
  N_OP += SZ;
  return;
}
void fop_30()
{
  T_Type s=0;
  int i,k;
  /*  printf("--- fop_30: Simple Loop triple: z[k] = x[k]*y[i] --- \n");  */
  for(k=0; k<SZ; k+=3) { z[k] = x[k]*y[99]; z[k+1] = x[k+1]*y[99]; z[k+2] = x[k+2]*y[99]; } 
  N_OP += SZ;
  return;
}

void fop_1()
{
  T_Type s=0;
  int i,j,k;
  printf("--- fop_1: Double Loop: z[k] = Somme_i(x[k]*y[i]) --- \n");
  for(k=0; k<SZ; k++) {
    for(i=0; i<SZ; i++)  s += x[k]*y[i];
    z[k] = s;
    N_OP += 2*SZ;
  }
  return;
}

void fop_2() 
{
  T_Type s=0;
  int i,j,k;
  printf("--- fop_2: Double Loop: z[k] = Somme_i(x[k]*y[i]+x[i]*y[k]) --- \n");
  for(k=0; k<SZ; k++) {
    for(i=0; i<SZ; i++)  s += x[k]*y[i]+x[i]*y[k];
    z[k] = s;
    N_OP += 4*SZ;
  }
  return;
}

void fop_3() 
{
  T_Type s=0;
  int i,j,k;
  printf("--- fop_3: Double Loop: z[k] = Somme_i(x[k]*y[i]+x[i]*y[k]-0.85*(y[k]+x[k])) --- \n");
  for(k=0; k<SZ; k++) {
    for(i=0; i<SZ; i++)  s += x[k]*y[i]+x[i]*y[k]-0.85*(y[k]+x[k]);
    z[k] = s;
    N_OP += 7*SZ;
  }
  return;
}

void fop_4() 
{
  T_Type s=0;
  int i,j,k;
  printf("--- fop_4: Double Loop: z[k] = Somme_i(x[k]*sin(y[i])+y[k]*cos(x[i])) --- \n");
  for(k=0; k<SZ; k++) {
    for(i=0; i<SZ; i++)  s += x[k]*sin(y[i])+y[k]*cos(x[i]);
    z[k] = s;
    N_OP += 70*SZ;  /* le facteur 70 est approximatif */
  }
  return;
}

void fop_5() 
{
  T_Type s=0;
  int i,j,k;
  printf("--- fop_5: Double Loop: z[k] = Somme_i(x[k]*sin(y[i])+log(fabs(x[i])+0.1)+y[k]) --- \n");
  for(k=0; k<SZ; k++) {
    for(i=0; i<SZ; i++)  s += x[k]*sin(y[i])+log(fabs(x[i])+0.1)+y[k];
    z[k] = s;
    N_OP += 100*SZ;  /* le facteur 100 est approximatif */
  }
  return;
}

/* Fonctions de timing (TCPU) - voir en fin de fichier */
void InitTim(void);                                                             
void PrtTim(const char * Comm);                                                        
double GetPartialCPUTime();
double GetTotalCPUTime();

int main(int narg, char* arg[])
{
  int maxnprt = 1000;
  int nloop = 1000;
  int i,nprt;
  double mflops;

  /*  SophyaInit();   */
  InitTim();   /* Initializing the CPU timer  */

  if (narg < 2) {
    printf("--- Programme cpupower: (Puissance de calcul) ---- \n");
    printf("  Usage cpupower Op=0/1/2/3/4/5 [Size=20000] [NLoop=1000] [CkPrt=0]\n");
    printf("  Op=0 Simple Loop: z[k] = x[k]*y[k]  effectue NLoop fois\n");
    printf("    Size=10^6 par defaut, NLoop=100 \n");
    printf("  Op=1 Double Loop: z[k] = Somme_i(x[k]*y[i])  \n");
    printf("    Op=1 *,+=2 op / tour de boucle \n");
    printf("  Op=2 Double Loop: z[k] = Somme_i(x[k]*y[i]+x[i]*y[k]) \n");
    printf("    Op=2 2x(*,+)=4 op /  tour de boucle \n");
    printf("  Op=3 Double Loop: z[k] = Somme_i(x[k]*y[i]+x[i]*y[k]-0.85*(y[k]+x[k])) \n");
    printf("    Op=3  7 op (*,+) /  tour de boucle \n");
    printf("  Op=4 Double Loop: z[k] = Somme_i(x[k]*sin(y[i])+y[k]*cos(x[i])) \n");
    printf("    Op=4   ~ 70 op float /  tour de boucle \n");
    printf("  Op=5 Double Loop: z[k] = Somme_i(x[k]*sin(y[i])+log(fabs(x[i])+0.1)+y[k]) \n");
    printf("    Op=5   ~ 100 op float /  tour de boucle \n");
    printf("    Test 4,5 Size=2000 par defaut  \n");
    return 1;
  } 
  OPE = atoi(arg[1]);
  SZ = 20000;
  if (OPE == 0) SZ = 1000000;
  if (OPE >= 4) SZ = 2000;
  if (narg > 2)  SZ = atoi(arg[2]);
  nloop = 100;
  if (narg > 3)  nloop = atoi(arg[3]);
  ckprt = 0;
  if (narg > 4)  ckprt = atoi(arg[4]);

  printf("::::::: cpupower: OPE=%d SZ= %d ::::::: \n", OPE,SZ);
  x = malloc(SZ*sizeof(T_Type));
  y = malloc(SZ*sizeof(T_Type));
  z = malloc(SZ*sizeof(T_Type));

  for(i=0; i<SZ; i++) {
    x[i] = (T_Type)((random()%10000)-5000);  
    y[i] = (T_Type)((random()%10000)-5000);  
    z[i] = (T_Type)(0);
  }
  PrtTim("--Fin malloc+init x,y,z ");
  N_OP = 0;
  if (OPE == 5)  fop_5();
  else if (OPE == 4)  fop_4();
  else if (OPE == 3)  fop_3();
  else if (OPE == 2)  fop_2();
  else if (OPE == 1)  fop_1();
  else if (OPE == 30)  for(i=0; i<nloop; i++) fop_30();
  else for(i=0; i<nloop; i++) fop_0();
  PrtTim("---Fin OpeDoubleBoucle ");
  mflops = N_OP/ GetPartialCPUTime()*1.e-6;
  printf("-> Nb Operations= %g  MFLOPS= %g \n",N_OP,mflops);
  if (ckprt > 0) {
    printf(" CheckPrint - ckprt= %d maxnprt= %d \n", ckprt, maxnprt);
    nprt = 0;
    for(i=0; i<SZ; i+= ckprt) {
      printf("i=%d x[i]=%g y[i]=%g z[i]=%g \n", i, (double)x[i], 
	(double)y[i], (double)z[i]);
      nprt++; 
      if (nprt > 1000) break;
    }
  }

  PrtTim("----Fin cpupower");
  printf(":::::: FIN cpupower N_OP= %g MFLOPS= %g ::::::: \n", N_OP,mflops);

  free(x);
  free(y);
  free(z);
  return 0;
}


static clock_t CPUT0, CPUT;
static time_t ELT0, ELT;
static double _totcpu, _partialcpu;

void InitTim(void)
{
CPUT0 = CPUT = clock();
ELT0 = ELT = time(NULL);
_totcpu = _partialcpu = 0.;
return;
}
double GetPartialCPUTime() { return _partialcpu; }
double GetTotalCPUTime() { return _totcpu; }

/* Nouvelle-Fonction */
void PrtTim(const char * Comm)
{
float tcal,tcalt;
clock_t cput;
time_t elt;
int etm,etmt;
 
elt = time(NULL);  cput = clock();
tcalt = ( (float)(cput) - (float)(CPUT0) ) / (float)(CLOCKS_PER_SEC);
tcal = ( (float)(cput) - (float)(CPUT) ) / (float)(CLOCKS_PER_SEC);
_totcpu = tcalt;
_partialcpu = tcal;

etm = elt - ELT;
etmt = elt - ELT0;
printf("%s CPUTime: Total= %g  (Partial= %g) Sec. \n",
       Comm, tcalt, tcal);
printf("ElapsedTime(hh:mm:ss): Total= %02d:%02d:%02d ",
       etmt/3600, (etmt%3600)/60, etmt%60);
printf(" (Partial= %02d:%02d:%02d)\n",
       etm/3600, (etm%3600)/60, etm%60);
 
ELT = elt;
CPUT = cput;
return;
}