TTherminator/Therminator/Hypersurface.cxx

   1 #include "Hypersurface.h"
   2
   3 Hypersurface::Hypersurface(void) {
   4 // Read FOHSI.txt file and init parameters
   5   FName = "fohsi.txt";
   6   HSFile = new ifstream;
   7   HSFile->open(FName);
   8   if(HSFile->is_open()) {
   9     HSFile->seekg (0, ios::beg);
  10 // phi
  11     HSFile->getline(buff,100);
  12     Np=atoi(buff);
  13     HSFile->getline(buff,100);
  14     ip=atof(buff);
  15     HSFile->getline(buff,100);
  16     fp=atof(buff);
  17 // zeta
  18     HSFile->getline(buff,100);
  19     Nz=atoi(buff);
  20     HSFile->getline(buff,100);
  21     iz=atof(buff);
  22     HSFile->getline(buff,100);
  23     fz=atof(buff);
  24 // freeze-out temperature
  25     HSFile->getline(buff,100);
  26     TFO=atof(buff);
  27 // hydro initial time
  28     HSFile->getline(buff,100);
  29     tau0=atof(buff);
  30
  31     HSFile->close();
  32     dp = (fp-ip)/(Np-1);
  33     dz = (fz-iz)/(Nz-1);
  34
  35 // create 2D array
  36     aArr   = new double* [Np];
  37     vArr   = new double* [Np];
  38     dArr   = new double* [Np];
  39     DpdArr = new double* [Np];
  40     DzdArr = new double* [Np];
  41     for(i=0;i<Np;i++) {
  42       aArr[i]   = new double [Nz];
  43       vArr[i]   = new double [Nz];
  44       dArr[i]   = new double [Nz];
  45       DpdArr[i] = new double [Nz];
  46       DzdArr[i] = new double [Nz];
  47     }
  48   } else {
  49     Np = 0; ip = 0.0; fp = 0.0; dp = 1.0;
  50     Nz = 0; iz = 0.0; fz = 0.0; dz = 1.0;
  51     TFO = 0.0;
  52     aArr   = NULL;
  53     vArr   = NULL;
  54     dArr   = NULL;
  55     DpdArr = NULL;
  56     DzdArr = NULL;
  57   }
  58   delete HSFile;
  59
  60   // Read FOHSa.txt file and fill array
  61   FName = "FOHSa.txt";
  62   HSFile = new ifstream;
  63   HSFile->open(FName);
  64   if(HSFile->is_open()) {
  65     HSFile->seekg (0, ios::beg);
  66     for(i=0;i<Np;i++) {
  67       for(j=0;j<Nz;j++) {
  68         HSFile->getline(buff,100);
  69         aArr[i][j]=atof(buff);
  70       }
  71     }
  72     HSFile->close();
  73   }
  74   delete HSFile;
  75
  76   // Read FOHSv.txt file and fill array
  77   FName = "FOHSv.txt";
  78   HSFile = new ifstream;
  79   HSFile->open(FName);
  80   if(HSFile->is_open()) {
  81     HSFile->seekg (0, ios::beg);
  82     for(i=0;i<Np;i++) {
  83       for(j=0;j<Nz;j++) {
  84         HSFile->getline(buff,100);
  85         vArr[i][j]=atof(buff);
  86       }
  87     }
  88     HSFile->close();
  89   }
  90   delete HSFile;
  91
  92   // Read FOHSd.txt file and fill array
  93   FName = "FOHSd.txt";
  94   HSFile = new ifstream;
  95   HSFile->open(FName);
  96   if(HSFile->is_open()) {
  97     HSFile->seekg (0, ios::beg);
  98     for(i=0;i<Np;i++) {
  99       for(j=0;j<Nz;j++) {
 100         HSFile->getline(buff,100);
 101         dArr[i][j]=atof(buff);
 102       }
 103     }
 104     HSFile->close();
 105   }
 106   delete HSFile;
 107
 108   // Read FOHSDpd.txt file and fill array
 109   FName = "FOHSDpd.txt";
 110   HSFile = new ifstream;
 111   HSFile->open(FName);
 112   if(HSFile->is_open()) {
 113     HSFile->seekg (0, ios::beg);
 114     for(i=0;i<Np;i++) {
 115       for(j=0;j<Nz;j++) {
 116         HSFile->getline(buff,100);
 117         DpdArr[i][j]=atof(buff);
 118       }
 119     }
 120     HSFile->close();
 121   }
 122   delete HSFile;
 123
 124   // Read FOHSDzd.txt file and fill array
 125   FName = "FOHSDzd.txt";
 126   HSFile = new ifstream;
 127   HSFile->open(FName);
 128   if(HSFile->is_open()) {
 129     HSFile->seekg (0, ios::beg);
 130     for(i=0;i<Np;i++) {
 131       for(j=0;j<Nz;j++) {
 132         HSFile->getline(buff,100);
 133         DzdArr[i][j]=atof(buff);
 134       }
 135     }
 136     HSFile->close();
 137   }
 138   delete HSFile;
 139 }
 140
 141 Hypersurface::~Hypersurface(void) {
 142   for(i=0;i<Np;i++) {
 143     delete[] aArr[i];
 144     delete[] vArr[i];
 145     delete[] dArr[i];
 146     delete[] DpdArr[i];
 147     delete[] DzdArr[i];
 148   }
 149   delete[] aArr;
 150   delete[] vArr;
 151   delete[] dArr;
 152   delete[] DpdArr;
 153   delete[] DzdArr;
 154 }
 155
 156 double Hypersurface::fahs(double p, double z) {
 157   i=int(p/dp);
 158   j=int(z/dz);
 159   if(i>=Np-1) i=Np-2;
 160   if(j>=Nz-1) j=Nz-2;
 161   return (aArr[i][j]   * (i+1-p/dp) + aArr[i+1][j]   * (p/dp-i)) * (j+1-z/dz) +
 162          (aArr[i][j+1] * (i+1-p/dp) + aArr[i+1][j+1] * (p/dp-i)) * (z/dz-j);
 163 }
 164
 165 double Hypersurface::fvhs(double p, double z) {
 166   i=int(p/dp);
 167   j=int(z/dz);
 168   if(i>=Np-1) i=Np-2;
 169   if(j>=Nz-1) j=Nz-2;
 170   return (vArr[i][j]   * (i+1-p/dp) + vArr[i+1][j]   * (p/dp-i)) * (j+1-z/dz) +
 171          (vArr[i][j+1] * (i+1-p/dp) + vArr[i+1][j+1] * (p/dp-i)) * (z/dz-j);
 172 }
 173
 174 double Hypersurface::fdhs(double p, double z) {
 175   i=int(p/dp);
 176   j=int(z/dz);
 177   if(i>=Np-1) i=Np-2;
 178   if(j>=Nz-1) j=Nz-2;
 179   return (dArr[i][j]   * (i+1-p/dp) + dArr[i+1][j]   * (p/dp-i)) * (j+1-z/dz) +
 180          (dArr[i][j+1] * (i+1-p/dp) + dArr[i+1][j+1] * (p/dp-i)) * (z/dz-j);
 181 }
 182
 183 double Hypersurface::fDpdhs(double p, double z) {
 184   i=int(p/dp);
 185   j=int(z/dz);
 186   if(i>=Np-1) i=Np-2;
 187   if(j>=Nz-1) j=Nz-2;
 188   return (DpdArr[i][j]   * (i+1-p/dp) + DpdArr[i+1][j]   * (p/dp-i)) * (j+1-z/dz) +
 189          (DpdArr[i][j+1] * (i+1-p/dp) + DpdArr[i+1][j+1] * (p/dp-i)) * (z/dz-j);
 190 }
 191
 192 double Hypersurface::fDzdhs(double p, double z) {
 193   i=int(p/dp);
 194   j=int(z/dz);
 195   if(i>=Np-1) i=Np-2;
 196   if(j>=Nz-1) j=Nz-2;
 197   return (DzdArr[i][j]   * (i+1-p/dp) + DzdArr[i+1][j]   * (p/dp-i)) * (j+1-z/dz) +
 198          (DzdArr[i][j+1] * (i+1-p/dp) + DzdArr[i+1][j+1] * (p/dp-i)) * (z/dz-j);
 199 }