CbmRoot
ThermalModelEVMF.h
Go to the documentation of this file.
1 //#pragma once
2 #ifndef THERMALMODELEVMF_H
3 #define THERMALMODELEVMF_H
4 #include "ThermalModelBase.h"
5 #include <cmath>
6 
7 #ifndef PI
8 #define PI 3.141592653589793
9 #endif
10 
11 
12 class ThermalModelEVMF : public ThermalModelBase {
13 public:
14  std::vector<double> densitiesid;
15  double fSuppression;
16  double fDensity;
17  double RHad;
18  int fMode;
19  ThermalModelEVMF(ThermalParticleSystem* TPS_,
20  double T = 0.16,
21  double muB = 0.23,
22  double muS = 0.1,
23  double muQ = -0.01,
24  double gammaS = 1.,
25  double V = 4000.,
26  double R = 1.,
27  int mode = 0)
28  : ThermalModelBase(TPS_, T, muB, muS, muQ, gammaS, V, R)
29  , densitiesid()
30  , fSuppression(1.)
31  , fDensity(0.)
32  , RHad(R)
33  , fMode(mode) {
34  densities.resize(TPS->fParticles.size());
35  densitiesid.resize(TPS->fParticles.size());
36  densitiestotal.resize(TPS->fParticles.size());
37  }
38  ThermalModelEVMF(ThermalParticleSystem* TPS_,
39  const ThermalModelParameters& params,
40  double RHad_ = 0.,
41  int mode = 0)
42  : ThermalModelBase(TPS_, params)
43  , densitiesid()
44  , fSuppression(1.)
45  , fDensity(0.)
46  , RHad(RHad_)
47  , fMode(mode) {
48  densities.resize(TPS->fParticles.size());
49  densitiesid.resize(TPS->fParticles.size());
50  densitiestotal.resize(TPS->fParticles.size());
51  if (TPS->fParticles.size() > 0) RHad = TPS->fParticles[0].fRadius;
52  fVolume = 4000.;
53  }
54  virtual ~ThermalModelEVMF(void);
55  virtual void SetParameters(double T,
56  double muB,
57  double muS,
58  double muQ,
59  double gammaS,
60  double V,
61  double R) {
62  Parameters.T = T;
63  Parameters.muB = muB;
64  Parameters.muS = muS;
65  Parameters.muQ = muQ;
66  Parameters.gammaS = gammaS;
67  Parameters.R = R;
68  Parameters.V = V;
69  fCalculated = false;
70  RHad = R;
71  fVolume = V;
72  }
73  void setRadius(double rad) {
74  Parameters.R = rad;
75  RHad = rad;
76  }
77  void setMode(int mode) { fMode = mode; }
78  virtual void SetParameters(const ThermalModelParameters& params) {
79  Parameters = params;
80  fCalculated = false;
81  RHad = params.R;
82  }
83  virtual void ChangeTPS(ThermalParticleSystem* TPS_) {
84  ThermalModelBase::ChangeTPS(TPS_);
85  densitiesid.resize(TPS->fParticles.size());
86  densities.resize(TPS->fParticles.size());
87  densitiestotal.resize(TPS->fParticles.size());
88  }
89  virtual void FixParameters();
90  virtual void FixParameters(double QB); // And zero net strangeness
91  virtual void CalculateDensities();
92  double Density(double n);
93  double UVdW(double n, double T, double vo) {
94  // van der Waals
95  if (fMode == 0) return T * vo * n / (1. - vo * n) - T * log(1. - vo * n);
96 
97  // Directly excluded
98  if (fMode == 1) return -T * log(1. - 2. * vo * n);
99  //return T * 2. * vo * n;
100 
101  // Carnahan-Starling
102  if (fMode == 2) {
103  double tmp = vo * n - 4.;
104  return -T * (3. + 16. * (-vo * n + 12.) / tmp / tmp / tmp);
105  }
106 
107  // Pade approximation for hard spheres
108  double x = vo * n;
109  double b1 = 0.063597, b2 = 0.017329, b3 = 0.561493, b4 = 0.081313;
110  double u = 1 + b1 * x + b2 * x * x;
111  double d = 1 - b3 * x + b4 * x * x;
112  double c1 = b2 / b4, c2 = (b1 * b4 + b2 * b3) / 2. / b4 / b4,
113  c3 =
114  1. + (b1 * b3 * b4 + b2 * b3 * b3 - 2. * b2 * b4) / 2. / b4 / b4,
115  c4 = 1. / sqrt(4. * b4 - b3 * b3);
116  return T
117  * (x * (c1 + u / d) + c2 * log(d)
118  + 2. * c3 * c4 * (atan(c4 * (2. * b4 * x - b3)) + atan(c4 * b3)));
119 
120  // Virial expansion for hard spheres
121  double b[8];
122  b[2] = 4 * vo, b[3] = 10. * vo * vo, b[4] = 18.365 * vo * vo * vo,
123  b[5] = 28.24 * vo * vo * vo * vo, b[6] = 39.5 * vo * vo * vo * vo * vo,
124  b[7] = 56.5 * vo * vo * vo * vo * vo * vo;
125  double ret = 0.;
126  for (int i = 7; i >= 2; --i) {
127  ret = ret * n + i * b[i] / (i - 1.) * n;
128  }
129  return ret * T;
130  }
131  double PVdW(double n, double T, double vo) {
132  // van der Waals
133  if (fMode == 0) return T * n * vo * n / (1. - vo * n);
134 
135  // Directly excluded
136  if (fMode == 1) return -T * (n + log(1. - 2. * vo * n) / vo / 2.);
137  //return T * vo * n * n;
138 
139  // Carnahan-Starling
140  if (fMode == 2) {
141  double tmp = (1 - vo * n / 4.);
142  return T * n * (vo * n - (vo * n * vo * n) / 8.) / tmp / tmp / tmp;
143  }
144 
145  // Pade approximation
146  double x = vo * n;
147  double b1 = 0.063597, b2 = 0.017329, b3 = 0.561493, b4 = 0.081313;
148  double u = 1 + b1 * x + b2 * x * x;
149  double d = 1 - b3 * x + b4 * x * x;
150  return T * x * x / vo * u / d;
151 
152  // Virial expansion
153  double b[8];
154  b[2] = 4 * vo, b[3] = 10. * vo * vo, b[4] = 18.365 * vo * vo * vo,
155  b[5] = 28.24 * vo * vo * vo * vo, b[6] = 39.5 * vo * vo * vo * vo * vo,
156  b[7] = 56.5 * vo * vo * vo * vo * vo * vo;
157  double ret = 0.;
158  for (int i = 7; i >= 2; --i) {
159  ret = ret * n + b[i] * n;
160  }
161  return ret * T * n;
162  }
163 
164  double GetEffectiveVO(double n, double T, double vo) {
165  //return 1./n * (1. - (1.-vo*n)*exp(-vo*n/(1.-vo*n)));
166  return 1. / n * (1. - exp(-UVdW(n, T, vo) / T));
167  }
168 
169  virtual double GetParticlePrimordialDensity(int part) {
170  if (!fCalculated) CalculateDensities();
171  if (part >= static_cast<int>(densities.size())) return 0.;
172  return densities[part];
173  }
174  virtual double GetParticleTotalDensity(int part) {
175  if (!fCalculated) CalculateDensities();
176  if (part >= static_cast<int>(densitiestotal.size())) return 0.;
177  return densitiestotal[part];
178  }
179 
180  virtual double CalculateHadronDensity() {
181  if (!fCalculated) CalculateDensities();
182  double ret = 0.;
183  for (unsigned int i = 0; i < TPS->fParticles.size(); ++i)
184  ret += densities[i];
185  if (fUseHagedorn) ret += fHagedornDensity;
186  return ret;
187  }
188  virtual double CalculateBaryonDensity() {
189  if (!fCalculated) CalculateDensities();
190  double ret = 0.;
191  for (unsigned int i = 0; i < TPS->fParticles.size(); ++i)
192  ret += TPS->fParticles[i].fB * densities[i];
193  return ret;
194  }
195  virtual double CalculateChargeDensity() {
196  if (!fCalculated) CalculateDensities();
197  double ret = 0.;
198  for (unsigned int i = 0; i < TPS->fParticles.size(); ++i)
199  ret += TPS->fParticles[i].fC * densities[i];
200  return ret;
201  }
202  virtual double CalculateStrangenessDensity() {
203  if (!fCalculated) CalculateDensities();
204  double ret = 0.;
205  for (unsigned int i = 0; i < TPS->fParticles.size(); ++i)
206  ret += TPS->fParticles[i].fS * densities[i];
207  return ret;
208  }
209  virtual double CalculateCharmDensity() {
210  if (!fCalculated) CalculateDensities();
211  double ret = 0.;
212  for (unsigned int i = 0; i < TPS->fParticles.size(); ++i)
213  ret += TPS->fParticles[i].fCharm * densities[i];
214  return ret;
215  }
216  virtual double CalculateAbsoluteStrangenessDensity() {
217  if (!fCalculated) CalculateDensities();
218  double ret = 0.;
219  for (unsigned int i = 0; i < TPS->fParticles.size(); ++i)
220  ret += TPS->fParticles[i].fAbsS * densities[i];
221  return ret;
222  }
223  virtual double CalculateAbsoluteCharmDensity() {
224  if (!fCalculated) CalculateDensities();
225  double ret = 0.;
226  for (unsigned int i = 0; i < TPS->fParticles.size(); ++i)
227  ret += TPS->fParticles[i].fAbsC * densities[i];
228  return ret;
229  }
230  virtual double CalculateEnergyDensity() {
231  double ret = 0.;
232  double vo = 0.;
233  vo = 4. * 4. * PI / 3. * RHad * RHad * RHad;
234  double dMu = -UVdW(fDensity, Parameters.T, vo);
235  for (unsigned int i = 0; i < TPS->fParticles.size(); ++i) {
236  ret += TPS->fParticles[i].CalculateDensity(Parameters.T,
237  Parameters.muB,
238  Parameters.muS,
239  Parameters.muQ,
240  Parameters.gammaS,
241  1,
242  fUseWidth,
243  dMu);
244  }
245  if (fUseHagedorn)
246  ret += fHag.CalculateDensity(Parameters.T, Parameters.muB, 1, dMu);
247  return ret; // / fSuppression;
248  }
249  virtual double CalculateEntropyDensity() {
250  return (1. / Parameters.T)
251  * (CalculateEnergyDensity() + CalculatePressure()
252  - Parameters.muB * CalculateBaryonDensity()
253  - Parameters.muS * CalculateStrangenessDensity()
254  - Parameters.muQ * CalculateChargeDensity());
255  }
256 
257  // Dummy
258  virtual double CalculateBaryonMatterEntropyDensity() {
259  double ret = 0.;
260  return ret;
261  }
262  virtual double CalculateMesonMatterEntropyDensity() {
263  double ret = 0.;
264  return ret;
265  }
266 
267  virtual double CalculatePressure() {
268  if (!fCalculated) CalculateDensities();
269  double ret = 0.;
270  double vo = 0.;
271  vo = 4. * 4. * PI / 3. * RHad * RHad * RHad;
272  double dMu = -UVdW(fDensity, Parameters.T, vo);
273  for (unsigned int i = 0; i < TPS->fParticles.size(); ++i) {
274  ret += TPS->fParticles[i].CalculateDensity(Parameters.T,
275  Parameters.muB,
276  Parameters.muS,
277  Parameters.muQ,
278  Parameters.gammaS,
279  3,
280  fUseWidth,
281  dMu);
282  }
283  if (fUseHagedorn)
284  ret += fHag.CalculateDensity(Parameters.T, Parameters.muB, 3, dMu);
285  return ret + PVdW(fDensity, Parameters.T, vo);
286  }
287 
288  virtual double CalculateShearViscosity();
289 };
290 
291 #endif
ThermalModelParameters::muS
double muS
Definition: ThermalModelBase.h:7
exp
friend F32vec4 exp(const F32vec4 &a)
Definition: L1/vectors/P4_F32vec4.h:134
ThermalModelEVMF::setMode
void setMode(int mode)
Definition: ThermalModelEVMF.h:77
ThermalParticleSystem::fParticles
std::vector< ThermalParticle > fParticles
Definition: ThermalParticleSystem.h:9
ThermalModelBase::densitiestotal
std::vector< double > densitiestotal
Definition: ThermalModelBase.h:33
ThermalModelEVMF::GetParticlePrimordialDensity
virtual double GetParticlePrimordialDensity(int part)
Definition: ThermalModelEVMF.h:169
ThermalModelEVMF::CalculateShearViscosity
virtual double CalculateShearViscosity()
Definition: ThermalModelEVMF.cxx:228
ThermalModelEVMF::fDensity
double fDensity
Definition: ThermalModelEVMF.h:16
sqrt
friend F32vec4 sqrt(const F32vec4 &a)
Definition: L1/vectors/P4_F32vec4.h:41
ThermalModelEVMF::CalculateBaryonMatterEntropyDensity
virtual double CalculateBaryonMatterEntropyDensity()
Definition: ThermalModelEVMF.h:258
ThermalModelBase.h
ThermalModelEVMF::PVdW
double PVdW(double n, double T, double vo)
Definition: ThermalModelEVMF.h:131
ThermalModelParameters::V
double V
Definition: ThermalModelBase.h:7
i
int i
Definition: L1/vectors/P4_F32vec4.h:25
ThermalModelEVMF::GetEffectiveVO
double GetEffectiveVO(double n, double T, double vo)
Definition: ThermalModelEVMF.h:164
HagedornSpectrum::CalculateDensity
double CalculateDensity(double T, double muB, int type=0, double dMu=0.) const
Definition: HagedornSpectrum.cxx:27
ThermalModelEVMF::RHad
double RHad
Definition: ThermalModelEVMF.h:17
ThermalModelEVMF::ThermalModelEVMF
ThermalModelEVMF(ThermalParticleSystem *TPS_, double T=0.16, double muB=0.23, double muS=0.1, double muQ=-0.01, double gammaS=1., double V=4000., double R=1., int mode=0)
Definition: ThermalModelEVMF.h:19
ThermalModelEVMF::fSuppression
double fSuppression
Definition: ThermalModelEVMF.h:15
ThermalModelParameters::gammaS
double gammaS
Definition: ThermalModelBase.h:7
ThermalModelEVMF::CalculatePressure
virtual double CalculatePressure()
Definition: ThermalModelEVMF.h:267
ThermalModelBase::fVolume
double fVolume
Definition: ThermalModelBase.h:28
ThermalModelBase
Definition: ThermalModelBase.h:18
ThermalModelEVMF::Density
double Density(double n)
Definition: ThermalModelEVMF.cxx:133
ThermalModelEVMF::~ThermalModelEVMF
virtual ~ThermalModelEVMF(void)
Definition: ThermalModelEVMF.cxx:98
ThermalModelEVMF::setRadius
void setRadius(double rad)
Definition: ThermalModelEVMF.h:73
ThermalModelBase::ChangeTPS
virtual void ChangeTPS(ThermalParticleSystem *TPS_)
Definition: ThermalModelBase.h:125
PI
#define PI
Definition: HagedornSpectrum.h:5
ThermalModelEVMF::CalculateAbsoluteStrangenessDensity
virtual double CalculateAbsoluteStrangenessDensity()
Definition: ThermalModelEVMF.h:216
ThermalModelEVMF::densitiesid
std::vector< double > densitiesid
Definition: ThermalModelEVMF.h:14
ThermalModelBase::TPS
ThermalParticleSystem * TPS
Definition: ThermalModelBase.h:20
d
double d
Definition: P4_F64vec2.h:24
ThermalModelEVMF::CalculateChargeDensity
virtual double CalculateChargeDensity()
Definition: ThermalModelEVMF.h:195
ThermalModelBase::fCalculated
bool fCalculated
Definition: ThermalModelBase.h:23
ThermalModelParameters::muB
double muB
Definition: ThermalModelBase.h:7
log
friend F32vec4 log(const F32vec4 &a)
Definition: L1/vectors/P4_F32vec4.h:135
ThermalModelBase::fUseWidth
bool fUseWidth
Definition: ThermalModelBase.h:22
ThermalModelParameters::T
double T
Definition: ThermalModelBase.h:7
ThermalModelEVMF::fMode
int fMode
Definition: ThermalModelEVMF.h:18
ThermalModelEVMF::CalculateBaryonDensity
virtual double CalculateBaryonDensity()
Definition: ThermalModelEVMF.h:188
ThermalModelBase::fHagedornDensity
double fHagedornDensity
Definition: ThermalModelBase.h:36
ThermalModelEVMF::CalculateMesonMatterEntropyDensity
virtual double CalculateMesonMatterEntropyDensity()
Definition: ThermalModelEVMF.h:262
ThermalModelEVMF::CalculateAbsoluteCharmDensity
virtual double CalculateAbsoluteCharmDensity()
Definition: ThermalModelEVMF.h:223
ThermalModelEVMF::UVdW
double UVdW(double n, double T, double vo)
Definition: ThermalModelEVMF.h:93
ThermalModelEVMF::FixParameters
virtual void FixParameters()
Definition: ThermalModelEVMF.cxx:100
ThermalModelEVMF::ThermalModelEVMF
ThermalModelEVMF(ThermalParticleSystem *TPS_, const ThermalModelParameters &params, double RHad_=0., int mode=0)
Definition: ThermalModelEVMF.h:38
ThermalModelEVMF::CalculateEntropyDensity
virtual double CalculateEntropyDensity()
Definition: ThermalModelEVMF.h:249
ThermalParticleSystem
Definition: ThermalParticleSystem.h:7
ThermalModelEVMF::SetParameters
virtual void SetParameters(const ThermalModelParameters &params)
Definition: ThermalModelEVMF.h:78
x
Double_t x
Definition: CbmMvdSensorDigiToHitTask.cxx:68
ThermalModelBase::fUseHagedorn
bool fUseHagedorn
Definition: ThermalModelBase.h:26
ThermalModelParameters
Definition: ThermalModelBase.h:6
ThermalModelBase::fHag
HagedornSpectrum fHag
Definition: ThermalModelBase.h:35
ThermalModelParameters::R
double R
Definition: ThermalModelBase.h:7
ThermalModelEVMF::CalculateDensities
virtual void CalculateDensities()
Definition: ThermalModelEVMF.cxx:159
ThermalModelParameters::muQ
double muQ
Definition: ThermalModelBase.h:7
ThermalModelEVMF::CalculateEnergyDensity
virtual double CalculateEnergyDensity()
Definition: ThermalModelEVMF.h:230
ThermalModelBase::densities
std::vector< double > densities
Definition: ThermalModelBase.h:32
ThermalModelEVMF::ChangeTPS
virtual void ChangeTPS(ThermalParticleSystem *TPS_)
Definition: ThermalModelEVMF.h:83
ThermalModelBase::Parameters
ThermalModelParameters Parameters
Definition: ThermalModelBase.h:21
ThermalModelEVMF::CalculateHadronDensity
virtual double CalculateHadronDensity()
Definition: ThermalModelEVMF.h:180
ThermalModelEVMF::CalculateCharmDensity
virtual double CalculateCharmDensity()
Definition: ThermalModelEVMF.h:209
ThermalModelEVMF::GetParticleTotalDensity
virtual double GetParticleTotalDensity(int part)
Definition: ThermalModelEVMF.h:174
ThermalModelEVMF::SetParameters
virtual void SetParameters(double T, double muB, double muS, double muQ, double gammaS, double V, double R)
Definition: ThermalModelEVMF.h:55
ThermalModelEVMF
Definition: ThermalModelEVMF.h:12
ThermalModelEVMF::CalculateStrangenessDensity
virtual double CalculateStrangenessDensity()
Definition: ThermalModelEVMF.h:202