Thermal-FIST  1.3
Package for hadron resonance gas model applications
Classes | Public Member Functions | Protected Member Functions | Protected Attributes | List of all members
thermalfist::ThermalModelEVDiagonal Class Reference

Class implementing the diagonal excluded-volume model. More...

#include <ThermalModelEVDiagonal.h>

Inheritance diagram for thermalfist::ThermalModelEVDiagonal:
thermalfist::ThermalModelBase

Public Member Functions

 ThermalModelEVDiagonal (ThermalParticleSystem *TPS, const ThermalModelParameters &params=ThermalModelParameters())
 Construct a new ThermalModelEVDiagonal object. More...
 
virtual ~ThermalModelEVDiagonal (void)
 Destroy the ThermalModelEVDiagonal object. More...
 
void FillVirialEV (const std::vector< double > &vi=std::vector< double >(0))
 Same as FillVirial() but uses the diagonal excluded-volume coefficients \( v_i \equiv b_{ii} \) as input instead of radii. More...
 
double ExcludedVolume (int i) const
 The excluded volume parameter of particle species i. More...
 
double CommonSuppressionFactor ()
 The density suppression factor \( (1 - \sum_i v_i n_i) \), common for all species. More...
 
void SetRadius (double rad)
 Set the same excluded volume radius parameter for all species. More...
 
void SetRadius (int i, double rad)
 Set the radius parameter for particle species i. More...
 
void FillVirial (const std::vector< double > &ri=std::vector< double >(0))
 Fills the vector of particle eigenvolume parameters. More...
 
virtual void ReadInteractionParameters (const std::string &filename)
 Read the set of eigenvolume parameters for all particles from an external file. More...
 
virtual void WriteInteractionParameters (const std::string &filename)
 Write the set of eigenvolume parameters for all particles to an external file. More...
 
virtual double CalculateEigenvolumeFraction ()
 
double VirialCoefficient (int i, int j) const
 Return the eigenvolume parameter of particle species i. More...
 
void SetVirial (int i, int j, double b)
 Sets the eigenvolume parameter of particle species i. More...
 
virtual void ChangeTPS (ThermalParticleSystem *TPS)
 Change the particle list. More...
 
virtual void CalculatePrimordialDensities ()
 Calculates the primordial densities of all species. More...
 
virtual void CalculateFluctuations ()
 Computes the fluctuation observables. More...
 
virtual void CalculateTwoParticleCorrelations ()
 Computes the fluctuations and correlations of the primordial particle numbers. More...
 
virtual std::vector< double > CalculateChargeFluctuations (const std::vector< double > &chgs, int order=4)
 Calculates fluctuations (diagonal susceptibilities) of an arbitrary "conserved" charge. More...
 
virtual double CalculatePressure ()
 Implementation of the equation of state functions. More...
 
virtual double CalculateEnergyDensity ()
 
virtual double CalculateEntropyDensity ()
 
virtual double CalculateBaryonMatterEntropyDensity ()
 The fraction of entropy carried by baryons (Ideal GCE only) More...
 
virtual double CalculateMesonMatterEntropyDensity ()
 The fraction of entropy carried by mesons (Ideal GCE only) More...
 
virtual double ParticleScaledVariance (int)
 Scaled variance of primordial particle number fluctuations for species i. More...
 
virtual double ParticleSkewness (int)
 Skewness of primordial particle number fluctuations for species i. More...
 
virtual double ParticleKurtosis (int)
 Kurtosis of primordial particle number fluctuations for species i. More...
 
virtual double ParticleScalarDensity (int part)
 The scalar density of the particle species i. More...
 
virtual void DisableMesonMesonVirial ()
 
virtual void DisableMesonMesonAttraction ()
 
virtual void DisableMesonBaryonVirial ()
 
virtual void DisableMesonBaryonAttraction ()
 
virtual void DisableBaryonBaryonVirial ()
 
virtual void DisableBaryonBaryonAttraction ()
 
virtual void DisableBaryonAntiBaryonVirial ()
 
virtual void DisableBaryonAntiBaryonAttraction ()
 
- Public Member Functions inherited from thermalfist::ThermalModelBase
 ThermalModelBase (ThermalParticleSystem *TPS, const ThermalModelParameters &params=ThermalModelParameters())
 Construct a new ThermalModelBase object. More...
 
virtual ~ThermalModelBase (void)
 
int ComponentsNumber () const
 Number of different particle species in the list. More...
 
bool UseWidth () const
 Whether finite resonance widths are considered. More...
 
void SetUseWidth (bool useWidth)
 Sets whether finite resonance widths are used. Deprecated. More...
 
void SetUseWidth (ThermalParticle::ResonanceWidthIntegration type)
 Sets the finite resonance widths scheme to use. More...
 
void SetOMP (bool openMP)
 OpenMP support. Currently not used. More...
 
void UpdateParameters ()
 Calls SetParameters() with current m_Parameters. More...
 
virtual void SetTemperature (double T)
 Set the temperature. More...
 
virtual void SetBaryonChemicalPotential (double muB)
 Set the baryon chemical potential. More...
 
virtual void SetElectricChemicalPotential (double muQ)
 Set the electric chemical potential. More...
 
virtual void SetStrangenessChemicalPotential (double muS)
 Set the strangeness chemical potential. More...
 
virtual void SetCharmChemicalPotential (double muC)
 Set the charm chemical potential. More...
 
virtual void SetGammaq (double gammaq)
 Set the light quark fugacity factor. More...
 
virtual void SetGammaS (double gammaS)
 Set the strange quark fugacity factor. More...
 
virtual void SetGammaC (double gammaC)
 Set the charm quark fugacity factor. More...
 
virtual void SetBaryonCharge (int B)
 Set the total baryon number (for canonical ensemble only) More...
 
virtual void SetElectricCharge (int Q)
 Set the total electric charge (for canonical ensemble only) More...
 
virtual void SetStrangeness (int S)
 Set the total strangeness (for canonical ensemble only) More...
 
virtual void SetCharm (int C)
 Set the total charm (for canonical ensemble only) More...
 
virtual void SetAttraction (int, int, double)
 Set the vdW mean field attraction coefficient \( a_{ij} \). More...
 
virtual double AttractionCoefficient (int, int) const
 QvdW mean field attraction coefficient \( a_{ij} \). More...
 
bool QuantumStatistics () const
 
virtual void SetStatistics (bool stats)
 
virtual void SetCalculationType (IdealGasFunctions::QStatsCalculationType type)
 Sets the CalculationType() method to evaluate quantum statistics. Calls the corresponding method in TPS(). More...
 
virtual void SetClusterExpansionOrder (int order)
 Set the number of terms in the cluster expansion method. Calls the corresponding method in TPS(). More...
 
void SetResonanceWidthShape (ThermalParticle::ResonanceWidthShape shape)
 Set the ThermalParticle::ResonanceWidthShape for all particles. Calls the corresponding method in TPS(). More...
 
void SetResonanceWidthIntegrationType (ThermalParticle::ResonanceWidthIntegration type)
 Set the ThermalParticle::ResonanceWidthIntegration scheme for all particles. Calls the corresponding method in TPS(). More...
 
virtual void FillChemicalPotentials ()
 Sets the chemical potentials of all particles. More...
 
virtual void SetChemicalPotentials (const std::vector< double > &chem=std::vector< double >(0))
 Sets the chemical potentials of all particles. More...
 
const std::vector< double > & ChemicalPotentials () const
 A vector of chemical potentials of all particles. More...
 
double ChemicalPotential (int i) const
 Chemical potential of particle species i. More...
 
virtual double FullIdealChemicalPotential (int i) const
 Chemical potential entering the ideal gas expressions of particle species i. More...
 
bool ConstrainMuB () const
 
void ConstrainMuB (bool constrain)
 
bool ConstrainMuQ () const
 
void ConstrainMuQ (bool constrain)
 
bool ConstrainMuS () const
 
void ConstrainMuS (bool constrain)
 
bool ConstrainMuC () const
 
void ConstrainMuC (bool constrain)
 
void UsePartialChemicalEquilibrium (bool usePCE)
 Sets whether partial chemical equilibrium with additional chemical potentials is used. More...
 
bool UsePartialChemicalEquilibrium ()
 Whether partial chemical equilibrium with additional chemical potentials is used. More...
 
void SetVolume (double Volume)
 Sets the system volume. More...
 
double Volume () const
 System volume (fm \(^3\)) More...
 
void SetVolumeRadius (double R)
 Sets the system radius. More...
 
double CanonicalVolume () const
 The canonical correlation volume V \(_c\) (fm \(^3\)) More...
 
void SetCanonicalVolume (double Volume)
 Set the canonical correlation volume V \(_c\). More...
 
void SetCanonicalVolumeRadius (double Rc)
 Set the canonical correlation system radius. More...
 
void ConstrainChemicalPotentials (bool resetInitialValues=true)
 Constrains the chemical potentials \( \mu_B,\,\mu_Q,\,\mu_S,\,\mu_C \) by the conservation laws imposed. More...
 
virtual void FixParameters ()
 Method which actually implements ConstrainChemicalPotentials() (for backward compatibility). More...
 
virtual void FixParametersNoReset ()
 Method which actually implements ConstrainChemicalPotentialsNoReset() (for backward compatibility). More...
 
virtual bool SolveChemicalPotentials (double totB=0., double totQ=0., double totS=0., double totC=0., double muBinit=0., double muQinit=0., double muSinit=0., double muCinit=0., bool ConstrMuB=true, bool ConstrMuQ=true, bool ConstrMuS=true, bool ConstrMuC=true)
 The procedure which calculates the chemical potentials \( \mu_B,\,\mu_Q,\,\mu_S,\,\mu_Q \) which reproduce the specified total baryon, electric, strangeness, and charm charges of the system. More...
 
virtual void CalculateDensities ()
 Calculates the primordial and total (after decays) densities of all species. More...
 
virtual void ValidateCalculation ()
 Checks whether issues have occured during the calculation of particle densities in the CalculateDensities() method. More...
 
std::string ValidityCheckLog () const
 All messaged which occured during the validation procedure in the ValidateCalculation() method. More...
 
virtual void CalculateDensitiesGCE ()
 Calculates the particle densities in a grand-canonical ensemble. More...
 
virtual void CalculateFeeddown ()
 Calculates the total densities which include feeddown contributions. More...
 
virtual void CalculateTwoParticleFluctuationsDecays ()
 Computes particle number correlations and fluctuations for all final-state particles which are marked stable. More...
 
virtual double TwoParticleSusceptibilityPrimordial (int i, int j) const
 Returns the computed primordial particle number (cross-)susceptibility \( \frac{1}{VT^3} \, \langle \Delta N_i \Delta N_j \rangle \) for particles with ids i and j. CalculateFluctuations() must be called beforehand. More...
 
virtual double TwoParticleSusceptibilityPrimordialByPdg (long long id1, long long id2)
 Returns the computed primordial particle number (cross-)susceptibility \( \frac{1}{VT^3} \, \langle \Delta N_i \Delta N_j \rangle \) for particles with pdg codes id1 and id2. CalculateFluctuations() must be called beforehand. More...
 
virtual double NetParticleSusceptibilityPrimordialByPdg (long long id1, long long id2)
 Returns the computed (cross-)susceptibility \( \frac{1}{VT^3} \, \langle \Delta N_i \Delta N_j \rangle \) between primordial net-particle numbers for pdg codes id1 and id2. CalculateFluctuations() must be called beforehand. More...
 
virtual double TwoParticleSusceptibilityFinal (int i, int j) const
 Returns the computed final particle number (cross-)susceptibility \( \frac{1}{VT^3} \, \langle \Delta N_i \Delta N_j \rangle \) for particles with ids i and j. CalculateFluctuations() must be called beforehand. Both particle species must be those marked stable. More...
 
virtual double TwoParticleSusceptibilityFinalByPdg (long long id1, long long id2)
 Returns the computed final particle number (cross-)susceptibility \( \frac{1}{VT^3} \, \langle \Delta N_i \Delta N_j \rangle \) for particles with pdg codes id1 and id2. CalculateFluctuations() must be called beforehand. Both particle species must be those marked stable. More...
 
virtual double NetParticleSusceptibilityFinalByPdg (long long id1, long long id2)
 Returns the computed (cross-)susceptibility \( \frac{1}{VT^3} \, \langle \Delta N_i \Delta N_j \rangle \) between final net-particle numbers for pdg codes id1 and id2. CalculateFluctuations() must be called beforehand. More...
 
virtual double PrimordialParticleChargeSusceptibility (int i, ConservedCharge::Name chg) const
 Returns the computed primordial particle vs conserved charge (cross-)susceptibility \( \frac{1}{VT^3} \, \langle \Delta N_i \Delta N_chg \rangle \) for particle with id i and conserved charge chg. CalculateFluctuations() must be called beforehand. More...
 
virtual double PrimordialParticleChargeSusceptibilityByPdg (long long id1, ConservedCharge::Name chg)
 Returns the computed primordial particle vs conserved charge (cross-)susceptibility \( \frac{1}{VT^3} \, \langle \Delta N_i \Delta N_chg \rangle \) for particle with pdg code id1 and conserved charge chg. CalculateFluctuations() must be called beforehand. More...
 
virtual double PrimordialNetParticleChargeSusceptibilityByPdg (long long id1, ConservedCharge::Name chg)
 Returns the computed primordial net particle vs conserved charge (cross-)susceptibility \( \frac{1}{VT^3} \, \langle \Delta N_i \Delta N_chg \rangle \) for particle with pdg code id1 and conserved charge chg. CalculateFluctuations() must be called beforehand. More...
 
virtual double FinalParticleChargeSusceptibility (int i, ConservedCharge::Name chg) const
 Returns the computed final (after decays) particle vs conserved charge (cross-)susceptibility \( \frac{1}{VT^3} \, \langle \Delta N_i \Delta N_chg \rangle \) for particle with id i and conserved charge chg. CalculateFluctuations() must be called beforehand. More...
 
virtual double FinalParticleChargeSusceptibilityByPdg (long long id1, ConservedCharge::Name chg)
 Returns the computed final (after decays) particle vs conserved charge (cross-)susceptibility \( \frac{1}{VT^3} \, \langle \Delta N_i \Delta N_chg \rangle \) for particle with pdg code id1 and conserved charge chg. CalculateFluctuations() must be called beforehand. More...
 
virtual double FinalNetParticleChargeSusceptibilityByPdg (long long id1, ConservedCharge::Name chg)
 Returns the computed final (after decays) net particle vs conserved charge (cross-)susceptibility \( \frac{1}{VT^3} \, \langle \Delta N_i \Delta N_chg \rangle \) for particle with pdg code id1 and conserved charge chg. CalculateFluctuations() must be called beforehand. More...
 
virtual void CalculateSusceptibilityMatrix ()
 Calculates the conserved charges susceptibility matrix. More...
 
virtual void CalculateProxySusceptibilityMatrix ()
 Calculates the susceptibility matrix of conserved charges proxies. More...
 
virtual void CalculateParticleChargeCorrelationMatrix ()
 Calculates the matrix of correlators between primordial (and also final) particle numbers and conserved charges. More...
 
double Pressure ()
 System pressure (GeV fm \(^{-3}\)) More...
 
double EnergyDensity ()
 Energy density (GeV fm \(^{-3}\)) More...
 
double EntropyDensity ()
 Entropy density (fm \(^{-3}\)) More...
 
double HadronDensity ()
 Total number density of all particle (fm \(^{-3}\)) More...
 
double BaryonDensity ()
 Net baryon density (fm \(^{-3}\)) More...
 
double ElectricChargeDensity ()
 Electric charge density (fm \(^{-3}\)) More...
 
double StrangenessDensity ()
 Net strangeness density (fm \(^{-3}\)) More...
 
double CharmDensity ()
 Net charm density (fm \(^{-3}\)) More...
 
double AbsoluteBaryonDensity ()
 Absolute baryon number density (baryons + antibaryons) (fm \(^{-3}\)) More...
 
double AbsoluteElectricChargeDensity ()
 Absolute electric charge density (Q+ + Q-) (fm \(^{-3}\)) More...
 
double AbsoluteStrangenessDensity ()
 Absolute strange quark content density (fm \(^{-3}\)) More...
 
double AbsoluteCharmDensity ()
 Absolute charm quark content density (fm \(^{-3}\)) More...
 
virtual double CalculateArbitraryChargeDensity ()
 Computes the density of the auxiliary ArbitraryCharge() More...
 
virtual double GetMaxDiff () const
 
virtual bool IsLastSolutionOK () const
 
double GetDensity (long long PDGID, int feeddown)
 Same as GetDensity(int,Feeddown::Type) More...
 
double GetDensity (long long PDGID, Feeddown::Type feeddown)
 Particle number density of species with a specified PDG ID and feeddown. More...
 
double GetYield (long long PDGID, Feeddown::Type feeddown)
 Particle number yield of species with a specified PDG ID and feeddown. More...
 
std::vector< double > GetIdealGasDensities () const
 
ThermalParticleSystemTPS ()
 
const std::vector< double > & Densities () const
 
const std::vector< double > & TotalDensities () const
 
const std::vector< std::vector< double > > & AllDensities () const
 
void SetTAG (const std::string &tag)
 Set the tag for this ThermalModelBase object. More...
 
const std::string & TAG () const
 The tag of this ThermalModelBase object. More...
 
void ResetCalculatedFlags ()
 Reset all flags which correspond to a calculation status. More...
 
bool IsCalculated () const
 
bool IsFluctuationsCalculated () const
 
bool IsGCECalculated () const
 
double ScaledVariancePrimordial (int id) const
 Scaled variance of primordial particle number fluctuations for particle species id. More...
 
double ScaledVarianceTotal (int id) const
 Scaled variance of final particle number fluctuations for particle species id. More...
 
double SkewnessPrimordial (int id) const
 Normalized skewness of primordial particle number fluctuations for particle species id. More...
 
double SkewnessTotal (int id) const
 Normalized skewness of final particle number fluctuations for particle species id. More...
 
double KurtosisPrimordial (int id) const
 Normalized excess kurtosis of primordial particle number fluctuations for particle species id. More...
 
double KurtosisTotal (int id) const
 Normalized excess kurtosis of final particle number fluctuations for particle species id. More...
 
double ConservedChargeDensity (ConservedCharge::Name chg)
 A density of a conserved charge (in fm^-3) More...
 
double Susc (ConservedCharge::Name i, ConservedCharge::Name j) const
 A 2nd order susceptibility of conserved charges. More...
 
double ProxySusc (ConservedCharge::Name i, ConservedCharge::Name j) const
 A 2nd order susceptibility of conserved charges proxies. More...
 
double ChargedMultiplicity (int type=0)
 Multiplicity of charged particles. More...
 
double ChargedScaledVariance (int type=0)
 Scaled variance of charged particles. More...
 
double ChargedMultiplicityFinal (int type=0)
 Multiplicity of charged particles including the feeddown contributions in accordance with the stability flags. More...
 
double ChargedScaledVarianceFinal (int type=0)
 Scaled variance of charged particles including the feeddown contributions in accordance with the stability flags. More...
 
ThermalModelEnsemble Ensemble ()
 The statistical ensemble of the current HRG model. More...
 
virtual bool IsConservedChargeCanonical (ConservedCharge::Name charge) const
 Whether the given conserved charge is treated canonically. More...
 
ThermalModelInteraction InteractionModel ()
 The interactions present in the current HRG model. More...
 
void SetNormBratio (bool normBratio)
 Whether branching ratios are renormalized to 100%. More...
 
bool NormBratio () const
 
virtual void SetParameters (const ThermalModelParameters &params)
 The thermal parameters. More...
 
const ThermalModelParametersParameters () const
 
void DisableMesonMesonRepulsion ()
 
void DisableMesonBaryonRepulsion ()
 
void DisableBaryonBaryonRepulsion ()
 
void DisableBaryonAntiBaryonRepulsion ()
 
double RepulsionCoefficient (int i, int j) const
 
void SetSoverB (double SB)
 The entropy per baryon ratio to be used to constrain the baryon chemical potential. More...
 
double SoverB () const
 
void SetQoverB (double QB)
 The electric-to-baryon charge ratio to be used to constrain the electric chemical potential. More...
 
double QoverB () const
 
virtual double CalculateHadronDensity ()
 
virtual double CalculateBaryonDensity ()
 
virtual double CalculateChargeDensity ()
 
virtual double CalculateStrangenessDensity ()
 
virtual double CalculateCharmDensity ()
 
virtual double CalculateAbsoluteBaryonDensity ()
 
virtual double CalculateAbsoluteChargeDensity ()
 
virtual double CalculateAbsoluteStrangenessDensity ()
 
virtual double CalculateAbsoluteCharmDensity ()
 

Protected Member Functions

virtual void SolvePressure ()
 Solves the transcdental equation for the pressure. More...
 
double Pressure (double P)
 Computes the l.h.s. of the transcdental equation for the pressure. More...
 
double DensityId (int i, double Pressure)
 Calculate the ideal gas density of particle species i for the given pressure value. More...
 
double PressureId (int i, double Pressure)
 Calculate the ideal gas pressure of particle species i for the given pressure value. More...
 
double ScaledVarianceId (int i, double Pressure)
 Calculate the ideal gas scaled variance of particle species i number fluctuations for the given pressure value. More...
 
virtual double MuShift (int i) const
 The shift in the chemical potential of particle species i due to the excluded volume interactions. More...
 

Protected Attributes

std::vector< double > m_densitiesid
 
std::vector< double > m_densitiesidnoshift
 
std::vector< double > m_v
 
double m_Suppression
 
double m_Pressure
 
double m_Densityid
 
double m_TotalDensity
 
EVSolution m_sol
 
- Protected Attributes inherited from thermalfist::ThermalModelBase
ThermalModelParameters m_Parameters
 
ThermalParticleSystemm_TPS
 
bool m_LastCalculationSuccessFlag
 
double m_MaxDiff
 
bool m_Calculated
 
bool m_FeeddownCalculated
 
bool m_FluctuationsCalculated
 
bool m_GCECalculated
 
bool m_UseWidth
 
bool m_NormBratio
 
bool m_QuantumStats
 
double m_QBgoal
 
double m_SBgoal
 
double m_Volume
 
bool m_ConstrainMuB
 
bool m_ConstrainMuQ
 
bool m_ConstrainMuS
 
bool m_ConstrainMuC
 
bool m_PCE
 
bool m_useOpenMP
 
std::vector< double > m_densities
 
std::vector< double > m_densitiestotal
 
std::vector< std::vector< double > > m_densitiesbyfeeddown
 
std::vector< double > m_Chem
 
std::vector< double > m_wprim
 
std::vector< double > m_wtot
 
std::vector< double > m_skewprim
 
std::vector< double > m_skewtot
 
std::vector< double > m_kurtprim
 
std::vector< double > m_kurttot
 
std::vector< std::vector< double > > m_PrimCorrel
 
std::vector< std::vector< double > > m_TotalCorrel
 
std::vector< std::vector< double > > m_PrimChargesCorrel
 
std::vector< std::vector< double > > m_FinalChargesCorrel
 
std::vector< std::vector< double > > m_Susc
 
std::vector< std::vector< double > > m_ProxySusc
 
std::string m_ValidityLog
 
double m_wnSum
 
std::string m_TAG
 
ThermalModelEnsemble m_Ensemble
 
ThermalModelInteraction m_InteractionModel
 

Additional Inherited Members

- Public Types inherited from thermalfist::ThermalModelBase
enum  ThermalModelEnsemble { GCE = 0, CE = 1, SCE = 2, CCE = 3 }
 The list of statistical ensembles. More...
 
enum  ThermalModelInteraction {
  Ideal = 0, DiagonalEV = 1, CrosstermsEV = 2, QvdW = 3,
  RealGas = 4, MeanField = 5
}
 Type of interactions included in the HRG model. More...
 

Detailed Description

Class implementing the diagonal excluded-volume model.

The model formulation can be found in G.D. Yen, M.I. Gorenstein, W. Greiner, S.-N. Yang, Phys. Rev. C 56 2210, (1997), http://arxiv.org/pdf/nucl-th/9711062.pdf

The transcendental equation for the pressure is solved using the Broyden's method.

Examples:
BagModelFit.cpp, cpc1-HRG-TDep.cpp, and cpc2-chi2-vs-T.cpp.

Definition at line 45 of file ThermalModelEVDiagonal.h.

Constructor & Destructor Documentation

thermalfist::ThermalModelEVDiagonal::ThermalModelEVDiagonal ( ThermalParticleSystem TPS,
const ThermalModelParameters params = ThermalModelParameters() 
)

Construct a new ThermalModelEVDiagonal object.

Parameters
TPSA pointer to the ThermalParticleSystem object containing the particle list
paramsThermalModelParameters object with current thermal parameters

Definition at line 32 of file ThermalModelEVDiagonal.cpp.

thermalfist::ThermalModelEVDiagonal::~ThermalModelEVDiagonal ( void  )
virtual

Destroy the ThermalModelEVDiagonal object.

Definition at line 45 of file ThermalModelEVDiagonal.cpp.

Member Function Documentation

virtual double thermalfist::ThermalModelEVDiagonal::CalculateBaryonMatterEntropyDensity ( )
inlinevirtual

The fraction of entropy carried by baryons (Ideal GCE only)

Reimplemented from thermalfist::ThermalModelBase.

Definition at line 176 of file ThermalModelEVDiagonal.h.

std::vector< double > thermalfist::ThermalModelEVDiagonal::CalculateChargeFluctuations ( const std::vector< double > &  chgs,
int  order = 4 
)
virtual

Calculates fluctuations (diagonal susceptibilities) of an arbitrary "conserved" charge.

Each particle specie is assumed to carry a conserved charge with a value provided by an input vector.

Restricted to the grand canonical ensemble.

Parameters
chgsA vector with conserved charge values for all species. 0-based indices of the vector must correspond to the 0-based indices of the particle list TPS()
orderUp to which order the susceptibilities are computed
Returns
std::vector<double> A vector with computed values of diagonal susceptibilities

Reimplemented from thermalfist::ThermalModelBase.

Definition at line 295 of file ThermalModelEVDiagonal.cpp.

double thermalfist::ThermalModelEVDiagonal::CalculateEigenvolumeFraction ( )
virtual

Fraction of the total volume occupied by the finite-sizes particles (Diagonal excluded volume model only)

Reimplemented from thermalfist::ThermalModelBase.

Definition at line 537 of file ThermalModelEVDiagonal.cpp.

double thermalfist::ThermalModelEVDiagonal::CalculateEnergyDensity ( )
virtual

Implements thermalfist::ThermalModelBase.

Definition at line 487 of file ThermalModelEVDiagonal.cpp.

double thermalfist::ThermalModelEVDiagonal::CalculateEntropyDensity ( )
virtual

Implements thermalfist::ThermalModelBase.

Definition at line 498 of file ThermalModelEVDiagonal.cpp.

void thermalfist::ThermalModelEVDiagonal::CalculateFluctuations ( )
virtual

Computes the fluctuation observables.

Includes the matrix of 2nd order susceptibilities of conserved charges, as well as particle number correlations and fluctuations, both for primordial yields and after decays.

Reimplemented from thermalfist::ThermalModelBase.

Definition at line 274 of file ThermalModelEVDiagonal.cpp.

virtual double thermalfist::ThermalModelEVDiagonal::CalculateMesonMatterEntropyDensity ( )
inlinevirtual

The fraction of entropy carried by mesons (Ideal GCE only)

Reimplemented from thermalfist::ThermalModelBase.

Definition at line 178 of file ThermalModelEVDiagonal.h.

double thermalfist::ThermalModelEVDiagonal::CalculatePressure ( )
virtual

Implementation of the equation of state functions.

Implements thermalfist::ThermalModelBase.

Definition at line 509 of file ThermalModelEVDiagonal.cpp.

void thermalfist::ThermalModelEVDiagonal::CalculatePrimordialDensities ( )
virtual

Calculates the primordial densities of all species.

Implements thermalfist::ThermalModelBase.

Definition at line 202 of file ThermalModelEVDiagonal.cpp.

void thermalfist::ThermalModelEVDiagonal::CalculateTwoParticleCorrelations ( )
virtual

Computes the fluctuations and correlations of the primordial particle numbers.

More specifically, computes the susceptibility matrix \( \frac{1}{VT^3} \, \langle \Delta N_i^* \Delta N_j^* \rangle \), where \( N_i^* \) is the primordial yield of species i.

Reimplemented from thermalfist::ThermalModelBase.

Definition at line 242 of file ThermalModelEVDiagonal.cpp.

void thermalfist::ThermalModelEVDiagonal::ChangeTPS ( ThermalParticleSystem TPS)
virtual

Change the particle list.

Parameters
TPSA pointer to new particle list.

Reimplemented from thermalfist::ThermalModelBase.

Definition at line 130 of file ThermalModelEVDiagonal.cpp.

double thermalfist::ThermalModelEVDiagonal::CommonSuppressionFactor ( )

The density suppression factor \( (1 - \sum_i v_i n_i) \), common for all species.

Returns
The density suppression factor

Definition at line 522 of file ThermalModelEVDiagonal.cpp.

double thermalfist::ThermalModelEVDiagonal::DensityId ( int  i,
double  Pressure 
)
protected

Calculate the ideal gas density of particle species i for the given pressure value.

Parameters
i0-based particle specie index
PressureInput pressure (GeV fm \(^{-3}\))
Returns
Ideal gas density (fm \(^{-3}\))

Definition at line 136 of file ThermalModelEVDiagonal.cpp.

virtual void thermalfist::ThermalModelEVDiagonal::DisableBaryonAntiBaryonAttraction ( )
inlinevirtual

Switches off QvdW attraction terms for all baryon-antibaryon pairs i.e. \( a_{ij} = 0 \) for baryon-antibaryon pairs

Reimplemented from thermalfist::ThermalModelBase.

Definition at line 102 of file ThermalModelEVDiagonal.h.

virtual void thermalfist::ThermalModelEVDiagonal::DisableBaryonAntiBaryonVirial ( )
inlinevirtual

Switches off eigenvolume terms for all baryon-antibaryon pairs i.e. \( \tilde{b}_{ij} = 0 \) for baryon-antibaryon pairs

Reimplemented from thermalfist::ThermalModelBase.

Definition at line 101 of file ThermalModelEVDiagonal.h.

virtual void thermalfist::ThermalModelEVDiagonal::DisableBaryonBaryonAttraction ( )
inlinevirtual

Switches off QvdW attraction terms for all baryon-baryon pairs i.e. \( a_{ij} = 0 \) for baryon-baryon pairs

Reimplemented from thermalfist::ThermalModelBase.

Definition at line 100 of file ThermalModelEVDiagonal.h.

virtual void thermalfist::ThermalModelEVDiagonal::DisableBaryonBaryonVirial ( )
inlinevirtual

Switches off excluded volume terms for all baryon-baryon pairs i.e. \( \tilde{b}_{ij} = 0 \) for baryon-baryon pairs

Reimplemented from thermalfist::ThermalModelBase.

Definition at line 99 of file ThermalModelEVDiagonal.h.

virtual void thermalfist::ThermalModelEVDiagonal::DisableMesonBaryonAttraction ( )
inlinevirtual

Switches off QvdW attraction terms for all meson-baryon pairs i.e. \( a_{ij} = 0 \) for meson-baryon pairs

Reimplemented from thermalfist::ThermalModelBase.

Definition at line 98 of file ThermalModelEVDiagonal.h.

virtual void thermalfist::ThermalModelEVDiagonal::DisableMesonBaryonVirial ( )
inlinevirtual

Switches off excluded volume terms for all meson-baryon pairs i.e. \( \tilde{b}_{ij} = 0 \) for meson-baryon pairs

Reimplemented from thermalfist::ThermalModelBase.

Definition at line 97 of file ThermalModelEVDiagonal.h.

virtual void thermalfist::ThermalModelEVDiagonal::DisableMesonMesonAttraction ( )
inlinevirtual

Switches off QvdW attraction terms for all meson-meson pairs i.e. \( a_{ij} = 0 \) for meson-meson pairs

Reimplemented from thermalfist::ThermalModelBase.

Definition at line 96 of file ThermalModelEVDiagonal.h.

virtual void thermalfist::ThermalModelEVDiagonal::DisableMesonMesonVirial ( )
inlinevirtual

Differential treatment of pair interactions is not applicable for the diagonal excluded volume model

Reimplemented from thermalfist::ThermalModelBase.

Definition at line 95 of file ThermalModelEVDiagonal.h.

double thermalfist::ThermalModelEVDiagonal::ExcludedVolume ( int  i) const

The excluded volume parameter of particle species i.

Parameters
i0-based particle species index.
Returns
The excluded volume parameter.

Definition at line 122 of file ThermalModelEVDiagonal.cpp.

void thermalfist::ThermalModelEVDiagonal::FillVirial ( const std::vector< double > &  ri = std::vector<double>(0))
virtual

Fills the vector of particle eigenvolume parameters.

Parameters
riA vector of radii parameters for all particles (in fm)

Reimplemented from thermalfist::ThermalModelBase.

Definition at line 58 of file ThermalModelEVDiagonal.cpp.

void thermalfist::ThermalModelEVDiagonal::FillVirialEV ( const std::vector< double > &  vi = std::vector<double>(0))

Same as FillVirial() but uses the diagonal excluded-volume coefficients \( v_i \equiv b_{ii} \) as input instead of radii.

Parameters
viA vector with diagonal excluded-volume coefficients for all species. 0-based indices of the vector must corresponds to the 0-based indices of the particle list TPS()

Definition at line 69 of file ThermalModelEVDiagonal.cpp.

double thermalfist::ThermalModelEVDiagonal::MuShift ( int  i) const
protectedvirtual

The shift in the chemical potential of particle species i due to the excluded volume interactions.

Equal to \(-v_i p\)

Parameters
i0-based particle specie index
Returns
The shift in the chemical potential

Reimplemented from thermalfist::ThermalModelBase.

Definition at line 529 of file ThermalModelEVDiagonal.cpp.

virtual double thermalfist::ThermalModelEVDiagonal::ParticleKurtosis ( int  )
inlinevirtual

Kurtosis of primordial particle number fluctuations for species i.

Reimplemented from thermalfist::ThermalModelBase.

Definition at line 187 of file ThermalModelEVDiagonal.h.

double thermalfist::ThermalModelEVDiagonal::ParticleScalarDensity ( int  i)
virtual

The scalar density of the particle species i.

Implements thermalfist::ThermalModelBase.

Definition at line 514 of file ThermalModelEVDiagonal.cpp.

virtual double thermalfist::ThermalModelEVDiagonal::ParticleScaledVariance ( int  )
inlinevirtual

Scaled variance of primordial particle number fluctuations for species i.

Reimplemented from thermalfist::ThermalModelBase.

Definition at line 181 of file ThermalModelEVDiagonal.h.

virtual double thermalfist::ThermalModelEVDiagonal::ParticleSkewness ( int  )
inlinevirtual

Skewness of primordial particle number fluctuations for species i.

Reimplemented from thermalfist::ThermalModelBase.

Definition at line 184 of file ThermalModelEVDiagonal.h.

double thermalfist::ThermalModelEVDiagonal::Pressure ( double  P)
protected

Computes the l.h.s. of the transcdental equation for the pressure.

Parameters
PInput pressure (GeV fm \(^{-3}\))
Returns
Computed pressure (GeV fm \(^{-3}\))

Definition at line 154 of file ThermalModelEVDiagonal.cpp.

double thermalfist::ThermalModelEVDiagonal::PressureId ( int  i,
double  Pressure 
)
protected

Calculate the ideal gas pressure of particle species i for the given pressure value.

Parameters
i0-based particle specie index
PressureInput pressure (GeV fm \(^{-3}\))
Returns
Ideal gas pressure (fm \(^{-3}\))

Definition at line 142 of file ThermalModelEVDiagonal.cpp.

void thermalfist::ThermalModelEVDiagonal::ReadInteractionParameters ( const std::string &  filename)
virtual

Read the set of eigenvolume parameters for all particles from an external file.

Parameters
filenameInput file name

Reimplemented from thermalfist::ThermalModelBase.

Definition at line 78 of file ThermalModelEVDiagonal.cpp.

double thermalfist::ThermalModelEVDiagonal::ScaledVarianceId ( int  i,
double  Pressure 
)
protected

Calculate the ideal gas scaled variance of particle species i number fluctuations for the given pressure value.

Parameters
i0-based particle specie index
PressureInput pressure (GeV fm \(^{-3}\))
Returns
Ideal gas scaled variance

Definition at line 148 of file ThermalModelEVDiagonal.cpp.

void thermalfist::ThermalModelEVDiagonal::SetRadius ( double  )
virtual

Set the same excluded volume radius parameter for all species.

Parameters
radRadius parameter (fm)

Reimplemented from thermalfist::ThermalModelBase.

Examples:
BagModelFit.cpp.

Definition at line 50 of file ThermalModelEVDiagonal.cpp.

void thermalfist::ThermalModelEVDiagonal::SetRadius ( int  ,
double   
)
virtual

Set the radius parameter for particle species i.

Parameters
i0-based index of particle species
radRadius parameter (fm)

Reimplemented from thermalfist::ThermalModelBase.

Definition at line 545 of file ThermalModelEVDiagonal.cpp.

void thermalfist::ThermalModelEVDiagonal::SetVirial ( int  i,
int  j,
double  b 
)
virtual

Sets the eigenvolume parameter of particle species i.

Parameters
i0-based particle specie index
jHere must be equal to i
bThe eigenvolume parameter (fm \(^{-3}\))

Reimplemented from thermalfist::ThermalModelBase.

Definition at line 556 of file ThermalModelEVDiagonal.cpp.

void thermalfist::ThermalModelEVDiagonal::SolvePressure ( )
protectedvirtual

Solves the transcdental equation for the pressure.

Solves \( p(T,\mu) = \sum_i p_i^{\rm id} (T, \mu_i - v_i p). \)

Definition at line 165 of file ThermalModelEVDiagonal.cpp.

double thermalfist::ThermalModelEVDiagonal::VirialCoefficient ( int  i,
int  j 
) const
virtual

Return the eigenvolume parameter of particle species i.

Parameters
i0-based particle specie index
jIrrelevant
Returns
The eigenvolume parameter (fm \(^{-3}\))

Reimplemented from thermalfist::ThermalModelBase.

Definition at line 551 of file ThermalModelEVDiagonal.cpp.

void thermalfist::ThermalModelEVDiagonal::WriteInteractionParameters ( const std::string &  filename)
virtual

Write the set of eigenvolume parameters for all particles to an external file.

One particle specie per line.

Parameters
filenameOutput file name

Reimplemented from thermalfist::ThermalModelBase.

Definition at line 102 of file ThermalModelEVDiagonal.cpp.

Member Data Documentation

std::vector<double> thermalfist::ThermalModelEVDiagonal::m_densitiesid
protected

Vector of ideal gas densities with shifted chemical potentials

Definition at line 253 of file ThermalModelEVDiagonal.h.

std::vector<double> thermalfist::ThermalModelEVDiagonal::m_densitiesidnoshift
protected

Vector of ideal gas densities without shifted chemical potentials

Definition at line 254 of file ThermalModelEVDiagonal.h.

double thermalfist::ThermalModelEVDiagonal::m_Densityid
protected

Definition at line 258 of file ThermalModelEVDiagonal.h.

double thermalfist::ThermalModelEVDiagonal::m_Pressure
protected

The solved pressure

Definition at line 257 of file ThermalModelEVDiagonal.h.

EVSolution thermalfist::ThermalModelEVDiagonal::m_sol
protected

Axuiliary

Definition at line 260 of file ThermalModelEVDiagonal.h.

double thermalfist::ThermalModelEVDiagonal::m_Suppression
protected

The common density suppression factor

Definition at line 256 of file ThermalModelEVDiagonal.h.

double thermalfist::ThermalModelEVDiagonal::m_TotalDensity
protected

Definition at line 259 of file ThermalModelEVDiagonal.h.

std::vector<double> thermalfist::ThermalModelEVDiagonal::m_v
protected

Vector of eigenvolumes of all hadrons

Definition at line 255 of file ThermalModelEVDiagonal.h.


The documentation for this class was generated from the following files: