1 Introduction

The propagation of perturbations in fluids have been widely used to unveil the thermodynamic state of the fluid. Depending on the magnitude of the perturbation, it has been categorized as linear [1, 2] or nonlinear. If the magnitude of the perturbation is small compared to the average value of the corresponding physical quantities then it is considered as linear. On the other hand, if the perturbation is comparable to the corresponding average value then it has to be treated as nonlinear. In the context of quark-gluon plasma (QGP), formed in relativistic heavy-ion collisions (RHICs) and Large Hadron Collider (LHC), the response of QGP to the perturbations may help in characterizing the medium. The locally equilibrated system of QGP responds to its high internal pressure by rapid expansion against the vacuum. The QGP formed in the violent collisions of nuclei cools due to expansion and subsequently transforms to hadronic matter via a cross over transition or through an intermediary coexisting phase (of QGP and hadrons) depending on the value of temperature (T) and baryonic chemical potential (\(\mu \)) attained in these collisions. When these hadrons cease to interact, their momenta get frozen to certain value and finally move freely to the detector. The momentum spectra of these detected hadrons are analyzed to understand the dynamics of the collisions. In addition to the hadrons, penetrating probes, such as photons and lepton pairs emanated from each space time point of the medium, are also detected to shed light on the early hot and dense phase of QGP [3].

The pressure gradient in QGP can be anisotropic and the degree of anisotropy depends on the value of impact parameter of the collisions. An anisotropic interacting system of quarks and gluons created after the collisions will develop a pressure gradient of different magnitudes along different directions. This anisotropic pressure gradient will result in different expansion rate which in turn will induce anisotropy in momentum distribution of the hadrons originating from the system on hadronization. For the small initial spatial anisotropy, the produced momentum anisotropy could be interpreted as linear response to the initial spatial anisotropy. However, nonlinear effects will be crucial if the anisotropy in initial geometry is large due to collisions at large impact parameter. Moreover, quarks and gluons produced in the early stage of collision with high momenta, do not equilibrate but propagate as jets through the thermal medium created due to the rescattering of the low momentum quarks and gluons. These jets, while propagating, deposit energy in the medium [4, 5]. If the magnitude of the energy density resulting from the jet-medium interaction is comparable to the energy density of the medium, then the effects of jet propagation on the medium needs to be treated as nonlinear. On the other hand, smaller value of the deposited energy density allows the study of the jet-medium interaction within the domain of linear response theory. The propagation of linear perturbations in QGP has been studied extensively by several authors [6,7,8,9,10,11,12,13]. However, the nonlinear aspects of the perturbations in QGP have been studied only by a few authors [14,15,16,17,18,19,20]. In Refs. [14, 15] the authors have considered the background medium to be an ideal fluid and found that the nature of perturbations is solitonic. Fogaca et al. in Refs. [18, 19] went further and introduced the shear viscous effects in the medium and considered the perturbations up to second order. They too, found the solitonic nature of the perturbations. In these studies, the local equilibrium is treated within the scope of Gibbs–Boltzmann (GB) extensive distribution. However, particle spectra from small systems produced in \(p+p\) collision are found to be well described by non-extensive equilibrium distribution suggested by Tsallis [21,22,23,24,25,26,27,28,29,30,31,32]. In fact, the data is better reproduced [31, 33,34,35,36] with Tsallis distribution than the GB distribution.

The appearance of the entropic index, \(q\ne 1\) (also called non-extensive parameter) in Tsallis distribution is attributed to the local fluctuations in the system [37]. It is important to recall at this point that, in the limit of \(q\rightarrow 1\), the Tsallis distribution approaches the GB distribution. In principle, q can be determined from the dynamics of the system. For some cases, like non-linear Fokker–Planck equation [38], Boltzmann lattice model [39], cold atoms in optical lattices [40] etc, the q is known analytically in terms of microscopic or mesoscopic quantities. A relation between q and the basic parameters of the QCD has been established in Refs. [41, 42]. Furthermore, it was shown that, this value of q can be used to reproduce the transverse momentum spectra of hadrons produced in relativistic \(p+p\) collisions [43]. The value of q for a many body system has been expressed in terms of the number of particles (N) as, \(q=1+(2/3)N^{-1}\) and it was shown that, the two particle correlation becomes weaker as \(N\rightarrow \infty \) i.e. when the system approaches the domain of molecular chaos [44]. This nicely depicts the importance of non-extensivity of a system with finite number of particles and has relevance for system formed in high energy collisions of protons and nuclei. It has been shown in Ref. [37] that, in a system with fluctuating temperature zones, q is related to the relative variance in temperature. The interplay between the parameter q and the QCD dynamics has been studied in Refs. [45, 46]. In relativistic nuclear collisions, the q is treated as a parameter to fit the transverse momentum spectra of hadrons, which can in principle vary from 0 to \(\infty \). However, for the phenomenological studies, the values of q are generally greater than one and from the thermodynamic considerations it has an upper bound of 4/3 [47]. The descriptions of different experimental observables, which carry the effect of fluctuations, require \(q\ne 1\) as argued in Refs. [37, 48,49,50,51,52,53,54]. The value of q varies with collision centrality, energy and system size [50, 51]. It has been reported that the system size alone is not sufficient to explain the value of non-extensive parameter [55].

Moreover, there are other observables advocated for choosing the non-extensive local thermal equilibrium. The system size dependence of quantum fluctuations has been studied theoretically in Refs. [56,57,58]. These fluctuations may be represented by the q parameter. The particle spectra for small systems depends on the \(q-\)parameter. A similar situation may arise in bigger systems also as the spatial inhomogeneity created in such systems (produced in nuclear collisions) may contain some smaller zones of size similar to that of the entire system produced in a \(p+p\) collision. This indicates that the local equilibrium of such systems may be considered to be non-extensive in nature. At this stage one should find out ways to distinguish whether non-extensive local equilibrium is the underlying nature or an extensive local equilibrium is the correct situation for the fluid. This requires investigations of the effect of non-extensivity on the response of the fluid to the perturbations.

To understand the distinguishing effects between the non-extensivity end extensivity in locally equilibrated fluid, one considers two aspects in this regard: (a) the equation of state to be of a non-extensive form, and (b) the fluid dynamic equations should also be of non-extensive type, i.e., to say that the fluid velocity must also be defined by taking into account the non-extensivity, which will clearly be different from the extensive one [59].

Recently, the study of the propagation of non-linear waves has considered only the first order deviation from the extensive one, where the non-linear wave equation for first-order perturbation has been derived for the ideal fluid background with fixed thermodynamic quantities [60]. In Ref. [60] breaking wave solution has been obtained for non-dissipative propagation of the nonlinear waves. However, it has been shown in Ref. [59] that the ideal non-extensive hydrodynamics (q-hydrodynamics) is equivalent to extensive dissipative hydrodynamics (d-hydrodynamics). The corresponding dissipation should result from the non-extensive nature of fluid velocity field [59]. Therefore, even in ideal non-extensive background, one expects dissipative propagations of perturbations. Moreover, the local equilibrium leads to inherent gradients in the dynamical fields of fluid which leads to transport phenomenon. Therefore, it is important to investigate the dissipative propagation of nonlinear waves in the non-extensive background to understand the role of non-extensivity on the propagation of perturbations.

In this work, we proceed to investigate the propagation of nonlinear waves in QGP with non-extensive equation of state (EoS) and also with viscous-coefficients appearing from non-extensive nature of fluid. The aim is to find how the conclusion changes with the change in background (either extensive or non-extensive) and how q affects the propagation of the nonlinear waves. The required equations for propagation (of nonlinear waves) are derived from Müller–Israel–Stewart (MIS) hydrodynamics [61] and the order of perturbations are considered up to the second order.

The paper is organized as follows: in the next section, i.e. in Sect. 2, we derive the required equations for nonlinear waves from MIS theory with all the relevant transport coefficients. In the Sect. 2.1 we discuss the MIS relativistic causal hydrodynamics. Section 2.2 contains the equation of motion of the fluid in (1+1) dimension. In Sect. 2.3, the non-linear equations governing the propagation of the perturbations through the fluid have been provided. The Equation of State (EoS) for the non-extensive fluid has been discussed in Sect. 3. We present our results in Sect. 4 and finally summarize the findings in Sect. 5. Some of the mathematical equations and expressions are presented in Appendices A, B and C. We have used natural units (\(\hbar =c=1, k_B=1\)) and the signature of the Minkowski metric as: \(g^{\mu \nu }=(1,-1,-1,-1)\).

2 Framework of relativistic viscous hydrodynamics

In the following subsections we derive the equation of motions of non-linear perturbations propagating through a relativistic viscous fluid.

2.1 Hydrodynamic equations

The non-extensive version of ideal and first order dissipative hydrodynamics have been investigated in Refs. [59] and [62] respectively. It is found in Ref. [62] that in the q-hydrodynamics, the main tensorial structure of energy momentum tensor remains same as the d-hydrodynamics whereas, changes appear only in thermodynamic quantities and EoS. Here we use this observation and derive the nonlinear wave equation from same tensorial structure of second order MIS theory, with thermodynamic quantities and equation of state from non-extensive thermodynamics. The fluid velocity field considered here has non-extensive nature.

Relativistic viscous hydrodynamics is a useful tool to describe the space-time evolution of the fluid created in relativistic heavy-ion collisions (RHICs). Depending on the magnitude of the expansion gradient, the order of the theory is being anticipated. Navier–Stokes (NS) theory deals with first-order gradients of hydrodynamic fields, usually designated as the first-order viscous hydrodynamics [63, 64]. The solution of the NS equation is acausal and unstable. However, there are recently developed theories which shows that the first-order hydrodynamics can be stable as well as causal [65,66,67,68,69,70].

The second-order hydrodynamics takes into account the second order gradient in hydrodynamic fields and the relaxation effects were introduced to attain causality and stability. Several prescriptions of the second-order hydrodynamics can be found in Refs. [61, 71,72,73,74,75]. In this work, the second-order MIS hydrodynamics [61] is used to obtain the required evolution equations for nonlinear waves. For a fluid without conserved quantum number, the choice of Landau frame [64] is very natural. However, Landau choice does not lead to a simple behaviour in the limit of low viscosities. But in the presence of high baryon number in a system with non-zero shear and bulk viscosities, we choose Eckart’s frame of reference [63]. Eckart frame represents a local rest frame where the net charge dissipation is vanishing but the energy dissipation is non-vanishing. Although we will use Eckart frame here, one may choose the Landau frame as well. In Eckart frame, the energy-momentum tensor (EMT) and particle current can be written as:

$$\begin{aligned} T^{\mu \nu }= & {} \epsilon u^{\mu }u^{\nu }-p \Delta ^{\mu \nu }+\Delta T^{\mu \nu }, \end{aligned}$$
(1)
$$\begin{aligned} J^{\mu }= & {} nu^{\mu }, \end{aligned}$$
(2)

where, \(\epsilon \) is the local energy density field, \(u^{\mu }\) is the fluid four velocity, p is the local pressure, and \(\Delta T^{\mu \nu }\) is the correction to the EMT due to dissipation. The \(u^\mu \) satisfies the condition \(u^{\mu }u_{\mu }=1\), implying, \(u^{\mu }\partial _{\nu }u_{\mu }=0\). The projection operator normal to \(u^{\mu }\) is defined as \(\Delta ^{\mu \nu }=g^{\mu \nu }-u^{\mu }u^{\nu }\) such that, \(\Delta ^{\mu \nu }u_{\nu }=0\) and \(\Delta ^{\mu \nu }\Delta _{\mu \nu }=3\). The symmetric traceless projection normal to \(u^{\mu }\) is defined as, \(\Delta ^{\mu \nu }_{\alpha \beta }= \frac{1}{2}(\Delta ^{\nu }_{\alpha }\Delta ^{\mu }_{\beta }+\Delta ^{\mu }_{\alpha }\Delta ^{\nu }_{\beta }- \frac{2}{3}\Delta ^{\mu \nu }\Delta _{\alpha \beta })~\).

The dissipative part of \(T^{\mu \nu }\) can be decomposed in terms of scalar, vector and tensor as:

$$\begin{aligned} \Delta T^{\mu \nu }=-\Pi \Delta ^{\mu \nu }+q^{\mu }u^{\nu }+q^{\nu }u^{\mu }+\pi ^{\mu \nu }. \end{aligned}$$
(3)

The energy momentum tensor in Eckart frame then reads as:

$$\begin{aligned} T^{\mu \nu }=\epsilon u^{\mu }u^{\nu }-(p+\Pi ) \Delta ^{\mu \nu }+q^{\mu }u^{\nu }+q^{\nu }u^{\mu }+\pi ^{\mu \nu }. \end{aligned}$$
(4)

In the fluid rest frame, the vector and tensor form of dissipations are considered to be non-existent that is, \(u_{\mu }q^{\mu }=0\), \(u_{\mu }\pi ^{\mu \nu }=0\). The scalar dissipation (\(\Pi \)) corresponds to the non-equilibrium pressure perturbation, which can not be related to the energy density through EoS. The system tend to achieve a new equilibrium state in response to non-equilibrium perturbation in volume through dissipative process which is related to \(\Pi \) through bulk viscosity.

The vector \(q^{\mu }\) and the tensor \(\pi ^{\mu \nu }\) correspond to dissipative energy flow and the shear stress tensor respectively. The relations among the energy density, pressure, dissipative fluxes and EMT are given by the following relations:

$$\begin{aligned} u_{\mu }T^{\mu \nu }u_{\nu }= & {} \epsilon ,\quad q_{\alpha }=u_{\mu }T^{\mu \nu }\Delta _{\nu \alpha }=u_{\mu }\Delta T^{\mu \nu }\Delta _{\nu \alpha }, \nonumber \\ p+\Pi= & {} -\frac{1}{3}\Delta _{\mu \nu }T^{\mu \nu },\nonumber \\ \Pi= & {} -\frac{1}{3}\Delta _{\mu \nu }\Delta T^{\mu \nu },\quad u_{\mu }\Delta T^{\mu \nu }=q^{\nu }. \end{aligned}$$
(5)

The equation of motion are obtained from the conservation of EMT and the net charge density (here the net baryon number) density as:

$$\begin{aligned} \partial _{\mu }T^{\mu \nu }= & {} 0, \end{aligned}$$
(6)
$$\begin{aligned} \partial _{\mu } J^{\mu }= & {} 0. \end{aligned}$$
(7)

The dissipative fluxes [61] are given by,

$$\begin{aligned} \Pi= & {} -\frac{1}{3}\zeta \Big [\partial _{\mu }u^{\mu }+\beta _{0}D\Pi -\tilde{\alpha _{0}}\partial _{\mu }q^{\mu }\Big ], \end{aligned}$$
(8)
$$\begin{aligned} \pi ^{\lambda \mu }= & {} -2\eta \Delta ^{\lambda \mu \alpha \beta }\Big [\partial _{\alpha }u_{\beta }+\beta _{2}D\pi _{\alpha \beta }-\tilde{\alpha _{1}}\partial _{\alpha }q_{\beta }\Big ], \end{aligned}$$
(9)
$$\begin{aligned} q^{\lambda }= & {} -\kappa T\Delta ^{\lambda \mu } \Big [\frac{1}{T}\partial _\mu T +Du_{\mu }+\tilde{\beta _1} D{q_\mu }-\tilde{\alpha _0}\partial _\mu \Pi \nonumber \\&-\tilde{\alpha _1}\partial _{\nu }\pi ^{\nu }_{\mu }~ \Big ] . \end{aligned}$$
(10)

Here, \(\eta \), \(\zeta \), \(\kappa \) are the coefficient of shear viscosity, bulk viscosity and thermal conductivity respectively. \(D\equiv u^\mu \partial _\mu \), is the co-moving derivative and in the local rest frame (LRF) it represents the time derivative, \(D\Pi ={\dot{\Pi }}\). The quantities, \(\beta _0,\tilde{\beta _1}, \beta _2\) are relaxation coefficients, \(\tilde{\alpha _0}\) and \(\tilde{\alpha _1}\) are coupling coefficients, and these coefficients in the Eckart frame (\(\tilde{\alpha _{0}}, \tilde{\alpha _{1}}, \tilde{\beta _{1}}\)) are connected to the corresponding coupling and relaxation coefficients in the Landau frame (\(\alpha _{0}\), \(\alpha _{1}\), \(\beta _{1}\)) by the following relation [61]

$$\begin{aligned} \tilde{\alpha _{0}}-\alpha _{0}=\tilde{\alpha _{1}}-\alpha _{1}=-(\tilde{\beta _{1}}-\beta _{1})=-[(\epsilon +p)]^{-1}. \end{aligned}$$
(11)

The \(\beta _0, {\beta _1}\) and \(\beta _2\) are also related to the relaxation time scale as [76, 77]:

$$\begin{aligned} \tau _{\Pi }=\zeta \beta _0, \quad \tau _q=\kappa T\beta _1,\quad \tau _{\pi }=2\eta \beta _2. \end{aligned}$$
(12)

And the coupling coefficients are related to the relaxation lengths, which couple to heat flux, and bulk pressure \((l_{q\Pi }, l_{\Pi q})\), the heat flux, and shear tensor \((l_{q\pi }, l_{\pi q})\) by the following relation

$$\begin{aligned} l_{\Pi q}=\zeta \alpha _0,\quad l_{q\Pi }=\kappa T \alpha _0, \quad l_{q\pi }=\kappa T\alpha _1,\quad l_{\pi q}=2\eta \alpha _1. \nonumber \\ \end{aligned}$$
(13)

The expressions for relaxation and coupling coefficients are given in Appendix A.

2.2 1-D flow equations

In this section the relevant equations are derived to study the propagation of nonlinear perturbation in a fluid in \((1+1)\) dimension. The evolution equations of energy density and fluid velocity are obtained from the Eqs. (6) and (7) by taking projection along the directions parallel and perpendicular to \(u^{\mu }\) as:

$$\begin{aligned}&D\epsilon =-(\epsilon +p)\theta +(\partial _{\mu } u_{\nu })\Delta T^{\nu \mu }-(\partial _{\mu } q^{\mu }), \end{aligned}$$
(14)
$$\begin{aligned}&(\epsilon +p)D u^{\alpha }=\Delta ^{\alpha \mu }\partial _{\mu }-\Delta ^{\alpha }_{\nu }\partial _{\mu } \Delta T^{\nu \mu }, \end{aligned}$$
(15)

where, \(u^{\mu }=(\gamma ,\gamma v, 0,0),\quad \gamma =1/\sqrt{1-v^{2}}\) and

$$\begin{aligned} \theta= & {} \partial _{\mu }u^{\mu }=v\gamma ^3\frac{\partial v}{\partial t}+\gamma ^3\frac{\partial v}{\partial x}, \end{aligned}$$
(16)
$$\begin{aligned} \partial _{\mu }q^{\mu }= & {} q^{x}\frac{\partial v}{\partial t}+v\frac{\partial q^{x}}{\partial t}+\frac{\partial q^{x}}{\partial x}. \end{aligned}$$
(17)

Here,

$$\begin{aligned}&u_{\mu }q^{\nu }=0,\quad u_{\mu }\pi ^{\mu \nu }=0,\quad g_{\mu \nu }\pi ^{\mu \nu }=0, \nonumber \\&\pi ^{xy}=\pi ^{zx}=\pi ^{yz}=0,\quad \pi ^{yy}=\pi ^{zz}=\pi ^{T} \end{aligned}$$
(18)

which leads to,

$$\begin{aligned} \pi ^{00}= & {} v^2\pi ^{xx},\quad \pi ^{xt}=\pi ^{tx}= v^2\pi ^{xx},\nonumber \\ \pi ^{T}= & {} \pi ^{tx}=\frac{( v^2-1)}{2}\pi ^{xx}, \end{aligned}$$
(19)

and

$$\begin{aligned} q^{t}=vq^{x}, \quad q^{y}=q^{z}= q^{T}. \end{aligned}$$
(20)

Now from the energy-momentum conservation equations i.e. from Eqs. (14) and (15), we get,

$$\begin{aligned}&\gamma v\Big (\frac{\partial \epsilon }{\partial t}+\frac{\partial \epsilon }{\partial x}\Big )\nonumber \\&\quad = \Big [ (h+\Pi )v\gamma ^3+v \gamma \pi ^{xx}+q^{x}\}\frac{\partial v}{\partial t} \Big ]\nonumber \\&\qquad +\Big [(h+\Pi )\gamma ^3+\gamma \pi ^{xx}\Big ] \frac{\partial v}{\partial x}+v\frac{\partial q^{x}}{\partial t}+\frac{\partial q^{x}}{\partial x}, \end{aligned}$$
(21)
$$\begin{aligned}&\Big [ (h+\Pi )v\gamma ^4+v \gamma ^2(1-2v^2)v \pi ^{xx}+\gamma ^3(1+v^2)v q^{x}\Big ]\frac{\partial v}{\partial t}\nonumber \\&\qquad +\Big [ (h+\Pi )v\gamma ^4-\gamma ^2 v^2\pi ^{xx} +\gamma v(1+2\gamma ^2)q^{x}\Big ] \frac{\partial v}{\partial x}\nonumber \\&\quad =-v^2\Big [\gamma ^2\frac{\partial p_{\zeta }}{\partial t}+\frac{\gamma ^2 }{v}\frac{\partial p_{\zeta }}{\partial x} +\frac{\partial \pi ^{xx}}{\partial t}+ \frac{1}{v}\frac{\partial \pi ^{xx}}{\partial x}\nonumber \\&\qquad +\frac{ \gamma }{v}\frac{\partial q^{x}}{\partial t}- \gamma \frac{\partial q^{x}}{\partial x}\Big ] \end{aligned}$$
(22)

where \(h=\epsilon +p\) is the enthalpy density. We define \(p+\Pi =p_{\zeta }\). The equation governing the net baryon number conservation reads,

$$\begin{aligned} \frac{\partial n}{\partial t}+ v\frac{\partial n}{\partial x} = -n \gamma ^2\Big [v\frac{\partial v}{\partial t}+\frac{\partial v}{\partial x}\Big ], \end{aligned}$$
(23)

The other three equations originating from dissipative fluxes are given by the Eqs. (8), (9), and (10), which can be written as:

$$\begin{aligned}&\Pi +\zeta \beta _{0}\Big [\gamma \frac{\partial \Pi }{\partial t}+v\frac{\partial \Pi }{\partial x}\Big ] =-\zeta \Big [\{v\gamma ^3 -\tilde{\alpha _{0}} q^{x}\}\frac{\partial v}{\partial t}\nonumber \\&\qquad +\{ \gamma ^3\} \frac{\partial v}{\partial x}-\tilde{\alpha _{0}}v\frac{\partial q^{x}}{\partial t}-\tilde{\alpha _{0}}\frac{\partial q^{x}}{\partial x}\Big ], \end{aligned}$$
(24)
$$\begin{aligned}&q^{x}+\kappa T \tilde{\beta _{1}}\Big [\gamma \frac{\partial q^{x}}{\partial t}+v\frac{\partial q^{x}}{\partial x}\Big ]\nonumber \\&\quad = \kappa T\Big [ \{-\gamma ^4+ \tilde{\beta _{1}}\gamma ^3 vq^{x}+\tilde{\alpha _{1}}(2-\gamma ^2) \pi ^{xx} \}\frac{\partial v}{\partial t}+\{-v\gamma ^4\nonumber \\&\qquad + \tilde{\beta _{1}}\gamma ^3 v^2 q^{x}- \tilde{\alpha _{1}}\gamma ^2 v \pi ^{xx} \} \frac{\partial v}{\partial x} +\tilde{\alpha _{1}}v \frac{\partial \pi ^{xx}}{\partial t}+\tilde{\alpha _{1}} \frac{\partial \pi ^{xx}}{\partial x}\nonumber \\&\qquad -\tilde{\alpha _{0}}v\gamma ^2 \frac{\partial \Pi }{\partial t}-\tilde{\alpha _{0}}\gamma ^2 \frac{\partial \Pi }{\partial x} +\frac{1}{T}v\gamma ^2 \frac{\partial T}{\partial t}+\frac{1}{T}\gamma ^2 \frac{\partial T}{\partial x}\Big ], \end{aligned}$$
(25)
$$\begin{aligned}&\pi ^{xx}+\frac{4}{3}\eta \beta _{2}\Big [\gamma \frac{\partial \pi ^{xx}}{\partial t}+v\frac{\partial \pi ^{xx}}{\partial x}\Big ]\nonumber \\&\quad = \frac{4}{3}\eta \Big [ \{\gamma ^5+\tilde{\alpha _{1}}\gamma ^4 v q^{x}+ 2\beta _{2}\gamma ^3 \pi ^{xx} \}\frac{\partial v}{\partial t}+\{\gamma ^5+ \tilde{\alpha _{1}}\gamma ^4 v q^{x}\nonumber \\&\qquad + 2\beta _{2}\gamma ^3v^2 \pi ^{xx} \} \frac{\partial v}{\partial x}-\tilde{\alpha _{1}}\gamma ^2 v\frac{\partial q^{x}}{\partial t}-\tilde{\alpha _{1}}\gamma ^2 \frac{\partial q^{x}}{\partial x}\Big ]. \end{aligned}$$
(26)

2.3 Derivation of nonlinear wave equations

The above hydrodynamic equations from MIS causal theory is used to derive the nonlinear equations governing the motion of the perturbations in the fluid. We have adapted Reductive Perturbative Method (RPM) [78,79,80,81]. To proceed further we need to define ‘stretched co-ordinates’ as,

$$\begin{aligned} X=\frac{\sigma ^{1/2}}{L}(x-c_s t), \quad \text {and} \quad Y=\frac{\sigma ^{3/2}}{L}(c_s t), \end{aligned}$$
(27)

where, L is the characteristic length and \(c_{s}\) is the speed of sound. Therefore, from the above equations we get

$$\begin{aligned} \frac{\partial }{\partial x}=\frac{\sigma ^{1/2}}{L}\frac{\partial }{\partial X},\quad \text {and} \quad \frac{\partial }{\partial t}=-c_{s}\frac{\sigma ^{1/2}}{L}\frac{\partial }{\partial X}+c_{s}\frac{\sigma ^{3/2}}{L}\frac{\partial }{\partial Y},\nonumber \\ \end{aligned}$$
(28)

where \(\sigma \) is the expansion parameter. The coordinate X is measured from the frame of propagating sound waves, whereas Y represents fast moving coordinate. The RPM technique is devised to preserve the structural form of the parent equation. Now we perform the simultaneous series expansion of hydrodynamic quantities in powers of \(\sigma \) to get,

$$\begin{aligned} {\hat{\epsilon }}= & {} \frac{\epsilon }{\epsilon _0}=1+\sigma \epsilon _1+\sigma ^2 \epsilon _2+\sigma ^3 \epsilon _3+\cdots ,\nonumber \\ {\hat{p}}= & {} \frac{p}{p_0}=1+\sigma p_1+\sigma ^2 p_2+\sigma ^3 p_3+\cdots ,\nonumber \\ {\hat{v}}= & {} \frac{v}{c_s}=\sigma v_1+\sigma ^2 v_2+\sigma ^3 v_3+\cdots ,\nonumber \\ {\hat{\Pi }}= & {} \frac{\Pi }{p_0}=\sigma \Pi _1+\sigma ^2 \Pi _2+\sigma ^3 \Pi _3+\cdots ,\nonumber \\ {\hat{q}}^x= & {} \frac{q^x}{\epsilon _0}=\sigma q^x_1+\sigma ^2 q^x_2+\sigma ^3 q^x_3+\cdots ,\nonumber \\ {\hat{\pi }}^{xx}= & {} \frac{\pi ^{xx}}{p_0}=\sigma \pi ^{xx}_1+\sigma ^2 \pi ^{xx}_2+\sigma ^3 \pi ^{xx}_3+\cdots , \end{aligned}$$
(29)

where \(\epsilon _0\) and \(p_0\) are the background energy density and pressure respectively (given by Eqs. (36) and (37) below) on which the perturbations propagate. By collecting terms with different order of \(\sigma \) from the perturbation series, different types of equation like Breaking wave equation, Burgers equation or KdV equation, etc can be obtained [2, 18, 60].

We rewrite the equations of motion (Eqs. (21), (22),(23), (24), (25),  (26)) using Eqs. (28), (29) and collect terms with different orders in \(\sigma \) and finally revert from \((X,Y) \rightarrow (t,x)\) to attain the evolution equations of different order of perturbations. It is found in Ref. [19] that for conformal background, the evolution equation of the first-order perturbation does not include the relaxation and coupling coefficients arising from the MIS theory. Therefore, it is required to go to the second order to find the role of these coefficients. To obtain the time evolution of second order perturbation, we need to collect terms up to third order in \(\sigma \) in the expansion, because the equations in n-th order in \(\sigma \) contains time derivatives of \({\hat{\epsilon }}_{n-1}\). Therefore, to get equations containing time derivative of \({\hat{\epsilon }}_{2}\) we go up to 3rd order in \(\sigma \). This leads to the following equations for the perturbation in the energy density (\({\hat{\epsilon }}\)) as:

$$\begin{aligned}&\frac{\partial \hat{\epsilon _{1}}}{\partial t} +\left[ 1+(1-c_{s}^2) \frac{\epsilon _{0}}{\epsilon _{0}+p_{0}}\hat{\epsilon _{1}}\right] c_{s}\frac{\partial \hat{\epsilon _{1}}}{\partial x}\nonumber \\&\qquad - \left[ \frac{1}{2(\epsilon _{0}+p_{0})}(\zeta +\frac{4}{3}\eta )\right] \frac{\partial ^{2} \hat{\epsilon _{1}}}{\partial x^{2}}=0, \end{aligned}$$
(30)

and

$$\begin{aligned} \frac{\partial \hat{\epsilon _{2}}}{\partial t}{+}{\mathcal {S}}_{1}\frac{\partial \hat{\epsilon _{2}}}{\partial x}{+}{\mathcal {S}}_{2}\frac{\partial \hat{\epsilon _{1}}}{\partial x}{+}{\mathcal {S}}_{3} \frac{\partial ^{2} \hat{\epsilon _{1}}}{\partial x^{2}}{+}{\mathcal {S}}_{4}\frac{\partial ^{3} \hat{\epsilon _{1}}}{\partial x^{3}}{+}{\mathcal {S}}_{5}\frac{\partial ^{2} \hat{\epsilon _{2}}}{\partial x^{2}}{=}0. \nonumber \\ \end{aligned}$$
(31)

where \(\hat{\epsilon _{i}}= \sigma ^i \epsilon _i\) for \(i=1,2,\ldots \). The coefficients \({\mathcal {S}}_i\)’s for \(i=1\) to 5 are given in Appendix B. Equations (30) and (31) have been solved to investigate the fate of the nonlinear waves propagating through a relativistic viscous and non-extensive fluid. The relevant EoS (i.e., relation between \(\epsilon _0\) and \(p_0\)) and the initial conditions required to solve these equations are discussed in Sects. 3 and 4 respectively.

3 Tsallis MIT bag equation of state

We assume that the quarks and gluons in the viscous plasma medium follow the quantum Tsallis Fermionic (F) and Tsallis Bosonic (B) single particle distributions given by [82],

$$\begin{aligned} f_{Fq}=\frac{1}{\left[ 1+(q-1)\frac{E_p-\mu }{T}\right] ^{\frac{q}{q-1}}+1}\nonumber \\ f_{Bq}=\frac{1}{\left[ 1+(q-1)\frac{E_p-\mu }{T}\right] ^{\frac{q}{q-1}}-1} , \end{aligned}$$
(32)

where, \(E_p=\sqrt{\vec {p}^2+m^2}\) is the single particle energy of a particle of mass m, \(\mu \) is the (baryonic) chemical potential (assumed to be zero in the present work), q is the Tsallis (non-extensivity) parameter and T is the Tsallis temperature. The \(\epsilon \) and p of a system of quarks and gluons in MIT bag model [83] is given by,

$$\begin{aligned} \epsilon _{\mathrm {}}= & {} {\mathcal {B}}+\epsilon _{B}+2\epsilon _F, \end{aligned}$$
(33)
$$\begin{aligned} p= & {} -{\mathcal {B}}+P_{B}+2P_F, \end{aligned}$$
(34)

where \({\mathcal {B}}\) is the bag constant (\({\mathcal {B}}^{1/4}= 200\) MeV). The \(\epsilon _{i}\) and \(P_{i}\) are given by,

$$\begin{aligned} \epsilon _i = g \int \frac{d^3p}{(2\pi )^3}E_p f_{iq},~ \ \ \ \ P_i = g \int \frac{d^3p}{(2\pi )^3}\frac{p^2}{3E_p} f_{iq}, \end{aligned}$$
(35)

where g is the degeneracy factor and \(i=B, F\).

The factor 2 in front of the Fermionic contributions in Eqs. (33) and (34) takes into account both quarks and anti-quarks. The net baryonic number density, n is assumed to be zero here. The Fermionic thermodynamic variables have been evaluated up to \({\mathcal {O}}(m^2T^2)\). It is assumed that the quark-gluon plasma consists of the up quark, down quark and gluons. It has been verified that for a wide range of the values of q and T relevant for phenomenological studies of high-energy collisions, the \({\mathcal {O}}(m^2T^2)\) approximated results overlap with the exact results when the mass is small (\(\sim \) 10 MeV). We write down the expressions for the energy density and pressure in the non-extensive bag model [84, 85] as:

$$\begin{aligned} \epsilon _{\mathrm {}}= & {} {\mathcal {B}} + \epsilon _{B,1}T^4 + 2 (\epsilon _{F,1}T^4+\epsilon _{F,2}m^2T^2)\nonumber \\= & {} {\mathcal {B}} + \epsilon _{\alpha }T^4 + \epsilon _{\beta }m^2T^2 , \end{aligned}$$
(36)
$$\begin{aligned} p_{\mathrm {}}= & {} -{\mathcal {B}} + P_{B,1}T^4 + 2 (P_{F,1}T^4+P_{F,2}m^2T^2) \nonumber \\= & {} -{\mathcal {B}} + P_{\alpha }T^4 + P_{\beta }m^2T^2. \end{aligned}$$
(37)

where,

$$\begin{aligned} \epsilon _{\alpha }= & {} \epsilon _{B,1} + 2\epsilon _{F,1}, ~P_{\alpha } = P_{B,1} + 2P_{F,1}, \nonumber \\ ~\epsilon _{\beta }= & {} 2\epsilon _{F,2},~P_{\beta } = 2P_{F,2}. \end{aligned}$$
(38)

The expressions for \(\epsilon _{i,\ell }\) and \(P_{i,\ell }\) (\(i=F, B\); \(\ell =1,2\)) are given in the Appendix C.

The expressions for pressure and energy density obtained in this section are used for getting \(\epsilon _0\) and \(p_0\) to solve the Eqs. (30) and (31), governing the evolution of the nonlinear perturbations.

4 Results and discussion

We present here the results on the fate of the nonlinear perturbations in the QGP fluid at zero baryonic chemical potential which leads to the vanishing thermal conductivity of the medium. Since, the obtained equations look like KdV equation (ideal background) whose analytic solution is a sec-hyperbolic function with a shifted argument. For better understanding the deviation of propagation from that of the ideal background, the initial profile of perturbations of both the orders are taken as of the same form as the sech function. We solve the equations describing the evolution of \(\hat{\epsilon _i}\) for the following initial profile of the perturbations:

$$\begin{aligned} {\hat{\epsilon }}_1= & {} A_{1}\bigg [\mathrm {sech} \frac{(x-x_{0})}{B_{1}}\bigg ]^{2} \nonumber \\ {\hat{\epsilon }}_2= & {} A_{2}\bigg [\mathrm {sech} \frac{(x-x_{0})}{B_{2}}\bigg ]^{2}, \end{aligned}$$
(39)

where \(x_0\) is the initial position of the peak of the perturbations. Since, here we are considering 1D propagation in \(+x\) direction in an extended and static background, any choice of value of \(x_0\) will provide the same salient features of the propagation. We have taken \(x_0=10\) fm here to demonstrate the propagation of the perturbation however, one may choose any other value of \(x_0\) as well. The solutions of Eqs. (30) and (31) possess solitonic nature (with sech\(^2\) dependence) similar to the solutions of KdV equations for small dissipations and this motivates us to consider initial profile given by Eq. (39). It may be mentioned here that for small spatial gradient the solitonic behaviour of the solution has also been demonstrated in Ref. [19] for Gaussian initial profiles in a conformally invariant extensive hydrodynamic background.

Fig. 1
figure 1

a \({\hat{\epsilon }}_1\) as a function of x at different time for \(A_{1}=0.65\) and \(B_{1}=1\) fm. b \({\hat{\epsilon }}_2\) as a function of x at different time for \(A_{2}=0.35\) and \(B_{2}=1\) fm. c \({\hat{\epsilon }}=1+{\hat{\epsilon }}_1+{\hat{\epsilon }}_2\) as a function of x at different time for \(A_{1}=0.65, B_{1}=1\) fm, \(A_{2}=0.35\) \(B_{2}=1\) fm. The background temperature, \(T=200\) MeV and \(q=1.08\) here

Full size image

Firstly, we investigate the effects of non-extensivity only through the EoS, keeping the value of the shear viscosity at the KSS bound [86], \(\eta /s=1/4\pi \). The bulk viscosity is taken as \(\zeta /s= 1/4\pi \) and \(\kappa = 0\). The value of T and q are taken as \(T=200\) MeV and \(q=1.08\). Figure 1a, b show the time evolution of the spatial shape of the first-order \((\hat{\epsilon _1})\) and the second-order \((\hat{\epsilon _2})\) perturbations respectively for the initial profiles given by the Eq. (39). In Fig. 1a, the values amplitude \((A_1)\) and width \((B_1)\) of the first-order perturbation are taken as \(A=0.65\), \(B=1\) fm (Fig. 1a), whereas the corresponding values for the second-order perturbation are taken as \(A_2=0.35\) and \(B_2=1\) fm (Fig. 1b). Comparison of results depicted in Fig. 1a, b indicate that the dissipation of the second-order perturbation is less as compared to first order. This can be understood from the Eq. (31), which indicates that \(\hat{\epsilon _2}\) gets contributions from \(\hat{\epsilon _1}\) which contains terms with third-order derivatives in space, known as dispersive term, responsible for preserving the shape of the propagating nonlinear wave. Figure 1c shows the evolution of the energy density up to second order scaled by the background energy density. Though the QGP system produced in heavy-ion collisions has a lifetime of the order of 10–15 fm, we have considered the time evolution up to 80 fm to have a better understanding of the late time behaviour. It is evident from the results that with time, amplitude decreases and width increases indicating towards the possibility of complete dissipation, unlike the case studied earlier in Ref. [60]. It may be noted that up to 15 fm the dissipation of the nonlinear wave is not very drastic. The change in amplitude is faster in the beginning than later times. At later time the propagation speed becomes smaller. The slowing down of the dissipation can be understood from a comparison of the results displayed in Figs. 1c, 2, 3a, and Fig. 3b. In Figs. 1c and 2, the amplitudes of the initial perturbations are taken to be different with the same value of the width, whereas in Fig. 1c, the amplitude is taken as \(A=2\) and in Fig. 2 it is considered to be \(A=4\). It is observed that with larger amplitude, the speed of propagation is larger, which is clearly seen in Fig. 2 (the result with \(t=80\) fm provides better visibility). This is due to amplitude dependent propagation of nonlinear waves, a well known feature of the nonlinear propagation. The nonlinear waves become slower as the amplitude reduces due to dissipation.

Fig. 2
figure 2

\({\hat{\epsilon }}\) as a function of x at different time for \(T= 200\) MeV, \(q= 1.08\), \(A_{1}=1.65, B_{1}=1, A_{2}=1.35, B_{2}=1\) (parameters are same as in Fig. 1c but with higher amplitude). We observe that with larger amplitude the speed is larger

Full size image
Fig. 3
figure 3

a \({\hat{\epsilon }}\) as a function of x at different time for \(T= 200\) MeV, \(q= 1.08\) \(A_{1}=0.65, B_{1}=2, A_{2}=0.35\) and \(B_{2}=2\) (parameters are same as in Fig. 1c) but with higher width). We observe that the wave propagates faster with lower dissipation as the width increases. b Same as Fig. 1c with \(A_{1}=0.65, B_{1}=0.5, A_{2}=0.35,\) and \(B_{2}=0.5\). In this case larger dissipation and smaller speed of propagation is observed

Full size image

In Fig. 3a, b, we compare the results for the widths taken as \(B=2\) fm and \(B=0.5\) fm respectively with the same value of the amplitude (\(A=2\)). It is observed that the increase in width of the initial perturbation causes less dissipation with the progression of time for the same background fluid conditions. The propagation speed is found to be greater for perturbation with larger width. This is because, the perturbation with larger width dissipates less, and thus due to amplitude dependent speed of propagation, it propagates faster in space. Comparison of results presented in Figs. 2 and 3 indicates that the increase in amplitude or in width reduces the damping but increases the speed of the peak position.

Now we move to discuss the role of the background of the fluid, inscribed through q and T, on the perturbations. We take the same values of \(A_1\), \(A_2\), \(B_1\) and \(B_2\) as in Fig. 1c to study the fate of the perturbations after \(t=10\) fm for three different values of the background temperatures and fixed value of \(q=1.08\). The result is displayed in Fig. 4. We find that the dissipation decreases with increase in the background temperature but the perturbation in the front becomes steep. However, for ideal background, the magnitude as well as the tilting remain unchanged with the variation of the background temperature. In the present work we consider the background energy density as a function of temperature in contrast to previous work [60] on ideal hydrodynamic background at fixed background energy density with varying temperature.

Fig. 4
figure 4

\({\hat{\epsilon }}\) as a function of x at time \(t= 10\) fm for different T. Other parameters are \(q= 1.08\), \(A_{1}=0.65, B_{1}=1, A_{2}=0.35, B_{2}=1\). We observe smaller dissipation at larger temperature

Full size image
Fig. 5
figure 5

\({\hat{\epsilon }}\) as a function of x at time \(t= 10\) fm and \(T= 200\) MeV for different values of q. The amplitudes and the widths of the perturbations are defined by Eq. (39) with \(A_{1}=0.65, B_{1}=1, A_{2}=0.35, B_{2}=1\). Smaller dissipation and tilting along x axis have been seen for larger q

Full size image
Fig. 6
figure 6

\({\hat{\epsilon }}\) as a function of x for different time for ideal fluid (\(\eta = \zeta = 0\)). We take \(T=200\) MeV, \(q=1.08\), \(A_{1}=0.65, B_{1}=1, A_{2}=0.35\) and \(B_{2}=1\). We observe that the waves tend to loose their localization with time

Full size image
Fig. 7
figure 7

The variation of \({\hat{\epsilon }}\) with spatial coordinate (x) at \(t=2\) fm for different values of q and viscous coefficients have been shown here. The temperature of the fluid, \(T=200\) MeV. The values of \(A_{1}, A_{2}\), \(B_{1}\), and \(B_{2}\) are same as Fig. 1c

Full size image

Figure 5 shows the fate of the same initial perturbation as that of displayed in Fig. 1c after a time of 10 fm for three different values of q for \(T=200\) MeV. In the presence of dissipation the q plays an interesting role. Namely for higher values of q the dissipation is lower and the stiffness is more. In the absence of dissipation (ideal background) we also find more stiffness for higher q [60]. Indicating that the effects of dissipation can be reduced by increasing q, or in other words systems with higher q imitate lower dissipation. It is also to be noted from the results shown in Fig. 5 that the propagation speed of the peak increases with q if the q dependence is considered only through EoS. The results of ideal q-hydrodynamics can be obtained from our results by setting the values of the transport coefficients to zero. In Fig. 6 the results for vanishingly small transport coefficients have been shown. In such cases we find that the waves loose their localization (breaking waves) and the shape preserving nature no longer survives.

So far we have discussed the effects of non-extensivity on the propagation of nonlinear waves only through the EoS. Since the EoS and transport coefficients control the bulk evolution of a system, therefore, for a thorough treatment of the propagation of nonlinear waves in a non-extensive background (q-background), one should consider the effects of non-extensivity, not only through the EoS (q-EoS) but also on the transport coefficients (q-viscosities) which are estimated for q-background (see Refs. [62, 87] for details). The results with the q-viscosities (q-dependence of viscous coefficients) along with the q-EoS (as discussed in Sect. 3) are presented in Fig. 7 (purple and blue dashed lines) for the initial profile shown by green line. The results obtained with the effects of q-EoS only is compared (black and red dashed lines). It is found that the dissipation of nonlinear waves is more when non-extensivity is introduced via the q-EoS and the q-viscosities together. Interestingly, here, the propagation speed of the peak of the perturbations decrease with the increasing value of q when the effects of q-viscosity is introduced. Along with the reduction of the dissipation the higher values of q also reduces the speed of the propagation. This means that the introduction of non-extensivity favours the non-dissipative (or anti-dissipative) slowing down of nonlinear perturbations. Considering the fact that for an extensive background, the slowing down is caused by the dissipative effects [20], the results displayed in Fig. 7 reveal that the propagation of non-linear wave can be useful to distinguish between the extensive and the non-extensive background. We also find that the dissipation is more for smaller q-values (compare with Fig. 5). With the q-viscosities along with q-EoS, the relative dissipation with q is more than that of the q-EoS only. This is in contrast to the expectation from the notion that the non-extensive distribution is treated as a non-equilibrium deviation from the extensive one. In such a scenario the dissipative correction is interpreted as a result of \(q\ne 1\). These findings are opposite to the normal expectations. Rather, it indicates that q-equilibrium should not always be interpreted as a dissipative non-equilibrium form of extensive statistics.

In a system with q-background, the higher values of q is accompanied by less dissipation in the system. This can be understood from the fact that the enhanced correlation in a system can be represented by higher q-values [44, 88]. The correlations enhance the fluctuations, which can also be represented by a value of \(q\ne 1\) [37]. The correlations opposes the entropy generation in a system and hence dissipation. Therefore, a higher value of q parameter representing a higher degree of correlation which is expected to results in less dissipation, this is manifested in the variation of viscosity to entropy ratio (\(\eta /s\)) with q [62]. The \(\eta /s\) decreases with increasing q due to greater rate of change of non-extensive entropy with q than of the bulk and shear viscosities, thereby reducing the dissipation.

The second and higher order space derivatives in the evolution equations contain the effects of dissipation. Without the presence of the viscous terms, the q alone (through q-EoS only) can not produce such dissipative nature, i.e., q can only have a role in dissipation in the presence of explicit dissipative terms. That is why, without considering the q-viscosities (or any dissipative effects arising from d-hydrodynamics), the authors in Ref. [60] do not find any dissipation with q-EoS alone. Also they have considered only the first order of energy perturbation and found that the equation to be of the breaking-wave type, causing a tilting (ultimately breaking) of nonlinear waves without any dissipation. Whereas, in this work we have obtained solution with behaviour similar to the solution of KdV-type equation by considering perturbations up to second-order with the inclusion of all the relevant transport coefficients and found to produce prominent dissipative effects without much tilting.

In such a scenario, one may be tempted to think that, perhaps, less dissipation is the unique signature of presence of finite q (increasing q-values). However, it may be noted that the less dissipation may be accounted through the smaller values of the transport coefficients also. So there may be an ambiguity behind the observed smaller dissipation in a system, whether it is a non-extensive system with larger values of q or it is an extensive system with lower viscosity. Such ambiguity can be removed through other effects which can distinguish these two separate situations. The role of higher values of q as the cause of less dissipation can be mimicked by reducing the values of the transport coefficients, but it may possess other distinct effects which will help in removing the ambiguity. In the case of nonlinear wave propagation, the nature of the nonlinear wave will provide the relevant distinct feature. With less dissipation in extensive background (i.e. without q-background) the nonlinear wave will not preserve the solitonic nature due to breaking of the wave (see Fig. 6), whereas the q makes the dissipation less with preserving the solitonic nature (see Fig. 7). Therefore, the causes for less dissipation with Tsallis background can be distinguished with the propagation of a nonlinear wave. This may pave the way for investigation to ascertain whether the viscosity of the QGP is really low or it has q-thermalized local equilibrium with not so low transport coefficients. Since \(p+p\) collision admits q-distributed particle spectra, therefore, one would argue for consideration of q-thermalized local equilibrium. In that case it is very important to remove the ambiguity about viscosities of QGP estimated by fitting the experimental results. Another important point to be noted is that, when the ideal q-hydrodynamics is interpreted as equivalent to d-hydrodynamics, one should not take q as a cause of dissipation in general sense. This is so because, q does not necessarily represents dissipation, especially when the background is in q equilibrium.

5 Summary and conclusion

In summary, we have derived the equations for the propagation of nonlinear waves in a dissipative fluid involving shear viscosity, bulk viscosity and thermal conductivity. The effect of the \(q-\)statistics is incorporated by taking into account the \(q-\)EoS and q-viscosities. We have obtained the equations of motion describing the perturbation up to second order in energy density. The solution of the second order equation has solitonic character similar to the solution of KdV type equation. It is observed that the dissipation of the wave is more in first order than the second order. For larger amplitude of the perturbations, we find that the speed of the propagation is larger, whereas for broader perturbations, the dissipation is comparatively less and propagation is faster. As far as the background temperature of the medium is concerned, dissipation is less for higher temperature. The effect of q is evident from the results shown in Fig. 5, where we see more dissipation of nonlinear waves for smaller values of q. This is in line with the interpretation of q as representative of correlations as well as fluctuations in a system. The larger values of q may imitate the situation with smaller values of the transport coefficients. This is reflected through the maintenance of solitonic nature of the nonlinear waves. Moreover, it is found that the non-extensivity favours the non-dissipative slowing down of nonlinear perturbations. This is in contrast to the extensive background where the slowing down comes with the dissipation [20]. This may be a useful criterion for making distinction between the extensive and non-extensive background. So the nature of propagation of a nonlinear wave can be potentially used to distinguish a non-extensive medium from an extensive one, which could be especially helpful in situations where momentum distribution of particles can not distinguish the nature of underlying statistics (Gibbs–Boltzmann or Tsallis). This indistinguishability may occur when there is an uncertainty on whether the equilibration is complete or not and when there are other sources of power law distributions such as the radiation from particles with high momentum or jets moving through the QGP, instead of the presence of correlations that causes non-extensivity of the medium.