TUHKMgen/UHKM/GrandCanonical.cxx

   1 /*
   2
   3         Nikolai Amelin, Ludmila Malinina, Timur Pocheptsov (C) JINR/Dubna
   4       amelin@sunhe.jinr.ru, malinina@sunhe.jinr.ru, pocheptsov@sunhe.jinr.ru
   5                            November. 2, 2005
   6
   7 */
   8 #include <iostream>
   9
  10 #include <TError.h>
  11 #include <TMath.h>
  12 #ifndef GRANDCANONICAL_INCLUDED
  13 #include "GrandCanonical.h"
  14 #endif
  15 #ifndef HANKELFUNCTION_INCLUDED
  16 #include "HankelFunction.h"
  17 #endif
  18 #ifndef UKUTILITY_INCLUDED
  19 #include "UKUtility.h"
  20 #endif
  21
  22 GrandCanonical::GrandCanonical():
  23   fTemperature(-1111),
  24   fBaryonPotential(-1111),
  25   fStrangePotential(-1111),
  26   fElectroPotential(-1111),
  27   fNMax(-1111),
  28   fInitialized(kFALSE)
  29 {
  30 }
  31
  32 GrandCanonical::GrandCanonical(Int_t nmax, Double_t temperature, Double_t baryonPotential, Double_t strangePotential, Double_t electroPotential):
  33   fTemperature(temperature),
  34   fBaryonPotential(baryonPotential),
  35   fStrangePotential(strangePotential),
  36   fElectroPotential(electroPotential),
  37   fNMax(nmax),
  38   fInitialized(kTRUE)
  39 {
  40 }
  41
  42 GrandCanonical::~GrandCanonical() {}
  43
  44
  45 void GrandCanonical::Temperature(Double_t value) {
  46   fTemperature = value;
  47   if(fNMax!=-1111 && fBaryonPotential!=-1111 && fStrangePotential!=-1111 && fElectroPotential!=-1111)
  48     fInitialized = kTRUE;
  49 }
  50
  51 void GrandCanonical::BaryonPotential(Double_t value) {
  52   fBaryonPotential = value;
  53   if(fNMax!=-1111 && fTemperature!=-1111 && fStrangePotential!=-1111 && fElectroPotential!=-1111)
  54     fInitialized = kTRUE;
  55 }
  56
  57 void GrandCanonical::StrangePotential(Double_t value) {
  58   fStrangePotential = value;
  59   if(fNMax!=-1111 && fTemperature!=-1111 && fBaryonPotential!=-1111 && fElectroPotential!=-1111)
  60     fInitialized = kTRUE;
  61 }
  62
  63 void GrandCanonical::ElectroPotential(Double_t value) {
  64   fElectroPotential = value;
  65   if(fNMax!=-1111 && fTemperature!=-1111 && fBaryonPotential!=-1111 && fStrangePotential!=-1111)
  66     fInitialized = kTRUE;
  67 }
  68
  69 void GrandCanonical::NMax(Int_t value) {
  70   fNMax = value;
  71   if(fTemperature!=-1111 && fBaryonPotential!=-1111 && fStrangePotential!=-1111 && fElectroPotential!=-1111)
  72     fInitialized = kTRUE;
  73 }
  74
  75 Double_t GrandCanonical::ParticleEnergyDensity(ParticlePDG* particle) {
  76   // Check if all the thermodinamic parameters are set
  77   if(!fInitialized)
  78     Fatal("GrandCanonical::ParticleEnergyDensity", "GrandCanonical object not fully initialized!!");
  79
  80   // Compute the particle energy density
  81   Double_t degFactor = 2.*particle->GetSpin() + 1.;                                    // degeneracy factor
  82   Double_t mass = particle->GetMass();                                                // PDG table mass
  83   Double_t d = Int_t(2.*particle->GetSpin()) & 1 ? 1. : -1;                                   // Bose-Einstein/Fermi-Dirac factor
  84   Double_t preFactor = (degFactor*mass*mass*fTemperature*fTemperature/hbarc/hbarc/hbarc)/(2.*TMath::Pi()*TMath::Pi());
  85
  86   Double_t postFactor = 0.;
  87   //compute chemical potential
  88   Double_t potential = fBaryonPotential * particle->GetBaryonNumber() +
  89                        fStrangePotential * particle->GetStrangeness() +
  90                        fElectroPotential * particle->GetElectricCharge();
  91
  92   for(Int_t n = 1; n <= fNMax; ++n) {
  93     postFactor += TMath::Power(-d, n+1)/(n*n) *
  94                   TMath::Exp(n*potential/fTemperature) *
  95                   (3.*HankelKn(2, n*mass/fTemperature) + n*mass/fTemperature*HankelK1(n*mass/fTemperature));
  96   }
  97   return preFactor * postFactor;
  98 }
  99
 100 Double_t GrandCanonical::ParticleNumberDensity(ParticlePDG* particle) {
 101   // Check if all the thermodinamic parameters are set
 102   if(!fInitialized)
 103     Fatal("GrandCanonical::ParticleNumberDensity", "GrandCanonical object not fully initialized!!");
 104
 105   Double_t degFactor = 2.*particle->GetSpin() + 1.;
 106   Double_t mass = particle->GetMass();
 107   Double_t d = Int_t(2*particle->GetSpin()) & 1 ? 1. : -1.;
 108   Double_t preFactor = (degFactor*mass*mass*fTemperature/hbarc/hbarc/hbarc)/(2.*TMath::Pi()*TMath::Pi());
 109
 110   Double_t postFactor = 0.;
 111   Double_t potential = fBaryonPotential * particle->GetBaryonNumber() +
 112                        fStrangePotential * particle->GetStrangeness() +
 113                        fElectroPotential * particle->GetElectricCharge();
 114
 115   for(Int_t n = 1; n <= fNMax; ++n) {
 116     postFactor += TMath::Power(-d, n+1)/n *
 117                   TMath::Exp(n*potential/fTemperature) *
 118                   HankelKn(2, n*mass/fTemperature);
 119   }
 120   return preFactor * postFactor;
 121 }
 122
 123
 124 Double_t GrandCanonical::EnergyDensity(DatabasePDG* database) {
 125   // Check if all the thermodinamic parameters are set
 126   if(!fInitialized)
 127     Fatal("GrandCanonical::EnergyDensity", "GrandCanonical object not fully initialized!!");
 128
 129   Double_t meanEnergyDensity = 0.;
 130
 131   for(Int_t currParticle = 0; currParticle<database->GetNParticles(); currParticle++) {
 132     ParticlePDG *particle = database->GetPDGParticleByIndex(currParticle);
 133     meanEnergyDensity += ParticleEnergyDensity(particle);
 134   }
 135
 136   return meanEnergyDensity;
 137 }
 138
 139 Double_t GrandCanonical::BaryonDensity(DatabasePDG* database) {
 140   // Check if all the thermodinamic parameters are set
 141   if(!fInitialized)
 142     Fatal("GrandCanonical::BaryonDensity", "GrandCanonical object not fully initialized!!");
 143
 144   Double_t meanBaryonDensity = 0.;
 145
 146   for(Int_t currParticle = 0; currParticle<database->GetNParticles(); currParticle++) {
 147     ParticlePDG *particle = database->GetPDGParticleByIndex(currParticle);
 148     meanBaryonDensity += ParticleNumberDensity(particle)*particle->GetBaryonNumber();
 149   }
 150   return meanBaryonDensity;
 151 }
 152
 153 Double_t GrandCanonical::StrangeDensity(DatabasePDG* database) {
 154   // Check if all the thermodinamic parameters are set
 155   if(!fInitialized)
 156     Fatal("GrandCanonical::StrangeDensity", "GrandCanonical object not fully initialized!!");
 157
 158   Double_t meanStrangeDensity = 0.;
 159
 160   for(Int_t currParticle = 0; currParticle<database->GetNParticles(); currParticle++) {
 161     ParticlePDG *particle = database->GetPDGParticleByIndex(currParticle);
 162     meanStrangeDensity += ParticleNumberDensity(particle)*particle->GetStrangeness();
 163   }
 164
 165   return meanStrangeDensity;
 166 }
 167
 168 Double_t GrandCanonical::ElectroDensity(DatabasePDG* database) {
 169   // Check if all the thermodinamic parameters are set
 170   if(!fInitialized)
 171     Fatal("GrandCanonical::ElectroDensity", "GrandCanonical object not fully initialized!!");
 172
 173   Double_t meanElectroDensity = 0.;
 174
 175   //hadrons
 176   for(Int_t currParticle = 0; currParticle<database->GetNParticles(); currParticle++) {
 177     ParticlePDG *particle = database->GetPDGParticleByIndex(currParticle);
 178     meanElectroDensity += ParticleNumberDensity(particle)*particle->GetElectricCharge();
 179   }
 180
 181   return meanElectroDensity;
 182 }