Merge pull request #16 from vlvovch/devel

Version 1.1

CMakeLists.txt

@@ -7,7 +7,7 @@ project (ThermalFIST)
 
 # The version number.
 set (ThermalFIST_VERSION_MAJOR 1)
-set (ThermalFIST_VERSION_MINOR 0)
+set (ThermalFIST_VERSION_MINOR 1)
 set (ThermalFIST_VERSION_DEVEL 0)
 
 # configure a header file to pass some of the CMake settings
@@ -42,18 +42,19 @@ endif(CMAKE_BUILD_TYPE STREQUAL "" AND NOT CMAKE_CXX_FLAGS MATCHES "-O[123]")
 include(CheckCXXCompilerFlag)
 CHECK_CXX_COMPILER_FLAG("-std=c++11" COMPILER_SUPPORTS_CXX11)
 CHECK_CXX_COMPILER_FLAG("-std=c++0x" COMPILER_SUPPORTS_CXX0X)
+#set(CMAKE_CXX_FLAGS "${CMAKE_CXX_FLAGS} -std=c99 -pedantic")
 if(COMPILER_SUPPORTS_CXX11)
     set(CMAKE_CXX_FLAGS "${CMAKE_CXX_FLAGS} -std=c++11")
 elseif(COMPILER_SUPPORTS_CXX0X)
     set(CMAKE_CXX_FLAGS "${CMAKE_CXX_FLAGS} -std=c++0x")
-#else()
+else()
 #        message(STATUS "The compiler ${CMAKE_CXX_COMPILER} has no C++11 support. Please use a different C++ compiler.")
 endif()
 # Here it is setting the Visual Studio warning level to 4
 if(MSVC)
 set (CMAKE_CXX_FLAGS "${CMAKE_CXX_FLAGS} /W4")
 else(MSVC)
-set (CMAKE_CXX_FLAGS "${CMAKE_CXX_FLAGS}")
+set (CMAKE_CXX_FLAGS "${CMAKE_CXX_FLAGS} -pedantic")
 endif(MSVC)
 
 
@@ -72,6 +73,15 @@ else(result)
   endif()
 endif(result)
 
+# Check whether the size of int is at least 32-bit
+INCLUDE (CheckTypeSize)
+check_type_size("int" SIZEOFINT BUILTIN_TYPES_ONLY)
+if (${SIZEOFINT} LESS "4")
+message(FATAL_ERROR "Size of int of at least 4 bytes is required. Your compiler has " ${SIZEOFINT} " bytes. Please get \
+a different compiler.")
+endif()
+message(STATUS "Size of int is " ${SIZEOFINT} " bytes")
+
 add_definitions(-DUSE_MINUIT)
 
 if(NOT STANDALONE_MINUIT)

README.md

@@ -28,7 +28,8 @@ For more details please see the [**Thermal-FIST documentation**](https://fias.un
 
 ### Notes
 Calculations may break down for excessive/overstressed parameters.
-This should be carefully monitored.
+This should be carefully monitored, and the accuracy of the output assessed.
+If in doubt, contact the program author(s).
 
 ### Third-party code used
 

docs/doxygen/Doxyfile

@@ -38,7 +38,7 @@ PROJECT_NAME           = Thermal-FIST
 # could be handy for archiving the generated documentation or if some version
 # control system is used.
 
-PROJECT_NUMBER         = 1.0
+PROJECT_NUMBER         = 1.1
 
 # Using the PROJECT_BRIEF tag one can provide an optional one line description
 # for a project that appears at the top of each page and should give viewer a

docs/quickstart.md

@@ -40,12 +40,14 @@ From the build folder, run with
 
 Using the GUI should be self-explanatory for most part. 
 
-The particle list is read from file.
+The particle list is read from a file.
 The default list is located in the `$(ThermalFIST)/input/list/PDG2014` folder, containing particles consisting of light and strange quarks in accordance with the PDG2014 compilation.
-There the `list-withnuclei.dat` file contains hadrons plus light nuclei, while the `list.dat` file contains hadrons only. In the same folder there is also the `decays.dat` where all the decay channels are listed.
-Additionally, the `$(ThermalFIST)/input/list/` folder also contains particle lists extracted from the open source THERMUS-2.3/3.0 ([link](http://www.phy.uct.ac.za/phy/people/academic/wheaton/research)) package, which can be used to compare the results between different lists and also compare Thermal-FIST to THERMUS.
+There the `list-withnuclei.dat` file contains hadrons plus light nuclei, while the `list.dat` file contains hadrons only. There are also two files with a `-withcharm` prefix which contain charm hadrons
+in addition to light and strange hadrons.
+In the same folder there is also the `decays.dat` file where all the decay channels are listed.
+Additionally, the `$(ThermalFIST)/input/list/` folder also contains particle lists extracted from the open source THERMUS-2.3/3.0 ([link](http://www.phy.uct.ac.za/phy/people/academic/wheaton/research)) package, which can be used to compare the results between different lists and also to compare Thermal-FIST and THERMUS.
 
-The **QtThermalFIST** GUI program consists of four tabs:
+The **QtThermalFIST** GUI program consists of five tabs:
 
 #### *Thermal model* tab
 
@@ -71,7 +73,7 @@ The model configuration is specified in the right-hand side panel.
 **NOTE**
 
 Analytic calculations within the strangeness-canonical ensemble in the presence of the EV/vdW interactions are implemented, but use approximatations. 
-More specifically, it is assumied that strange particles are a very small part of the whole system. This is appropriate for low energies (e.g. HADES), but not for high energies (like LHC). Use with caution!
+More specifically, it is assumed that strange particles form a very small part of the whole system. This is appropriate for low energies (e.g. HADES), but not for high energies (like LHC). Use with caution!
 
 ---
 
@@ -80,13 +82,13 @@ More specifically, it is assumied that strange particles are a very small part o
 The possibilities are:
 
 1. Zero-width approximation
-2. Energy-independent Breit-Wigner
+2. Energy-independent Breit-Wigner, restricted to a ±2Γ interval around the pole mass
 3. Energy-dependent Breit-Wigner (eBW)
-4. Energy-dependent Breit-Wigner but with constant branching ratios when evaluating feeddown.
+4. Energy-dependent Breit-Wigner but with constant branching ratios when evaluating feeddown
 
 **Statistics.** Maxwell-Boltzmann or quantum statistics. 
 
-For quantum statistics the default method (*Use quadratures* checkbox checked) is to use Gauss-Legendre quadratures to compute the relevant integrals numerically. A faster method is the series expansion in Bessel functions. This method however does not work well if chemical potential is close to or larger than particle's mass. 
+For quantum statistics the default method (*Use quadratures* checkbox checked) is to use Gauss-Legendre quadratures to compute the relevant integrals numerically. A faster method is the series (cluster) expansion in Bessel functions. This method however does not work well if chemical potential is close to or larger than particle's mass. 
 
 It is possible to include quantum statistics for pions or for mesons only. This can speed up significantly calculations in the canonical ensemble.
 
@@ -99,58 +101,60 @@ Chemical potentials in the GCE can be constrained by various conservation laws.
 4. Charm chemical potential can be constrained from the charm neutrality condition.
 
 **Excluded volume/van der Waals interactions.** 
-Set the interaction parameters for excluded volume or van der Waals HRG model. The are two possibilities.
+Set the interaction parameters for an excluded volume or a van der Waals HRG model. The are two possibilities.
 
-1. Set parameters manually, by inputting the values of a and b. Interpreting the input parameters as those of protons, it is possible to assign he same parameters to other particles, impose a mass-proportional (bag model), or baryon conent proportional scaling of parameters. For Crossterms EV or QvdW models it is possible to switch off interactions for various hadron pairs such as meson-meson, meson-baryon, baryon-baryon, and baryon-antibaryon.
-2. The interaction parameters can be read from external file. 
+1. Set parameters manually, by providing the values of the a and b parameters. Interpreting these input parameters as those of protons, it is possible to assign the same parameters to all other particles, or impose a mass-proportional (bag model) or baryon conent proportional scaling of parameters. For Crossterms EV or QvdW models it is also possible to switch off interactions for various hadron pairs such as meson-meson, meson-baryon baryon-baryon, and baryon-antibaryon.
+2. The interaction parameters can be read from an external file. 
    * For *Diagonal EV* model the external file should contain a list of rows (one for each species) with two columns: 1 - PDGID, 2 - v_i. 
    * For *Crossterms EV*, one row for each pair of species, three columns: 1 - PDGID1, 2 - PDGID2, 3 - b_ij. 
    * For *QvdW-HRG*, one row for each pair of species, four columns: 1 - PDGID1, 2 - PDGID2, 3 - b_ij, 4 - a_ij. 
-   * If some species/pair of species is not found in the input file, then the parameters are zero for that. Some sample input files are provided in the `$(ThermalFIST)/input/list/interaction` folder.
+   * If some species/pair of species is not found in the input file, then the parameters are assumed to be zero for those. A number of sample input files is provided in the `$(ThermalFIST)/input/list/interaction` folder.
 
-The button *Calculate* invokes calculation of all the hadron yields, equation of state properties, and, optionally, fluctuations, for the specified values of thermal parameters.
+The button *Calculate* invokes a calculation of all the hadron yields, equation of state properties, and, optionally, fluctuations, for the specified values of thermal parameters.
 
-The button *Equation of state...* shows the equation of state properties corresponding to a given calculation, as well as some extra quantities, such as the second-order fluctuations incuding the contributions from probabilistic decays.
+The button *Equation of state...* shows the equation of state properties corresponding to a given calculation, as well as some extra quantities, such as the second-order fluctuations incuding contributions from probabilistic decays.
 
-Double-clicking on a particle in the list will open a window with additional useful information about the particle, including e.g. all the feeddown contributions from different resonance decays.
+Double-clicking on a particle in the list will open a window with additional useful information about the particle, including, e.g., all the feeddown contributions from different resonance decays.
 
 
 #### *Thermal fits* tab
 
 Thermal fitting. The HRG model specification is done in the same manner as in the previous tab.
 
 It is possible to specify which parameters should be fitted, their initial/fixed values, and the value range where the fit is performed.
-For fits in the (partially) canonical picture it also possible to fix the correlation volume $V_c$ relative to the total volume $V$.
+For fits in the (partially) canonical picture it also possible to fix the value of the correlation volume $V_c$ relative to the total volume $V$.
 
 The experimental data can be loaded from an external file. Some samples are provided in the [$(ThermalFIST)/input/data](../input/data) folder. Possibility to input the data directly in the GUI is provided as well (double-click on the yield in the table to edit it).
-It is possible to view the fit results directly within the GUI in a form of the common thermal fit plots.
+It is possible to view the fit results directly within the GUI in a form of some of the common thermal fit plots.
 Another useful feature is the analysis of the $χ^2$ profiles (the *Chi2 profile...* button), which can be done after a global fit was performed.
 
 #### *Equation of state* tab
 
-Offers a possibility to study the temperature dependence at fixed $μ_B$ (or at fixed S/B) 
+Offers a possibility to study the temperature dependence at fixed μ_B (or at fixed S/B) 
 of some common equation of state observables, conserved charges susceptibilities, and number densities (primordial or with feeddown). It is also possible to consider ratios of any pair of these observables.
-At finite $μ_B$ or S/B. 
-Chemical potentials $μ_Q$, $μ_S$, and/or $μ_C$ can be fixed by conservation laws, or set to zero otherwise.
+For a finite μ_B or S/B value, 
+the chemical potentials μ_Q, μ_S, and/or μ_C can be fixed by conservation laws, or set to zero otherwise.
 The thermal model specification is done similarly to the previous two tabs, except that here only the grand-canonical ensemble is considered.
-For a number of observables at $μu_B = 0$ the published lattice QCD data of the Wuppertal-Budapest and/or HotQCD collaborations is plotted along with the calculation results for convenience.
+For a number of observables at $μ_B = 0$ the published lattice QCD data of the Wuppertal-Budapest and/or HotQCD collaborations is plotted along with the calculation results for convenience.
 
 #### *Event generator* tab
 
 A Monte Carlo generator which generates events with particle numbers distributed according to the corresponding multiplicity distribution in a HRG model, while momenta are distributed in accordance with the Blast-Wave model (spherically or cylindrically symmetric). Currently restricted to the case of the Maxwell-Boltzmann statistics.
-In theory provides the same yield, fluctuations etc. as analytic calculations (in case of EV/vdW interactions a sufficiently large volume is necessary), also applicable to study combined effect of EV/vdW interactions and exact charge conservation.
-Can be useful in detector response simulations.
+In theory provides the same yield, fluctuations etc. as analytic calculations (in case of EV/vdW interactions a sufficiently large system volume is necessary), also applicable to study combined effect of EV/vdW interactions and exact charge conservation.
+Can be useful, e.g., in detector response simulations.
 
 #### *Particle list editor* tab
 
 Allows to edit the particle list and decay channels on the fly.
-Note that antiparticles are not listed, but are generated autmatically from particles if at least one of the quantum numbers (baryon number, electric charge, strangeness, charm) is non-zero.
-All changes can be saved to file.
+Note that antiparticles are not listed, but are generated automatically from particles if at least one of the quantum numbers (baryon number, electric charge, strangeness, charm) is non-zero.
+All changes can be saved to a file.
 
 ### Other notes
 
 - Sample C++ macros which use Thermal-FIST as a library can be found in the [$(ThermalFIST)/src/examples](../src/examples) folder
+
 - An example of using the Thermal-FIST library as a submodule is provided by [this repository](https://github.com/vlvovch/1807.02079)
+
 - The canonical ensemble calculations use the method from [this paper](https://arxiv.org/abs/nucl-th/0112021) to perform analytically the integration over the baryon number fugacity whenever this is possible. This allows to significantly speed up the calculations. This method cannot be used if quantum statisics is considered for baryons, or if there are multi-baryons (light nuclei) in the particle list. In this case direct numerical integration is performed, which can make calculations very slow.
 
 ### Common issues
@@ -167,4 +171,6 @@ All changes can be saved to file.
 
 - Program outputs warnings about the Bose-Einstein condensation. This occurs if chemical potentials of bosons at some point exceed their mass. This sometimes happens in the process of constaining electric/strange chemical potentials through conservation laws. If these warnings eventually disappear and final values of chemical potentials are reasonable, things should be fine.
 
+- Accuracy of the calculations in the canonical ensemble will eventually break down as the system volume is increased to a large value, such as that corresponding to central heavy-ion collisions. Accuracy of a canonical ensemble calculation can be verified by comparing the values of the exactly conserved charges computed in the model with the given ones.
+
 - The GUI program appears to hang when the cylindrically-symemetric event generator is invoked. This event generator has a fairly long initialization time (approx. 20-30 sec) after which the program proceeds smoothly.

include/HRGBase/IdealGasFunctions.h

@@ -37,6 +37,9 @@ namespace thermalfist {
      */
     enum QStatsCalculationType { ClusterExpansion, Quadratures };
 
+    /// \brief Whether \mu > m Bose-Einstein condensation issue was encountered for a Bose gas
+    extern bool calculationHadBECIssue;
+
     /**
      * \brief Computes the particle number density of a Maxwell-Boltzmann gas.
      * 

include/HRGBase/NumericalIntegration.h

@@ -207,7 +207,7 @@ namespace thermalfist {
      * \param [out] xleg Gauss-Legendre nodes for the variable y. 
      * \param [out] wleg Gauss-Legendre weights for the variable y. 
      */
-    void GetCoefs2DLegendre32Legendre32(double ax, double bx, double ay, double by, std::vector<double> *xlag, std::vector<double> *wlag, std::vector<double> *xleg, std::vector<double> *wleg);
+    void GetCoefs2DLegendre32Legendre32(double ax, double bx, double ay, double by, std::vector<double> *xleg1, std::vector<double> *wleg1, std::vector<double> *xleg2, std::vector<double> *wleg2);
 
     /**
      * Populates the nodes and weights for integrating

include/HRGBase/ParticleDecay.h

@@ -71,9 +71,9 @@ namespace thermalfist {
    * \brief Structure containing information about a single decay channel of a particle.
    */
   struct ParticleDecayChannel {
-    double mBratio;                /**< Branching ratio */
-    std::vector<int> mDaughters;   /**< PDGID numbers of daughter particles */
-    double mM0;                    /**< Sum of masses of decay products */
+    double mBratio;                      /**< Branching ratio */
+    std::vector<long long> mDaughters;   /**< PDGID numbers of daughter particles */
+    double mM0;                          /**< Sum of masses of decay products */
 
     double mPole;
     double mL;                      /**< Orbital angular momentum for decay, used in the eBW scheme */
@@ -89,7 +89,7 @@ namespace thermalfist {
      * \param bratio Branching ratio of the decay. Between 0 and 1.
      * \param daughters A vector of PDG ID numbers of all daughter products.
      */
-    ParticleDecayChannel(double bratio = 0., const std::vector<int> &daughters = std::vector<int>(0)) :
+    ParticleDecayChannel(double bratio = 0., const std::vector<long long> &daughters = std::vector<long long>(0)) :
       mBratio(bratio), mDaughters(daughters), mM0(0.), mPole(0.), mL(0.),
       mBratioVsM(std::vector<double>(0)), mBratioAverage(bratio), mChannelName("decay1") {
     }