Nonextensive quasiparticle description of QCD matter

The dynamics of QCD matter is often described using effective mean field (MF) models based on Boltzmann-Gibbs (BG) extensive statistics. However, such matter is normally produced in small packets and in violent collisions where the usual conditions justifying the use of BG statistics are not fulfilled and the systems produced are not extensive. This can be accounted for either by enriching the original dynamics or by replacing the BG statistics by its nonextensive counterpart described by a nonextensivity parameter $q\neq 1$ (for $q \to 1$ one returns to the extensive situation). In this work we investigate the interplay between the effects of dynamics and nonextensivity. Since the complexity of the nonextensive MF models prevents their simple visualization, we instead use some simple quasi-particle description of QCD matter in which the interaction is modelled phenomenologically by some effective fugacities, $z$. Embedding such a model in a nonextensive environment allows for a well-defined separation of the dynamics (represented by $z$) and the nonextensivity (represented by $q$) and a better understanding of their relationship.


Introduction
Dense hadronic matter is usually described using relativistic mean field (MF) theory models (like, for example, the Walecka model for nucleons [1][2][3] or the Nambu-Jona-Lasinio model (NJL) for quarks [4][5][6][7]). All of them use the Boltzmann-Gibbs (BG) statistics, which means that they assume a homogeneous and infinite heat bath and in their original versions they do not account for any intrinsic fluctuations or long range correlations. However, this kind of matter is typically produced in violent collision processes and in rather small packets, which rapidly evolve in a highly nonhomogeneous way and whose spatial configurations (like the correlations between quarks located in different nucleons in NJL models) remain far from being uniform (in fact, there is no global equilibrium established, cf. [8][9][10] and references therein). As a result, some quantities become non-extensive and develop power-law tailed rather than exponential distrubutions, making application of the usual BG statistics questionable (cf., [11,12] and references therein). The remedy is either to supplement the BG statistics by some additional dynamical input or, when it is not known, to use some form of nonextensive statistics generalizing the BG one, for example Tsallis statistics [13,14]. The latter is characterized by a nonextensivity parameter q = 1 (for q = 1 one recovers the usual BG statistics). In fact, such an approach has already been investigated some time ago and the q-versions of essentially all types of MF models were formulated (see [15][16][17][18] and references therein). In the meantime the validity of the nonextensive q-thermodynamics used in such cases was also confirmed [19][20][21][22] and the conditions for its thermodynamical consistency were established [23][24][25][26].
In the nonextensive approach one investigates the way in which some selected observables change when one departs from the extensive statistics with value of q = 1. The goal is to disclose how, and to what extent, these changes are correlated with the possible modifications of the dynamics governing the model considered or with the possible influence of some external factors caused by the surroundings in which formation of dense QCD matter takes place and which is not accounted for in the usual extensive approach. In fact, it is expected that when these factors are gradually identified and their impact is accounted for by a suitable modification of the original model, the value of |q − 1| obtained from comparison with experiment gradually diminishes and |q − 1| = 0 signals that our improved dynamical model fully reproduces all aspects of the process considered [27].
The investigation of the interplay between these two factors is the subject of our work. However, in the case of the MF models such a procedure is not transparent because of the complexity of the dynamics of MF models (for example, as shown in our nonextensive NJL model [16], particles acquire dynamical masses which implicitly depend on the nonextensivity parameter). This prevents a clear interpretation of the role played by the parameter q and its interplay with the dynamics. There is therefore a need to simplify the dynamics, for example by reducing it to a number of well defined (temperature dependent) parameters. Such a possibility is offered by quasi-particle models (QPM) in which the interacting particles (quarks and gluons) are replaced by free quasi-particles. They can be formulated in a number of ways, the most popular approaches are: the model encoding the interaction in the effective masses [28,29], the model using the Polyakov loop concept [30,31] and the model based on the Landau theory of Fermi liquids where the effects of the interaction are modelled by some temperature dependent factors called effective fugacities, z (i) (T ), which distort the original Bose-Eistein or Fermi-Dirac distributions [32][33][34][35][36][37]. We will continue to use this model, and call it the z-QPM (note that there are also quite a number of other works on the QPM, cf., for example, [38][39][40][41][42]). This choice is motivated by the fact that in z-QPM the masses of quasi-particles are not modified by the interaction (they do not depend on the fugacities z (i) ) what allows us to avoid problems encountered in other approaches.
In the z-QPM the effective fugacities z (i) (T ) ( z (i) ≤ 1, the z (i) = 1 correspond to a noninteracting gas of gluons and quarks) are obtained from fits to lattice QCD results [43][44][45] which serve as a kind of experimental input [32][33][34][35]. Note that the effective fugacities have nothing to do with the usually used fugacities corresponding to the observation of particle number and are therefore not related to the chemical potential; they just encode the effects of the interactions between quarks and gluons. Because there are problems with allowing for a nonzero chemical potential in lattice simulations [46,47], the z-QPM was initially formulated assuming a vanishing chemical potential, µ = 0. Starting from [36] a small amount of non-vanishing chemical potential µ was introduced in the matter sector (to reproduce a realistic equation of state of the QGP), and assumed to be a constant whose value varies between µ = 0 and 100 MeV (depending on the circumstances, but such that µ/T << 1). However, so far in all fits to lattice data used by the z-QPM the chemical potential µ was neglected.
The aim of this paper (which is an extension of our previous work [48]) is two-fold. Firstly, after embedding the z-QPM in a nonextensive environment characterised by a nonextensive parameter q, we investigate the qz-QPM created in this way in terms of the changes in the effective fugacities, z (i) → z (i) q , necessary to fit the same lattice data. Secondly, we use our qz-QPM but retain the same effective fugacities z (i) as in the z-QMP model and calculate the changes in the densities and pressure induced only by the changes in the nonextensivity q. This parallels, in a sense, our nonextensive q-NJL model [16] with its dynamics replaced by a phenomenological parametrization in terms of fugacities z. However, unlike in the q-NJL model, in both cases our investigations are limited to T above the critical temperature T c because only such are considered in lattice simulations.
Please note that, in terms of dynamics, we do not introduce here any new model. We have just adapted for our purposes the widely known z-QPM [32][33][34][35], accepting its physical motivation which, when combined with its transparency and simplicity, makes this model especially useful for our purposes. However, this also means that the conclusions of this work have, at most, the same level of credibility as those of the z-QPM.
The paper is organized as follows. In Section 2 we provide a short reminder of z-QPM. Section 3 contains a formulation of the qz-QPM. Our results are presented in Section 4 and Section 5 concludes and summarises our work. Technical details are placed in the Appendices A, B, C and D.

A short reminder of the z-QPM
We start with a short reminder of the z-QPM proposed and used in [32][33][34][35][36][37]. It is based on the following effective equilibrium distribution function for quasi-partons (i = q, s, g for, respectively, u and d quarks, strange quarks and gluons): Here e(x) = exp(x), ξ = +1 for bosons and −1 for fermions and β = 1/T . In the z-QPM u and d quarks are assumed massless, E i=q = p, and strange quarks have mass m, E s = m 2 + p 2 ; for gluons E g = p. The z (i) ≤ 1 denote the effective fugacity describing the interactions, they are assumed to depend only on the scaled temperature, τ = T /T c where T c is the temperature of transition to the deconfined phase of QCD). The dynamics described by the lattice QCD data is encoded in z (i) . For z (i) = 1 one has free particles. The appearance of a chemical potential needs some comment. In the equation of state the fugacity z, which is connected with the interactions between particles, changes the pressure P and is therefore connected with the change of the chemical potential µ. It reflects the evolution of the system from some initial state, described by µ 0 and P 0 , to a state described by µ and P with ∆(µ) = µ − µ 0 = T ln (P/P 0 ), which can be derived from the equation of state for constant temperature T . For a noninteracting gas where the relative pressure (P/P 0 ) → 1, this correction vanishes, ln(P/P 0 ) → 0. In the z-QPM [32][33][34][35][36][37] one considers a gas of quarks and gluons above the critical temperature, T > T c , and assumes a quasi-particle description of the lattice QCD equation of state, which in the limit of high temperature (T → ∞) is given by a noninteracting gas of quarks and gluons. The correction ∆(µ) is replaced here by the fugacity z = exp[−∆(µ)] multiplying distribution function. By analogy to a perfect gas the effective pressure becomes unity in limit of the large T and z(T → ∞) → 1. Consequently, in the isothermal evolution of a hadron gas for finite temperatures, the chemical potential, or a single particle energy, are corrected by ∆(µ) = T ln(z) . Note that whereas usually the chemical potential µ enters together with the energy E, cf. Eq. (2), it can also be associated with the fugacity x(τ ) modifying it by an exponential, temperature dependent, factor: The effective fugacity,z (i) q , obtained this way combines the action of the original effective fugacity and that of the chemical potential.
Some remarks concerning the way of the effective fugacities are obtained from the lattice data used in z-QPM [32][33][34][35][36][37] are in order here. The QCD thermodynamics at high temperature can be described in terms of a grand canonical ensemble which can be expressed in terms of the distribution functions which, in turn, depend on the fugacities, cf. Eq. (1). One of the most important quantities calculated on the lattice is pressure. The pressures of the gluons and quarks (expressed as functions of the fugacities) were therefore compared with the corresponding pressures obtained from the lattice data; in this way one gets effective fugacities as functions of scaled temperature, z(τ ) (τ = T /T cr with T cr being the critical temperature). Because it turns out there is no single universal functional form describing the lattice QCD data over the whole range of τ , the low and high τ domains were therefore described by different functional forms with the cross-over points at τ g = 1.68 for gluons and τ q = 1.7 for quarks and were chosen as: They were then used to describe the QCD lattice data [43][44][45] with the parameters listed in Table 1.
To formulate the qz-QPM one has to replace the previous extensive effective distribution function for quasi-partons by its nonextensive equivalent, where and Thermodynamical consistency demands that the n q (x) obtained in this way must be replaced by [25,26,49] (this requirement follows from the proper theoretical formulation of the nonextensive thermodynamics provided in [15,17], cf. Eqs(10) below). A comment on the conditions of validity of the qz-QPM is in order here. The tacit assumption of the z-QPM is that both x and (−x) remain positive, i.e., that z (i) (τ ) ≤ 1 [32]. However, immersing our system in a nonextensive environment means that some part of the dynamics is now modelled by the parameter q, therefore the above constraints are not sufficient because e q (x) and e 2−q (x) must always be nonnegative real valued and the allowed range of x is given by the condition that [1+(q −1)x] ≥ 0 which must be satisfied and which can limit the available phase space [16]. Referring for details to [16,49] we say only that out of three possibilities of introducing nonextensivity discussed in [16], only two (one for particles and one for antiparticles) limiting appropriately the available phase space are applicable for our purpose. The third method, which does not limit the available phase space (and which was discussed in detail in [17]), introduces some novel dynamical effects, not observed in dense nuclear matter; therefore we shall not use it here [16]. Both the form of n q (x) and the fact that it effectively emerges as n q q (x) can be derived from the formulation of the nonextensive thermodynamics in which one starts from the nonextensive partition function Ξ q (the meaning of the index i and the parameter ξ is the same as in Eqs. (1), (2) and (5)) taken as [15,17] (V denotes the volume, i = g, q, s for, respectively, gluons, light quarks (u and d) and strange quarks, and ν i are the corresponding degeneracy factors which we take the same as in [32]: ν g = 16, ν q = 24 and ν s = 12): Integrating by parts, and noting that one arrives at the following alternative expression for the nonextensive partition function, with effective distribution functions equal now [n q (x)] q . A note of caution is necessary here. After closer inspection one realizes that the definition of e q (x) used in [15], when used together with the duality relation (7), leads to [n 2−q ] 2−q in Eq. (10), instead of n q q presented in [15] (cf., their Eq. (35)). The nonextensive versions of the particle density ρ q and the energy density, ε q , are defined, respectively, as (we use Eqs. (9), (37), (39) and (40) Note that the energy density in our qz-QPM depends explicitly on the nonextensivity via the nonextensive particle density and implicitly via a possible q-dependence of the effective fugacities mentioned previously. In the extensive limit, q → 1, Eq. (12) becomes equal to the corresponding equation (10) from the z-QPM [32]. The physical significance of the effective nonextensive fugacities is best seen when looking at the corresponding nonextensive dispersion relations defined as (cf., Eq. (12)) Note that the masses of the quasiparticles remain intact and the single quasiparticle energies ω (i) q are modified only by the action of the effective fugacities, z (i) q (T ). In both extensive and nonextensive cases this results in some additional contributions to the quasiparticle energies which can be interpreted as coming from the collective excitations. They occur because of the temperature dependence of the effective fugacities (deduced from the lattice calculations) which can be interpreted as representing the action of the gap equation in [16] taken at constant energy E i .
Because in the version of q-thermodynamics used here all thermodynamic relations are preserved, the pressure P q is given by the usual thermodynamic relation, From Eqs. (12) and (14) one gets an expression for the trace anomaly (β = 1/T )

Results
We shall now calculate the nonextensive effective fugacities, z q (τ ), which, for a given value of the nonextensivity parameter q, reproduce the original z-QPM results [32]. These results, in turn, were obtained from a comparison with the lattice QCD simulations from [43] using relation (14) to match the pressures in the z-QPM and in the lattice QCD simulations. Note that this procedure assumes in fact that the trace anomaly in the z-QPM (Eq. (15) with q = 1) is the same as that resulting from the QCD lattice data [32]. We adopt the same procedure and use Eq. (14) to match the pressures calculated, respectively, for q = 1 (as in in [32]) and for q = 1, (it is tacitly assumed that in both extensive and nonextensive environments the temperature T remains the same). This means that in our case the trace anomaly remains the same as in the z-QPM (and as in the lattice data) and does not depend on the nonextensivity. To this end the following conditions must be satisfied: for gluons (with ν g = 16) and for quarks. They give us the τ and q-dependent relations between the extensive fugacities obtained in [32]), z (i) (τ ) (which are our input), and the nonextensive fugacities, z (i) q (τ ) (which are our results). The function Θ(q, q) defines the allowed phase space; its details are presented in Appendix A.
Figs. 1 shows the resulting effective fugacities, z q = z (q) q (τ ) and z g = z (g) q (τ ), as functions of the scaled temperature, τ = T /T c . They can be fitted using the same parametrization as before, i.e., Eq. (4), with the parameters displayed in Table II. Since the values of z (i) q=1 obtained in [32] were obtained assuming µ = 0, the same assumption was used in obtaining our z (i) q here. Note that for the nonextensivites q used here the changes in the fugacities are small, q (τ ) and z (g) q (τ ) as a function of the scaled temperature τ = T /T c (calculated for µ = 0). Lower panels: As above but shown in more detail and with an enlarged range of the nonextensivity parameter q. and can be approximated (with very good accuracy of a few percent) by (cf. Appendix B), . (20) Eq. (20), together with Fig. 1, allows for a better understanding of interrelation between the dynamics (represented by the fugacities z) and the nonextensivity described by q. The central point is that all the z q must describe the lattice QCD data (directly for q = 1 in z-QPM and indirectly for q = 1 in qz-QPM, where they are forced to reproduce the results of z-QPM). The τ -dependence of z(τ ) starts from small values (corresponding to strong attraction) towards z = 1 (corresponding to free, noninteracting particles). The case of z > 1 would formally mean the emergence of repulsive forces and is not allowed in z-QPM, therefore we shall also keep this limitation in our qz-QPM. The replacement of extensive media by not extensive means adding some repulsive interaction (in the case of q < 1) or an attractive one (for q > 1). Therefore, in the first case it must be compensated by an increase in z (i.e., δz q > 0) and in the second case by a decrease (i.e., δz q < 0). Note now that whereas in the latter case we have z q (τ ) < z q=1 (τ ) < 1, in the former there is limiting value of τ = τ lim (q), depending on q, for which z q (τ lim ) = 1. This means that for τ > τ lim (q) the attraction represented by z(τ ) is already too weak to compensate the repulsion introduced by q < 1. The value of τ lim diminishes with the increase Table 2. Numerical values of the coefficients a (i) , b (i) , a ′ (i) and b ′ (i) (i = q, g) in Eq. (4) when used for different values of q resulting in the curves displayed in the lower panels of Fig. 1. of this repulsion (i.e., with the increase of |q − 1|). Not wanting to introduce the problem of repulsion we limit our considerations to τ > τ lim only. So far, results for z(τ ) and z q (τ ) have been obtained with µ = 0. The formal introduction of the chemical potential µ in z-QPM [36] makes z-QPM more flexible and applicable to possible future lattice QCD data with the chemical potential accounted for. Following this new development in z-QPM we have also formally introduced µ into our z-QPM. We can therefore check what would be the value of our z q in the case when part of the dynamics is shifted from fugacity z to the chemical potential µ. Eq.
(3) shows the effective fugacity with the chemical potential included. It is visualized in Fig. 2 where we plot a number of results for different values of the chemical potential µ and for two values of the nonextensivity parameter: q = 0.99 and q = 1.01. As one can see, nonzero µ diminishes the real values   21) and (22)), respectively) in the nonextensive environment as functions of the scaled temperature, τ = T /T c for q = 0.96, 1, 1.04, calculated for µ = 0. Right panel: The same as above but shown in more detail and with an enlarged range of the nonextensivity parameter q. The introduction of the chemical potential µ also changes the q-dependence of the relative density of the quarks, R ρ = ρq ρ q=1 , where ρ q is given be Eq. (11). As can be seen in the left panel of Fig. 3, in a nonextensive environment one observes a clear separation of the situations with R ρ > 1 and R ρ < 1. The first occurs for q < 1 and the observed increase of density is consistent with lowering of the entropy which, in turn, is connected with the tighter packing of the quarks in this case [16]. The second occurs for q > 1 and the picture is reversed; it is consistent with an increase of the entropy and with looser packing of the quarks in this case. Note that this behaviour of R ρ = R ρ (q) is fully consistent with the behaviour of the nonextensive fugacities presented in Fig. 1. Essentially the same result can be obtained using the linear approximation of the n q q in (q − 1) as given by Eqs. (56) and (57). Using now the same values of z q but adding some amount of the chemical potential µ yields the results shown in the right panel of Fig. 3. We observe some increase of the relative density with µ with a possible trace of a small upper bending.
We shall now calculate the modifications of partonic charges in a hot QCD medium embedded in a nonextensive environment calculating the corresponding Debye mass, M qD . Following [32] we use for the extensive Debye mass the expression derived in semiclassical transport theory in which M D is given in terms of equilibrium parton distribution functions (N c denotes the number of colors): In the nonextensive environment described by the nonextensivity parameter q we simply replace In Fig. 4, following [32], we present the ratio of M D /M I D where M I D denotes the Debye mass for the ideal EOS case (i.e., with z g = 1 and z q = 1) which, following [32], equals Note that because the Debye mass is essentially a combination of the densities of quarks and gluons the above results resemble those for the effective fugacities and all previous remarks also apply here. We now proceed to the second part of our work in which we keep the original dynamics of the qz-QPM intact using the same effective fugacities as in [32] (i.e., we assume that z q (i) → z (i) as given by Eq. (4) with the parameters listed in Table I). This parallels to some extent our approach in the nonextensive q-NJL model [16] (but now the dynamics is simplified and represented by the temperature dependent fugacities z(τ ) reproducing the lattice QCD results) and allows us to investigate the sensitivity of the selected observables to the nonextensive environment only. The only drawback is the limitations in the temperatures allowed because the fugacities are only defined for τ > 1, i.e., above the critical temperature T c .
To start with we present in Fig. 5 the corresponding relative pressure P q /P q=1 and relative density ρ q /ρ q=1 as functions of the nonextensivity parameter q at fixed temperature (left panel) and their dependencies on T for some fixed nonextensivity q (right panel). Note that whereas before the pressure was assumed to be the same for extensive and nonextensive environments, P q /P q=1 = 1, it now increases linearly with q in the same way as in the q-NJL model [16]. The relative density also increases with the nonextensivity q, contrary to its previous behaviour (demonstrated in the left panel of Fig. 3) where it decreased with q.
We now proceed to the trace anomaly, T q , Eq. (15). Note that now it acquires some explicit dependence on the nonextensive parameter q. Using the definitions of energy density, ε q , Eq. (12) and pressure P q (Eq. (14) (and additionally Eqs. (10) and (2)), we obtain that Figure 5. (Color online) Dependencies of relative pressure P q /P q=1 and density ρ q /ρ q=1 on the nonextensivity parameter q at fixed temperature (left panel) and on the temperature T for fixed q (right panel). The change in the trace anomaly generated by the nonextensivity q is given by (see Eqs. (56) and (57))    21)) and (22) as a function of τ for some selected values of q (left panel) and the same but for some selected values of the chemical potential µ (right panel). The dynamics is the same as in the z-QPM, i.e., z q = z q=1 = z. Unlike the results presented in Fig. 4, this time they are caused solely by the action of the nonextensive environment with z q = z q=1 = z as used in the z-QPM.

Summary and conclusions
In this work we investigate the interrelation between nonextensive statistics and the effects of dynamics in dense QCD matter. It continues our previous analysis of this problem on the example of the nonextensive version of the NJL model, the q-NJL. However, the complexity of its dynamics does not allow for a clear separation of the purely dynamical effects from the nonextensive ones. Therefore, in this work, following some specific quasi-particle models (z-QPM) [32][33][34][35][36][37], we used simplified dynamics reduced to a number of well defined parameters, the effective fugacities, z ∈ (0, 1). In this kind of QPM the masses of the quasi-particles are not modified by the interaction, which enables the problems and inconsistencies encountered in other approaches to be avoided. The fugacities z increase with temperature T from very small values in the vicinity of the critical temperature, T cr (which corresponds to strong interactions between quarks and gluons), towards unity (which corresponds to a free gas of quarks and gluons). They modify only the argument of the exponent in the corresponding Bose-Einstein or Fermi-Dirac distributions: e(x) → e(x − ln z). The action of nonextensivity is different, it changes the functional form of the exponent, e(x) → e q (x), leaving the argument x unchanged. This means that the actions of nonextensivity and dynamics are complementary and cannot be replaced by each other (although sometimes they describe the same, or comparable, situations). The fugacity z therefore models phenomenologically the dynamics of the mean field theory in the extensive environment and does not account for intrinsic correlations and fluctuations present in the system, while these are most naturally described phenomenologically by the nonextensivity q. Phenomenologically, both approaches nicely complement each other in what concerns the description of the dense QCD system. If we wanted to replace the action of nonextensivity, q, by the respective action of dynamics, z, (or vice versa) then either z or q would have to acquire energy dependence, which we consider as untenable.
Note that, contrary to the q-NJL model [16], the qz-QPM model is formulated in such a way as to reproduce the effective fugacities of the original z-QPM [32], which, in turn, describes the lattice QCD results [46,47]. This means that the qz-QPM also describes them; in fact they serve as a kind of experimental data. Such constraints were not present in the q-NJL model. Therefore, our conclusions are more reliable than those presented in [16]. The interplay between dynamics and nonextensivity is best seen in Fig. 3 (left panel) which shows results for the relative densities, R ρ = ρ q /ρ q=1 , in the nonextensive environment. For q < 1 (which corresponds to a lowering of the entropy) one observes R ρ > 1, which can be interpreted as caused by some positive (attractive) correlations in the system and may be connected with a tighter packing of quarks. The opposite is observed for q > 1 (corresponding to an increase of the entropy) where R ρ < 1. This can be interpreted as resulting from the repulsion of the quarks and fluctuations developing in the system. Both of these correlations and fluctuations are imposed on the effects of the interaction described by the fugacities z q . This is the clearest example of dynamical effects introduced by the nonextensive environment and characterized by the nonextensivity parameter q.
Let us now look more closely at the results on z q (τ ) presented in Fig. 1. Note that for q < 1, when, according to the left panel of Fig. 3 our system becomes more dense, one observes that δz q = z q<1 − z q=1 > 0 and increases with |q − 1| (i.e., increases with density). This means that the interaction represented by z q becomes weaker. As a result, the upper limit of z q = 1 (corresponding to a noninteracting gas of quarks and gluons) is reached for smaller temperature T , the more so the bigger |q − 1| (i.e., the smaller q). This means that to obtain the same pressure in the system one needs a weaker interaction described by fugacity; the increasing part of it is caused by the effect of the nonextensivity q. In other words: the change of statistics from extensive (q = 1) to nonextensive with q < 1, allows the attainment of the limit of the ideal gas with weaker correlations between quarks and gluons caused by the fugacity z. For the q > 1 case our system becomes, according to the left-panel of Fig. 3, less dense; the correction term needed to obtain the same pressure as in the extensive case is now negative, δz q = z q>1 − z q=1 < 0, and |δz q | grows only very slowly with increasing q (i.e., with decreasing density) becoming constant for higher T ; the limit z q = 1 is never reached for finite temperature T . This is because for q > 1 one expects some intrinsic fluctuations (for example temperature T fluctuations) which work against the dynamical interactions represented by z. Therefore, these interactions cannot cease and z q cannot grow too fast. In fact, with increasing T they seem to become constant and one observes a kind of equilibrium between dynamics and nonextensivity.
Because the z-QPM [32][33][34][35][36][37] uses the lattice QCD results [43][44][45] as its input and because there are problems with nonzero chemical potential µ in the lattice calculations [46,47], the z-QPM was initially formulated for zero chemical potential, µ = 0, which substantially limits its applications. However, anticipating the possibility of the emergence of some new lattice QCD results with the chemical potential included (if only partially), starting from [36] some small amount of non-vanishing µ in the matter sector was introduced. We have therefore also allowed for some nonvanishing µ. In Fig. 2 we show how nonzero µ influences the extracted z q for q < 1 and q > 1. Fig. 3 (right panel) shows that the relative density (both for q < 1 and q > 1) increases (almost) linearly with the chemical potential. Note that the possible introduction of the chemical potential in the lattice QCD calculations will change profoundly the z-QPM (and the qz-QPM); it will therefore become a third phenomenological parameter modelling the interaction. Our results shows in what direction these changes will proceed and in Appendix D we provide a scheme of expansion of the pressure in the chemical potential to allow for the possible further application of our qz-QPM should similar results occur in the lattice calculations [46,47]. Fig. 4 presents the results for the Debye mass in a nonextensive environment. Note that because it is essentially a combination of densities of quarks and gluons the results therefore resemble those for the effective fugacities. Calculations of more involved quantities, like, for example, dissipative effects would be much more involved because they would demand the use of the nonextensive version of the transport or hydrodynamic equations, which is beyond the scope of this work and will be presented elsewhere. It would also be desirable to be able to compare directly the results of qz-QPM with some future nonextensive lattice QCD simulations, which seems to be gaining some interest recently [50,51].
Finally, Figs. 5 -7 present some selected results on the, respectively, relative pressure and density, trace anomaly and Debye mass obtained when we use the qz-QPM with the same effective fugacities z (i) as in the z-QMP and calculate changes in densities and pressure induced only by changes in the nonextensivity q. Because in this case the original dynamics represented by z (i) remains intact, all changes in the results are caused only by the nonextensivity, i.e., by the fact that |q−1| = 0. These results correspond, in a sense, to the results obtained in our q-NJL model [16], with the proviso that now our investigations are limited to T above the critical temperature T c (i.e., to τ = T /T c > 1 because only such are considered in z-QMP). In both models we observe similar dependencies of the pressure and density on the nonextensivity parameter q while maintaining all dynamical parameters for a given temperature T the same. They are reduced for q < 1 and enhanced for q > 1. This means that when changing the amount of nonextensivity one cannot keep the same pressure P in the system without changing the dynamical parameters (or their temperature dependencies). Our qz-QPM is therefore a simple example of such changes needed to achieve equalization of the pressure in extensive and nonextensive systems. Note that now, as a result of the pressure equalization in extensive and nonextensive systems, the relative densities, R ρ = ρ q /ρ q=1 , change with q in opposite ways, becoming higher for q < 1 and lower for q < 1 than the density in an extensive system. This observation could be important for some new version of the q-NJL model, in which one could insist on keeping the same pressure for different nonextensivities and looking for the corresponding changes in its dynamical parameters (which would become q-dependent).
Such an approach could have its further application in investigations of the EoS of dense matter.
The authors declare no conflict of interest.

A. Limitations of the allowed phase space in the nonextensive approach
The functions Θ(p; i) (i = g for gluons, i = q for light quarks and q = s for strange quarks) provide the limitations of the allowed phase space resulting from the condition In the q < 1 case, for gluons (with zero mass and zero chemical potential) we have that Because z (g) q < 1, our integral is non-vanishing (i.e., p > 0) only for Stronger interactions (corresponding to smaller values of the fugacity) are in this case not allowed for the q used here. In the case of quarks (with chemical potential µ > 0 and with mass m for strange quarks) condition (27) results in the following limitation Now the phase space is open if: In the first case the (µ − m) choice is more restrictive and results in the condition that which for µ = 0 and m = 0 coincides with the corresponding condition for gluons. In the second case the choice (µ + m) is the more restrictive, for which For nonzero mass (strange quarks) it is more restrictive than condition (31).
In the case of q > 1 we have for gluons that In our case it is always satisfied and there are no limitations on z However, for z g < 1 it does not introduce any further limitations.

B. Approximate calculation of z q
Let us denote z q = z + δ (where z = z(τ ) are the fugacities obtained in [32] from lattice QCD and δ = z q − z is the change in fugacity emerging from the nonextensive environment). We shall now calculate δ/z for the case of small δ, |δ/z| << 1. We start by expanding L q (x) from Eq. (9), in δ and keeping only linear terms: Denoting one can write that (cf., Eqs. (6) and (5)) obtaining (note that for δ = 0 x q → x − ln z) The integrals of the type presented in Eqs. (17) and (18) can therefore be rewritten as integrals over Because z q is obtained from the condition that ∆L q = 0, in the first approximation the correction term δ is equal to Θ(p) represents the possible limitation of the phase space caused by the nonextensivity (cf. Appendix A). It depends on the nonextensivity parameter q and on the type of particle considered (gluons, light quarks or strange quarks). Formula (43) can be further approximated by expanding it in q − 1 and retaining only the linear terms in (q − 1). Following Eqs. (65) and (66) one obtains that Similarly, following Eqs. (56) and (57), one has that n q q ≃ n(x; z) + (q − 1)n(x; z)∆ q [n(x; x)] where ∆ q [n(x; z)] = ln n(x; z) + Therefore where B(x; z) = ξ 2 ∞ 0 dpp 2 ∆I(x; z)Θ(p), C(x; z) = ∞ 0 dpp 2 n(x; z)Θ(p), D(x; z) = ∞ 0 dpp 2 Θ(p)n(x; z) · ln n(x; z) + 1 2 [1 + ξn(x; z)]x 2 .
More involved expressions.  1 We do not address the question of the applicability of such an approach, assuming its validity for the range of variables used here (cf. [52]).
The corresponding q-logarithm and (2 − q)-logarithm functions to be used in what follows are connected with the q-exponential function e q (x) and its dual e 2−q (x): =⇒ ln X and ln 2−q X = ln q [e q (X)] = X, and ln 2−q [e 2−q (X)] = X, From them one gets that: Finally, the generalization of the relation n q (x) + n 2−q (−x) = 1 to the case where the effective particle densities are given not by n q but by n q q is approximately given by D. Expansion of pressure in chemical potential µ In the case when we allow for a chemical potential µ, in some applications we need to know the expansion of the pressure P (as given by Eq. (14)) in the chemical potential µ (in fact inμ = µ/T = βµ < 1). We present below the two first terms of such an expansion, where ln q (Ξ q ) is given by Eq. (9) and ∂ ln q (Ξ q ) ∂µ We have used here Eq. (11) and Eqs. (68) -(73) (with n q [x q (µ)], e q [x q (µ)] and x q (µ) are defined by Eqs. (6) and (5)) and {n q [x q (µ)]} q = N q (µ)): ∂N q ∂n q = qn q−1 q ,