Susceptibilities with multiquark interactions in PNJL model
Abstract
We have investigated the fluctuations and the higher order susceptibilities of quark number, isospin number, electric charge and strangeness at vanishing chemical potential for 2+1 flavor Polyakov loop extended Nambu–JonaLasinio model. The calculations are performed for the bound effective potential in the quark sector requiring up to eight quark interaction terms. These have been contrasted to the lattice results which currently have somewhat heavier quarks in the light flavor sector. The results show sufficient qualitative agreement. For comparison we also present the results obtained with the conventional effective potential containing upto six quark interaction terms.
pacs:
12.38.Aw, 12.38.Mh, 12.39.xI Introduction
Confinement and chiral symmetry breaking are the most fundamental properties of strong interaction physics at low temperature and density, where the physics is mainly governed by the nonperturbative QCD. In principle, the deconfinement phase transition and chiral phase transition are defined in two extreme limits of current quark mass. Deconfinement phase transition and its order parameter is well defined for infinite current quark mass and chiral phase transition is exact for zero quark mass. But in real world with a finite value of quark masses, the nature of these two phase transitions is an open question. Another important feature of the QCD phase diagram is the existence of the critical end point (CEP), where first order phase transition, from hadronic phase to quarkgluonplasma (QGP) phase, ends. However, the exact location of CEP is still unknown. Investigation of these properties of strongly interacting matter are necessary to understand the various astrophysical and cosmological scenario.
The experimental explorations to understand the properties of strongly interacting matter has been studied at Relativistic Heavyioncollider (RHIC), Brookhaven adcox where the heavy nuclei collide with each other at relativistic energies to form hot and dense strongly interacting matter. More data are expected from LHC and FAIR in future. But to analyze the data from the experiment, we need a thorough understanding of the theory of strong interaction physics.
Due to our limited knowledge of the nonperturbative physics, QCD, which is the theory of strong interaction, can not be used to study the phase transition picture. In this regard Lattice QCD (LQCD) provides the most direct approach to study QCD at high temperature Boyd ; engels ; fodor1 ; fodor2 ; allton1 ; allton2 ; allton3 ; forcrand ; aoki1 ; aoki2 ; arriola1 ; arriola2 ; arriola3 ; arriola4 . However LQCD has its own restrictions due to the discretization of spacetime. Furthermore, at finite chemical potential, LQCD faces the well known sign problem.
Another approach to study low energy limit of QCD and the QCD phase transition is to use effective theories of QCD. Polyakov loop extended Nambu–JonaLasinio model (PNJL) is one of the successful approach which combines the confinement and chiral symmetry breaking properties in a simple formalism. There have been series of work to study the thermodynamic properties of strongly interacting matter using PNJL model for both 2 flavor and 2+1 flavor fuku1 ; ratti1 ; pisarski ; fuku2 ; ratti2 ; gatto ; ghosh ; rayvdm ; deb1 . These studies suggest that this model reproduces the zero density lattice data quite successfully. However, the vacuum of the NJL part of the PNJL model seems to to be unbound in a 2+1 flavor scenario. A plausible solution of this problem has been proposed by Osipov et.al using eightquark interaction term osipov1 ; osipov2 ; osipov3 ; osipov4 . Also the 2 flavor PNJL model have been studied in Ref. kashiwa1 ; kashiwa2 with eightquark interaction. In our previous work, we developed 2+1 flavor PNJL model with eightquark interaction terms in the Lagrangian with threemomentum cutoff regularisation deb2 .
The thermodynamic aspect of the phase transition from hadronic phase to QGP phase can be understood properly if we study the thermodynamic variables like quark number susceptibility (QNS), isospin number susceptibility (INS), specific heat and speed of sound () etc. Susceptibilities are related to fluctuations via the fluctuationdissipation theorem. A measure of the intrinsic statistical fluctuations in a system close to thermal equilibrium is provided by the corresponding susceptibilities. At zero chemical potential, charge fluctuations are sensitive indicators of the transition from hadronic matter to QGP. Also the existence of the CEP can be signalled by the divergent fluctuations. For the small net baryon number, which can be met at different experiments, the transition from hadronic to QGP phase is continuous and the fluctuations are not expected to lead any singular behavior. Recently, the computations on the lattice have been performed for many of these susceptibilities at zero chemical potentials gottlieb ; gavai ; bernard1 ; bernard2 . It was shown that at vanishing chemical potential the susceptibilities rise rapidly around the continuous crossover transition region.
The study of higher order moments of fluctuations are also necessary to locate the transition point more accurately. In 2 flavor QCD, it has been shown that the quark number and isospin number fluctuations increase with temperature and their fourth moments start to show pronounced peaks in the transition region from low to high temperature allton3 ; ejiri . In fact the higher order coefficients become increasingly sensitive in the vicinity of phase transition. Fluctuations are computed with respect to the quark chemical potential in 2 flavor PNJL model with threemomentum cutoff regularisation ghosh ; roessner1 ; friman ; ray2 . Also QNS at finite density has been estimated in some works within 2 flavor PNJL model roessner2 . Recently the idea of the Taylor expansion in terms of chemical potential for PNJL model has been computed within the constraint that the net strange quark density is zero, which is the case in ultra relativistic heavy ion collision fuku3 . There have been recent calculations towards the fluctuations in 2+1 flavor Lattice QCD cheng1 ; cheng ; fodor10 and also in the 2+1 flavor PNJL model taking upto sixquark determinant interaction terms fuku2 ; wu1 . Similar calculations have been carried out in Polyakov loop coupled quarkmeson (PQM) model schaefer1 ; schaefer2 ; schaefer3 ; schaefer4 and its renormalization group improved version skokov .
In this paper we have investigated the susceptibilities within 2+1 flavor PNJL model framework with two different kinds of NJL interaction. In one case we take the conventional form which takes into account upto sixquark determinant interaction, which we will call as model A. In the other case an extra eightquark interaction term is added to make the effective potential bound, which will be hereafter denoted by model B. We have also studied the specific heat and the speed of sound. Specific heat is related to the eventbyevent temperature fluctuations stod and mean transverse momentum fluctuations korus in heavyion reactions. These fluctuations show diverging behavior near the critical end point (CEP). The speed of sound determines the flow properties in heavyion reactions sorge ; kolb .
Our paper is organized as follows: In Sec. II, the basic formalism of the PNJL model and also the calculation of fluctuations have been discussed. Various thermodynamic quantities such as specific heat, speed of sound are also defined in this section. In the next section we describe our results and compared with the recent lattice data and also with other models like Polyakov extended quarkmeson model (PQM). In the last section we conclude.
Ii Formalism
ii.1 Thermodynamic Potential
The PNJL model was formulated to study the chiral properties and the confinement physics of the QCD phase transition at finite temperature and density. In this model quark dynamics is studied with a background gauge field having only the temporal component. For a detailed review of the PNJL model with 2 flavor and 2+1 flavor see Ref. fuku1 ; ghosh ; gatto ; roessner1 ; deb1 ; friman ; roessner2 ; skokov ; ray2 ; fuku2 ; fuku3 ; wu1 . Osipov et.al pointed out that the effective potential of the NJL part is seems to be unbound in the conventional form and they have introduced the eightquark interaction terms in NJL model to stabilize the vacuum osipov1 ; osipov2 ; osipov3 ; osipov4 . Also the 2 flavor PNJL model have been studied in Ref. kashiwa1 ; kashiwa2 with eightquark interaction. In our previous paper we developed the 2+1 PNJL model with eightquark interactions within threemomentum cutoff scheme deb2 for the first time. We have reproduced the thermodynamic properties of QCD, calculated on Lattice, at zero baryon density quite satisfactorily. The thermodynamic potential for the multifermion interaction in the mean field approximation (MFA) of the PNJL model can be written as deb2 ,
(1) 
and are the fourquark and sixquark coupling constants respectively and and are the eightquark coupling constants. Here denotes chiral condensate of the quark with flavor and is the single quasiparticle energy. Here, constituent mass of flavor is given by the selfconsistent gap equation;
where , and take the labels of flavor , and in cyclic order. So, when then and and so on. In the above integrals, the vacuum integral has a cutoff whereas the medium dependent integrals have been extended to infinity. Polyakov loop being the normalized trace of the Wilson line L, which is an SU(3) matrix, should lie in the range 01. But earlier studies of PNJL model ghosh ; ray2 ; deb1 show that the becomes greater than 1 above 2. To solve this problem one has to take a proper Jacobian of transformation from the matrix valued field L to the complex valued field . This will constrain the value of within 1. Thus one has to modify the Polyakov loop potential by introducing Vandermonde (VdM) term. The necessity of the modification was reflected through an excellent agreement of flavour mixing effects in PNJL model with lattice data as shown in Ref. rayvdm . There are a few extra terms in in this prescription as compared to other recent works fuku2 ; roessner1 . These extra terms are put on the basis of the global Z(3) symmetry of the Polyakov loop potential and therefore their presence seem to be quite natural though may not be absolutely necessary as shown by fuku2 ; roessner1 . The modified potential can be expressed as ,
(2) 
where is the LandauGinzburg type potential given by ratti1 ,
(3) 
with,
(4) 
and , being constants. in eqn. (2) is known as VdM determinant rayvdm , is given by,
is a phenomenological constant. Polyakov loop and its charge conjugate is defined as,
The parameter is taken as , whereas the lattice determines its value to be for pure gauge theory. The reason to take a lower value of is to get the crossover temperature () consistent with the lattice data. In this work we have taken the parameter set obtained in our previous paper deb2 . For fixing the parameters we have used the following physical conditions;
and is kept fixed at 5.5 MeV. The parameters are given in table 1 for both PNJL model A and PNJL model B.
Interaction  






For the Polyakov loop potential we choose the parameters which reproduces the lattice data of pure gauge thermodynamics Boyd . It was shown in Ref. ratti1 that, pure gauge Lattice QCD data of scaled pressure, entropy and energy density are reproduced extremely well in Polyakov loop model using the ansatz (3) and (4) with parameters summarized below,
ii.2 Taylor expansion of pressure
The pressure of the strongly interacting matter can be written as,
(5) 
where is the temperature, is the quark chemical potential, is the charge chemical potential and is the strangeness chemical potential. From the usual thermodynamic relations we can show that the first derivative of pressure with respect to gives the quark number density and the second derivative is the quark number susceptibility (QNS).
Our first job is to minimize the thermodynamic potential numerically with respect to the fields , , , and . The values of the fields can then be used to evaluate the pressure using the equation (5). Then we can expand the scaled pressure at a given temperature in a Taylor series for the chemical potentials , , as,
(6) 
where,
(7) 
The flavor chemical potentials , , are related to , , by,
(8) 
It should be mentioned that one can also choose the independent chemical potentials as , , . Then eqn.(8) becomes,
(9) 
where, is the isospin chemical potential. Eq.(6) is a general expression of Taylor expansion of pressure for different chemical potentials. Since in this paper we are only concerned with the diagonal terms of the expansion, we can write the above equation in simpler way as,
(10) 
where,
(11) 
Where X is , or and . Here we will use the expansion around , where the odd terms vanish due to CP symmetry. In this work we evaluate the expansion coefficients up to eighth order. To obtain the Taylor coefficients, first the pressure is obtained as a function of for each value of T, then fitted to a polynomial about . All orders of derivatives are then obtained from the coefficients of the polynomial extracted from the fit. For the stability of the fit we have checked the values of least squares.
ii.3 Specific heat and speed of sound
We have studied the specific heat , which is important to find the location of CEP, the speed of sound, which determines the flow properties in heavyion reactions. The energy density is obtained from the thermodynamic potential as,
The specific heat is defined as the rate of change of energy density with temperature at constant volume, which is given by,
For a continuous phase transition near CEP, it is expected that shows a diverging behavior, which will translate into highly enhanced transverse momentum fluctuations or highly suppressed temperature fluctuations. The square of speed of sound at constant entropy is given by,
(12) 
Divergence in specific heat near CEP means vanishing speed of sound.
Iii Result
We now present the coefficients of the Taylor expansion of pressure for 2+1 flavor PNJL model with model A and model B Lagrangians and make a comparative study between the quark number susceptibility (QNS), isospin number susceptibility (INS), charge and strangeness susceptibility and their higher order derivatives. We then compare the specific heat and the speed of sound for both Lagrangians. We also compare our results with the recent lattice data available for 2+1 flavor with cheng and also with that of Polyakov extended quarkmeson model (PQM) schaefer1 ; schaefer2 ; schaefer3 ; schaefer4 .
iii.1 Coefficients of the Taylor expansion
The pressure is fitted to a polynomial in using the “gnuplot”gnu program at different values of temperature. Here we consider to take maximum eighth order term in the polynomial in . We restrict our expansion range to above which the diquark physics is expected to become important. Also the pion condensation and kaon condensation takes place in NJL model for and respectively. So we restrict our range within and below . However, above , approximate restoration of chiral symmetry implies that the chiral condensates become almost zero. So above , we have extended the range of and for the better fit of the coefficients. Near the (which is same as the least square here) of the fit varies rapidly with the variation of range of over which the fit was done. So near we have fitted the pressure for gap of temperature and the data points are spaced by of chemical potential for all temperature values. The leastsquares of all the fits came out to be or less.
Now we study the behavior of the coefficients , , and for three sets of chemical potentials for model A and model B. In figure (1) we show the variation of , , and with for for both models and lattice data. It can be seen that QNS () shows an order parameter like behavior. At low temperature there are small differences between model A and model B and model A is much closer to the lattice data. At high temperature for model A reaches almost 98% of its ideal gas value whereas for model B it reaches almost 99% of its ideal gas value. However lattice data for at high temperature reaches almost the StefanBoltzmann (SB) limit. The fourth order derivative can be thought of as the susceptibility of . The figure shows a peak near . Near the model B shows much higher peak than the model A and the peak of the eightquark interaction is closer to the lattice data. At higher temperature both cases match very well with the lattice data. But near and above lattice value converges to the SB limit, however both cases of PNJL model are slightly away from the SB limit. Note that both and have only fermionic contribution in the SB limit. Since, the coupling strength is large enough for , a sufficient amount of interaction is present in the system. So, it is expected that will not converge exactly to the SB limit within . The higher order coefficients and show interesting behavior near . Although at very low and high temperatures both of them converge to zero. Near , shows sharp peaks for both cases. However for model B the peak is much sharper. Similar behavior can be observed for , which shows more peaks near . The reason behind the peaks near the transition temperature may be due to the increase in fluctuation near . The sharper peaks of model B is probably due to the introduction of enhanced repulsive interaction through eightquark term. This was also reflected through the increase of scaled pressure in case of model B Lagrangian over the model A deb2 . The number of peaks increases near for higher order coefficients.
In figure (2) the variation of susceptibility and the higher order coefficients for the charge chemical potential is shown. The nature of all the coefficients are same as the quark chemical potential. At high temperature the fluctuation for model B is closer to the SB limit compared to the model A. However lattice data is slightly above the SB limit for . For the case of , our data (for both cases) show a better convergence towards SB limit, unlike . At low temperature the behavior of model A is closer to the lattice data compared to the model B. The quartic fluctuations show a peak near . The peak for model B is sharper than model A and the plot for model B matches well with the lattice result. The higher order coefficients show similar behavior as the quark chemical potential case.
The figure (3) shows the variation of susceptibility and the higher order coefficients for the strangeness chemical potential with . The for both the models are slightly different from the lattice data. Both the plots are almost 98% of the SB limit at high temperature, however the lattice data coincides with the SB limit. The has a similar behavior as the . However the peak is not near in this case. Near we can see a small bump for both type of Lagrangian, but the peak in both case is at higher temperature. This is due to the fact that during the chiral crossover the strange quark (the only element which carries strangeness) is sufficiently heavy and the corresponding condensate melts at much higher temperature than . The maxima of and coincide at the temperature where peaks are observed. This behavior is quite consistent with the lattice results which also indicates two peaks. However in case of lattice, the peak at is higher than the second peak. So, one can not really pin down the cause of the double peak structure. It may be a model artifact. At high temperature both the curves are above the SB limit. But the lattice data is slightly below the SB limit. The higher order derivatives also show a sharp peak at the transition temperature followed by a broader peak at higher temperature for both model A and Model B.
For the sake of completeness we have also plotted different moments of pressure for isospin chemical potential in fig. (4). It is clearly seen that, has also an order parameter like behavior like all other susceptibilities both models. At high temperature, the plot for model B reaches almost 99.5% of the SB limit, whereas the plot for model A reaches close to 98% . For , the plot for the model B shows higher peak than the case of model A at and both of them converges very well towards the SB limit at high temperature. The and shows very rapidly fluctuating peak structure near and goes to zero in both high and low temperature regime.
In figure (5) we show a comparative study for , , and with for , , and for both type of Lagrangian. In all cases we can see that the quadratic fluctuations rise rapidly in the transition region where the quartic fluctuations show a peak. This peak is most pronounced in case of . The generic form of this temperature dependence, a smooth crossover for quadratic fluctuations and a peak in quartic fluctuations, is in fact expected to occur in the vicinity of the chiral phase transition of QCD. For all of the coefficients the fluctuation in direction is strongest followed by the isospin direction fluctuations and the strangeness fluctuation. The charge fluctuation is the least pronounced. The quartic fluctuation for shows a peak at much higher temperature than the transition temperature. However other three , and show peaks at . For higher order coefficients the strangeness and charge fluctuations are negligible near compared to the quark number and isospin fluctuations.
In fig (6) we have plotted the kurtosis the ratio of for both type of potential (where ) and compared with the lattice data. The plot for for model B shows more fluctuation near than model A. However lattice data shows higher fluctuation near than the model study. At higher temperature both models coincide with the lattice data and converges well with the SB limit. Our results are qualitatively similar with PQM result schaefer3 . In case of the ratio the model B shows more fluctuation than the model A and as well as the lattice data. The model B shows almost 99% convergence with the SB limit at high temperature. For the strangeness fluctuation, we can see two peaks in curve for the both models. First peak occurs at chiral transition for light flavors and second peak occurs when chiral transition occurs in strange sector. At intermediate temperatures PNJL model overestimates the ratio than LQCD result. This feature is also observed in PQM model schaefer3 . In PQM model the ratio approaches the SB limit at high temperature and near the transition region sharp peaks appear. In our case the ratio shows sharper peak near the transition temperature and a broader peak at a higher temperature than . However lattice data shows much higher fluctuations than the plots for both models near . At high temperature the model A is closer to the SB limit than the model B.
iii.2 Specific heat and the speed of sound
We now discuss the thermodynamic quantities like specific heat () and the speed of sound (). In Fig. (7) we have plotted with for both models. From the plot we can see that grows with increasing temperature and reaches a peak at for both models. However model B shows a sharper peak at compared to the model A. Just above both the plots decrease sharply for a short range of temperature. Thereafter it gradually converges to a value slightly lower than the ideal gas value at high temperature. However the convergence towards the ideal gas value is better in case of model B. For comparison we have also plotted the values of , at which the specific heat is expected to coincide for a conformal gas. From the graph we can see that the specific heat converges very well with the at high temperature for both the models. Both the plots of and show similar behavior as PQM model schaefer3 . In PQM model grows with temperature and shows a sharp peak at . After the peak a broad bump is found around 1.2 and then goes gradually to the ideal gas value.
We now consider the speed of sound and for model A and model B in figure (7). We can see that is slightly below the ideal gas value at temperature for both cases. Our result is quite consistent with the lattice data for 2+1 flavor staggered fermions reported in szabo . A similar behavior as ours is observed in PQM model schaefer3 . We get the minimum of just below the similar to lattice data szabo and the softest point of the equation of state is found to be for model A and for model B. Model B gives better agreement with lattice data, which has its softest point of equation of state as . We have compared two sets of lattice data of Ref. cheng1 and Ref. fodor10 with our model study. Our result shows a better agreement with that of Ref. cheng1 . The softest point of equation of state of cheng1 is at , whereas the softest point of Ref. fodor10 is at much higher value . In case of PQM model the softest point of the equation of state is found to be around schaefer3 which is closer to our value.
Iv Discussion
We have studied the various fluctuations and some of the thermodynamic quantities using two diferent versions of PNJL model to understand the properties of the strongly interacting matter. In fact it is expected that the susceptibilities and the higher order fluctuations might provide the direct evidence of the order of the QCD phase transition. A pronounced peak in the susceptibility can depict the crossover transition and the sharp diverging behavior would indicate the existence of a phase transition.
We have also obtained the susceptibilities and the higher order derivatives by the Taylor expansion of pressure for two kinds of PNJL model near , where or and . In all cases the second derivative of pressure which is known as the susceptibility, show a steep rise near the transition region, which indicates near the transition region the fluctuation increases. However at higher temperature almost saturates and almost converges to the ideal gas value. This result is quite consistent with the lattice data. The higher order fluctuation shows a peak near for both models and the result matches with the lattice data. The finite height of the peak confirms the crossover nature of transition at . Both and show rapid variation around for all cases.
We have also calculated the specific heat, speed of sound for both kinds of PNJL model. The plot for , after showing a peak at , converges very well to curve at high . At high temperature almost reaches its ideal gas value and the softest point of the equation of state has a better agreement with lattice result for model B where eightquark interaction is taken into account.
In our formalism, we do not include pion condensate, kaon condensate and diquark condensate which may play an important role for higher values of chemical potential. Inclusion of those degrees of freedom may improve our result. Involved though, such studies are in progress.
V Acknowledgement
P.D. and A.L. would like to thank CSIR for financial support. A.B. thanks CSIR and UGC (UPE and DRS) for support. The authors thank Saumen Datta for some useful discussions.
References
 (1) K. Adcox et al, Nucl. Phys. A 757, 184 (2005).
 (2) G. Boyd, J. Engels, F. karsch, E. Laermann, C. legeland, M. Lugermeier and B. Peterson, Nucl. Phys. B 469, 419 (1996).
 (3) J. Engels, O. Kaczmarek, F. Karsch, and E. Laermann, Nucl. Phys. B 558, 307 (1999).
 (4) Z. Fodor, and S. D. Katz, Phys. Lett. B 534, 87 (2002).
 (5) Z. Fodor, S. D. Katz and K. K. Szabo, Phys. Lett. B 568, 73 (2003).
 (6) C. R. Allton, S. Ejiri, S. J. Hands, O. Kaczmarek, F. Karsch, E. Laermann, Ch. Schmidt, and L. Scorzato, Phys. Rev. D 66, 074507 (2002).
 (7) C. R. Allton, S. Ejiri, S. J. Hands, O. Kaczmarek, F. Karsch, E. Laermann, and Ch. Schmidt, Phys. Rev. D 68, 014507 (2003).
 (8) C. R. Allton, M. Doring, S. Ejiri, S. J. Hands, O. Kaczmarek, F. Karsch, E. Laermann, and K. Redlich, Phys. Rev. D 71, 054508 (2005).
 (9) P. de Forcrand, and O. Philipsen, Nucl. Phys. B 642, 290 (2002); Nucl. Phys. B 673, 170 (2003).
 (10) Y. Aoki, Z. Fodor, S. D. katz, and K. K. Szabo, Phys. Lett. B 643, 46 (2006).
 (11) Y. Aoki, G. Endrodi, Z. Fodor, S. D. Katz and K. K. Szabo, Nature 443, 675 (2006).
 (12) E. Megias, E. Ruiz Arriola, L. L. Salcedo, Rom. Rep. Phys. 58, 081086, (2006)
 (13) E. Megias, E. R. Arriola, L. L. Salcedo, Pos JHW2005, 025, (2006).
 (14) E. Megias, E. R. Arriola, L. L. Salcedo, Nucl. Phys. Proc. Suppl. 186, 256259, (2009).
 (15) E. Megias, E. R. Arriola, L. L. Salcedo, Phys. Rev. D 81, 096009, (2010).
 (16) K. Fukushima, Phys. Lett. B 591, 277 (2004).
 (17) C. Ratti, M. A. Thaler and W. Weise, Phys. Rev. D 73, 014019 (2006).
 (18) R. D. Pisarski, Phys. Rev. D 62, 111501 (2000); A. Dumitru and R. D. Pisarski, Phys. lett. B 504, 282 (2001), Phys. lett. B 525, 95 (2002), Phys. Rev. D 66, 096003 (2002).
 (19) K. Fukushima, Phys. Rev. D 77, 114028, (2008).
 (20) H. Hansen, W. M. Alberico, A. Beraudo, A. Molinari, M. Nardi and C. Ratti, Phys. Rev. D 75, 065004 (2007).
 (21) M. Ciminale, R. Gatto, N. D. Ippolito, G. Nardulli and M. Ruggieri, Phys. Rev. D 77, 054023 (2008).
 (22) S. K. Ghosh, T. K. Mukherjee, M. G. Mustafa, and R. Ray, Phys. Rev D 73, 114007 (2006).
 (23) S. K. Ghosh, T. K. Mukherjee, M. G. Mustafa, and R. Ray, Phys. Rev D 77, 094024 (2008).
 (24) P. Deb, A. Bhattacharyya, S. Datta and S. K. Ghosh, Phys. Rev. C 79, 055208 (2009).
 (25) A. A. Osipov, B. Hiller, J. da Providencia, Phys. Lett. B 634, 4854, (2006).
 (26) A. A. Osipov, B. Hiller, V. Bernard and A. H. Blin, Annals of Physics 321, 25042534, (2006).
 (27) A. A. Osipov, B. Hiller, A . H. Blin and J. da Providencia, Annals of Physics 322, 20212054, (2007).
 (28) B. Hiller, J. Moreira, A. A. Osipov and A. H. Blin, Phys. Rev. D 81, 116005, (2010).
 (29) K. Kashiwa, H. Kouno, T. Sakagauchi, M. Matsuzaki and M. Yahiro, Phys. Lett. B 647, 446451, (2007).
 (30) K. Kashiwa, H. Kouno, M. Matsuzaki and M. Yahiro, Phys. Lett. B 662, 2632, (2008).
 (31) A. Bhattacharyya, P. Deb, S. K. Ghosh and R. Ray, Phys. Rev. D 82, 014021 (2010)
 (32) S. A. Gottlieb et al., Phys. Rev Lett. 59, 2247 (1987).
 (33) R. V. Gavai, S. Gupta and P. Majumdar, Phys. Rev. D 65, 054506 (2002).
 (34) C. Bernard et al., Phys. Rev. D 71, 034504 (2005).
 (35) C. Bernard et al., Phys. Rev. D 77, 014503 (2008).
 (36) S. Ejiri, F. Karsch and K. Redlich, Phys. Lett. B 633, 275 (2006).
 (37) S. Roessner, C. Ratti and W. Weise, Phys. Rev. D 75, 034007, (2007).
 (38) C. Sasaki, B. Friman and K. Redlich, Phys. Rev. D 75, 074013, (2007).
 (39) S. Mukherjee, M. G. Mustafa and R. Ray, Phys. Rev. D 75, 094015 (2007).
 (40) C. Ratti, S. Roessner and W. Weise, Phys. Lett. B 649, 5760, (2007).
 (41) K. Fukushima, Phys. Rev. D 79, 074015, (2009).
 (42) M. Cheng et al., Phys. Rev. D 77, 014511, (2008).
 (43) M. Cheng et al., Phys. Rev. D 79, 074505, (2009).
 (44) Borsanyi. S et al., heplat/1007.2580v1.
 (45) W. J. Fu, Y. X. Liu and Y. L. Wu, Phys. Rev. D 81, 014028, (2010).
 (46) B. J. Schaefer and J. Wambach, Phys. Rev. D 75, 085015, (2007).
 (47) B. J. Schaefer, J. M. Pawlowski and J. Wambach, Phys. Rev. D 76, 074023, (2007).
 (48) B. J. Schaefer, M. Wagner and J. Wambach, Phys. Rev. D 81, 074013, (2010).
 (49) J. Wambach, B. J. Schaefer and M. Wagner, Acta Phys. Polon. Supp. 3 691, 2010.
 (50) V. Skokov, B. Friman, E. Nakano, K. Redlich and B. J. Schaefer, Phys. Rev. D 82 034029 (2010).
 (51) L. Stodolsky, Phys. Rev. Lett. 75, 1044 (1995).
 (52) R. Korus et al., Phys. Rev. C 64, 054908 (2001).
 (53) H. Sorge, Phys. Rev. Lett. 82, 2048 (1999).
 (54) P. F. Kolb et al., Phys. Lett. B 459, 667 (1999).
 (55) G. Boyd et al., Nucl. Phys. B 69, 419 (1996).
 (56) Y. Aoki, Z. Fodor, S. D. Katz and K. K. Szabo, JHEP, 0601, 089 (2006).
 (57) http://www.gnuplot.info/.