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