Nonlinear coupling of flow harmonics: Hexagonal flow and beyond
Higher Fourier harmonics of anisotropic flow ( and beyond) get large contributions induced by elliptic and triangular flow through nonlinear response. We present a general framework of nonlinear hydrodynamic response which encompasses the existing one, and allows to take into account the mutual correlation between the nonlinear couplings affecting Fourier harmonics of any order. Using Large Hadron Collider data on Pb+Pb collisions at TeV, we perform an application of our formalism to hexagonal flow, , a coefficient affected by several nonlinear contributions which are of the same order of magnitude. We obtain the first experimental measure of the coefficient , which couples to and . This is achieved by putting together the information from several analyses: event-plane correlations, symmetric cumulants, as well as new higher-order moments recently analyzed by the ALICE collaboration. The value of extracted from data is in fair agreement with hydrodynamic calculations, although with large error bars, which would be dramatically reduced by a dedicated analysis. We argue that within our formalism the nonlinear structure of a given higher harmonic can be determined more accurately than the harmonic itself, and we emphasize potential applications to future measurements of and .
Anisotropic flow () in heavy-ion collisions (1) has been measured up to the sixth Fourier harmonic, (2); (3); (4), and preliminary results on were recently reported (5). In ultra-central collisions, is to a good extent determined by linear response to the initial-state anisotropy in the harmonic (6); (7). In less central collisions, however, higher-order harmonics () get important contributions induced by and , through non-linear couplings (8). The magnitude of these non-linear couplings is encoded in the so-called response coefficients, which are largely insensitive to the initial state and directly probe the hydrodynamic behavior (9); (10). As a consequence, a generic prediction of hydrodynamics is that these coefficients depend weakly on both the collision centrality and the details of hydrodynamic calculations (11); (12); (13). Therefore, nonlinear response coefficients are robust probes of hydrodynamic behavior: Any disagreement between the calculated values and experimental data cannot easily be fixed via a tuning of the parameters.
While it can be easily argued that there is only one leading nonlinear contribution to and , several nonlinear couplings need to be considered in the decomposition of harmonics (14) and higher. Existing theoretical (9); (13); (10); (15); (16) and experimental (4) analyses of hexagonal flow isolate the various contributions by assuming that they are pairwise uncorrelated. For instance, they neglect the modest event-plane correlation between elliptic flow and triangular flow, which is measured (17). In this article, we improve the existing formalism by relaxing this assumption. We show in Sec. II that even if the nonlinear terms are strongly correlated, they can still be separated by means of a simple matrix inversion. In Sec. III, we explain how the corresponding matrix elements, which are moments (18), can be obtained from existing data. The values of the nonlinear response coefficients involving obtained from experimental data are presented in Sec. IV, and they are compared to simple hydrodynamic calculations in Sec. V. Eventually, in Sec. VI we stress the importance of using our formalism in the characterization of the nonlinear structure of harmonics beyond hexagonal flow, and .
Ii Improved formalism of nonlinear coupling
Let us first recall how the nonlinear coupling is defined in the simple case of quadrangular flow, (9).
In a hydrodynamic calculation, anisotropic flow is given in every event by , where curly brackets denote an average value taken with the single-particle distribution at freeze-out (19); (20).
The transformation of under an azimuthal rotation is .
In this way, and both get the same factor , so that azimuthal symmetry allows a coupling between and .
Therefore, one can separate into a contribution proportional to , and a remaining part, which we dub
where the nonlinear response coefficient is the same for all events in a centrality class. This decomposition is uniquely determined if one imposes the condition that the two components and are uncorrelated, that is, , where angular brackets denote an average over events in the centrality class. This condition, together with Eq. (1), implies
This equation defines uniquely the response coefficient , and was recently employed in experimental analyses (4). In an experiment, though, is not measured in every event due to finite multiplicity fluctuations, but the averaged quantities appearing in Eq. (2) can be measured accurately (18); (21). If is proportional to in every event, then as defined by Eq. (2) reconstructs the proportionality coefficient, even if fluctuates event to event (9). In this case, vanishes, so that can generally be interpreted as the part of which is not induced by .
Similarly, can be decomposed as
which is also measured (4).
As for hexagonal flow, , azimuthal symmetry allows several nonlinear terms (13), so that it can be written as:
We write the third nonlinear term as rather than so as to avoid double counting with the term in the case where is proportional to . Note, however, that the structure of the decomposition is unchanged if one replaces with . Using Eq. (1), one can indeed rewrite Eq. (5) as:
We show now that the nonlinear response coefficients appearing in Eq. (5) are uniquely determined as soon as one imposes that all the nonlinear terms are uncorrelated with the last term, . In previous works (9); (13); (10); (15); (4), though, such construction was supplemented by a stronger assumption, namely, that all terms in the right-hand side of Eq. (5) are pairwise independent. As we shall see, this turns out to be a reasonable approximation, but an unnecessary one.
To simplify the notation, let us rewrite a decomposition such as (6), in the generic form
where is the number of nonlinear terms ( in the case of ), are products of lower-order harmonics (, , for Eq.(6)), and denote the corresponding coupling constants. As done for , we define in Eq. (7) by the condition that it is linearly uncorrelated with all the nonlinear contributions:
The left-hand side is a moment involving the higher harmonic , while the terms are a set of moments involving lower-order harmonics. Note that since each is itself nonlinear, these moments are at least of order 4. Eq. (9) is a linear system of equations for the coupling constants, .
We define a matrix, , by:
It is hermitian by construction, and real if one neglects parity violation (22). We denote by the -vector whose components are the moments , and by the -vector whose components are the response coefficients . With these notations, the system (9) can be rewritten in matrix form:
It is solved by inverting the matrix:
Assuming that all nonlinear terms are mutually independent, as done in previous works, amounts to assuming that is diagonal. As we shall see in Sec. III, the off-diagonal elements of can all be extracted from existing data, hence, Eq. (12) can be applied directly without any approximation.
Note that these equations involves the higher harmonic only through the moments which are linear in the higher harmonic. By constrast, standard measures of , say, , are quadratic. Since the magnitude of decreases rapidly with the order , observables linear in a high-order harmonic are typically measured more accurately than quadratic observables. Hence, the nonlinear couplings of higher-order harmonics can be determined more precisely than these harmonics themselves (4); (5).
Iii Extracting the matrix elements from experimental data
Note that the diagonal elements of are moments involving only the magnitude of anisotropic flow, . On the other hand, the off-diagonal elements of , as well as the components of , involve relative phases between different Fourier harmonics, and are related to the so-called event-plane correlations (24); (17).
In principle, all these moments can be measured using the same experimental setup involving two subevents separated by a rapidity gap (18).
This method has recently been implemented by the ALICE collaboration (4).
The values of , , , , can be directly obtained from these data through simple algebraic manipulations.
For the remaining moments, we need to combine information from different analyses. could be extracted from the quantity dubbed in the first ALICE analysis of triangular flow (25). We instead choose to obtain it through the event-plane correlation measured by the ATLAS collaboration (17) via
where the left-hand side is in ATLAS notation. and are related to three-plane correlations through:
We extract the values of and from these equations using ATLAS data on event-plane correlations
Finally, involves the correlation between the magnitude of different harmonics, and . It is related to the so-called symmetric cumulant recently measured by the ALICE collaboration (27):
Figure 1–(a) displays the elements of as a function of the centrality percentile.
where is measured by the ALICE collaboration (4). The transformed elements are displayed in Fig. 1–(b). When written in terms of , the off-diagonal elements (full symbols) are much smaller than the diagonal elements (open symbols), which validates the approximations made in previous analyses (4) where they were neglected. However, they are not compatible with zero. In particular, is negative. It is therefore important to check to what extent the ALICE measurements, carried out under the assumption of negligible off-diagonal terms, are modified when the full pattern of correlations is taken into account.
As a by-product of our analysis, we can test whether the two components in the decomposition of , Eq. (1), are independent. We have imposed that they are uncorrelated, which is a weaker assumption. If they are independent, it implies in addition that the matrix element (open squares in Fig.1–(b)) factorizes into the product . This can be directly tested using ALICE data, as discussed in Appendix A.
Iv Nonlinear coefficients of from data
We now present our results for the nonlinear response coefficients of extracted from experimental data using Eq. (12). The coefficient has already been calculated in hydrodynamics (28); (10) but its experimental value is shown here for the first time. and have already been measured by ALICE under the approximation that the nonlinear terms are independent, while our new analysis takes into account the full correlation matrix.
Our results for the coefficients and are shown as full symbols in Fig. 2, and we compare them to the previous ALICE results, shown as open symbols. Comparison between the two sets of points shows that mutual correlations between nonlinear terms in Eq. (5) only have small effects. When they are taken into account, however, the centrality dependence of the response coefficient is somewhat flatter. As will be discussed in Sec. V, this generically improves agreement with hydrodynamics. Taking into account the error bars, one cannot exclude that both response coefficients are independent of centrality. Note that, although within error bars, extracted from the full correlation matrix appears to be systematically larger that the previous ALICE data. This is a signature of the largest non-diagonal term in , namely, .
Figure 3 displays the first experimental result for . It is larger than and for all centralities. It also has a stronger centrality dependence, but the large errors prevent any definite conclusion. We estimate these errors by taking into account only the error on the event-plane correlation (second line of Eq. (17)), which is the largest error. If a dedicated analysis of was carried out, however, the error bar would likely be as small as that on and . The reason is that a dedicated analysis would measure directly , which is linear in , while we extract it by combining the event-plane correlation and , which are both quadratic in , and therefore have a larger error. In Fig. 3 (open symbols) we present as well the coefficient extracted from data in absence of mutual correlations between the nonlinear terms. It is given by the following expression
Figure 4 displays the values of all the terms appearing in this equation as a function of centrality. is calculated as the difference between and the sum of the other terms. The largest nonlinear contribution to is that due to the coupling with and , while terms proportional to and are subleading for all centralities.
V Hydrodynamic calculation
Nonlinear response coefficients are unique probes of hydrodynamic behavior, because they are essentially independent of the initial state (9).
The uncertainty on the initial state is the bottleneck in traditional hydro-to-data comparisons (29), because traditional flow observables such as or are driven by initial anisotropies in the corresponding harmonics, which are poorly constrained.
The dependence on these initial anisotropies cancels in the ratios defining the nonlinear response coefficients.
This explains why all nonlinear response coefficients depend weakly on centrality in hydrodynamic calculations (13); (10).
This is in sharp contrast with the steep centrality dependence of .
Two different models of initial conditions with different and also return the same value for the nonlinear response coefficients.
Because of the aforementioned weak dependence on initial state, we evaluate these coefficients in a single collision event with a smooth initial density profile. The density profile is constructed by a smooth deformation of a symmetric 2-dimensional Gaussian, following Ref. (12). If one imprints a small elliptic deformation to the initial profile (a small asymmetry between the width along and ), the subsequent expansion generates all even Fourier harmonics, while odd harmonics (such as ) vanish because of symmetry. In this situation, is solely generated by , so that vanishes. Only the first nonlinear term remains in Eq. (5), therefore, . Similarly, one evaluates by imprinting a small triangular deformation to a radially symmetric profile, in which case . With this choice of initial conditions, the interpretation of and is transparent: they directly quantify the hexagonal flow produced by elliptic flow and triangular flow, respectively.
In the case of , one needs to deform the initial profile in harmonics 2 and 4. We carry out two hydrodynamic calculations, one with an asymmetric Gaussian profile, labeled (A) and one with a small quadrangular deformation (with the same orientation) on top of the asymmetric Gaussian, labeled (B). We compute , and at the end of the hydrodynamic evolution for both initial conditions. The initial conditions (A) and (B) differ only in the fourth Fourier harmonic, hence the values of are almost identical, . We first evaluate the change of induced by the quadrangular deformation, which is denoted by in Eq. (1), and given by . Then, the response coefficient defined by Eq. (5) is given by the increase of driven by the quadrangular deformation, i.e.:
We carry out this calculation for both ideal hydrodynamics and viscous hydrodynamics with (30).
We compute pion spectra at a freeze-out temperature of MeV.
This simple setup is justified by previous studies which have shown (13) that a more elaborate calculation taking into account the full hadron spectrum and strong decays returns almost identical nonlinear response coefficients.
The resulting nonlinear response coefficients are plotted as lines in Figs. 2 and 3.
Interestingly, the values for and are very close to those given by a full event-by-event hydrodynamic calculation (13), implementing MC Glauber or MC KLN initial conditions.
In our hydrodynamic calculation, the centrality enters only through the initial transverse radius, which we estimate in a Glauber model. As the centrality percentile increases, this radius decreases and off-equilibrium effects become larger (earlier freeze out, and larger dissipative corrections both during the hydrodynamic expansion and at freeze out). This explains the mild centrality dependence of nonlinear response coefficients. This mild dependence is a characteristic of hydrodynamic models in general. A steep variation of any nonlinear response coefficient would indicate a failure of the hydrodynamic picture. Experimental results for and are so far compatible with hydrodynamics. Our results on , on the other hand, seem not capture the centrality dependence of the extracted experimental values. This underlines the necessity of a dedicated analysis for reducing the error bars on this quantity.
Our results for and are mildly affected by adding shear viscosity to the hydrodynamic calculation. Viscosity results in a modest reduction of response coefficients, the effect being largest in peripheral collisions. For , we find a larger dependence on viscosity, and data are only compatible with viscous results. This large depenence is not observed in event-by-event calculations (13).
Vi Extension to higher harmonics
The data-driven analysis carried out in the previous sections show that the assumptions made in the literature about are reasonable: The matrix is essentially diagonal. In this section we argue, though, that it will be crucial to take into account the full pattern of correlations in harmonics of higher-order, such as (5), and potentially .
For heptagonal flow, there are also three leading nonlinear terms:
Note that involves four different Fourier harmonics, and is related to a 4-plane correlation (24). Hydrodynamic calculations of coefficients involving are already on the market (13); (15); (10), and always implicitly assume that the matrix is diagonal. It will be interesting to see if the hydrodynamic results for , and are modified once the full correlation structure is taken into account.
For completeness, let us also provide the correlation matrix of , a coefficient which is likely to be accessible to experimental analyses thank to the massive statistics of data collected in Pb+Pb collisions at LHC2. We decompose as follows:
which, in terms of the harmonics , reads
This leads to the following correlation matrix:
Let us stress once more, then, that the extraction of the nonlinear coefficients in our framework involves only moments which are linear in , and, therefore, experimentally easier to achieve than typical observables such as .
We have proposed a new framework which allows to systematically isolate the various nonlinear contributions to a given higher-order harmonic, and measure the nonlinear coupling coefficients. It can be applied to both experimental data and event-by-event hydrodynamic calculations. The main improvement over previous analyses is that we take into account the mutual correlations between nonlinear contributions. We have applied this new framework to using existing data. When mutual correlations are properly taken into account, the centrality dependence of the response coefficients and becomes somewhat flatter, thus improving agreement with hydrodynamic predictions. We have provided the first experimental determination of the coefficient coupling to and . It is in fair agreement with hydrodynamic predictions, though with large error bars. The corresponding nonlinear term is the largest of the three nonlinear contributions to . With the advent of large statistics Pb+Pb LHC2 data, we expect that this new framework will enable detailed analyses of higher-order flow coefficients, which will provide precision tests of hydrodynamic behavior.
Appendix A Correlation between and
In the definition of quadrangular flow in Eq. (1), one chooses and to be uncorrelated by construction, i.e.,
which implies statistical independence between the two terms. Statistical independence is a much stronger constraint than just requiring the two terms to be uncorrelated.
which vanishes if and are independent. All the terms appearing in the definition of are available from ALICE data. The quantity was recently measured (4), and is given by in Eq. (22). The resulting cumulant is displayed as black circles in Fig. 5.
This result illustrates that is indeed small in magnitude, so that Eq. (34) is a good approximation, although error bars, driven by the uncertainty on SC(4,2), are large in non-central collisions.
We ask now whether this quantity carries any useful information about the initial state of the hydrodynamic evolution. In heavy-ion collisions, elliptic flow is to a good approximation proportional to the second eccentricity harmonic of the initial state, (33). As for the fourth harmonic, it is often argued that may scale linearly with the fourth cumulant eccentricity of the initial medium (15), which is customarily taken from Ref. (12):
where corresponds to the moment-defined eccentricity harmonic of order (33).
Now, dubbing , if scales linearly with , then the following quantity
should match to a good extent observed in ALICE data, at least in central collisions where dominates. To check this, we compute from initial-state simulations of Pb+Pb collisions at TeV. We perform this by employing two different setups of the TENTo model of initial conditions (34). We use TENTo with , corresponding to a Glauber Monte Carlo model, and , which provides eccentricities in agreement with models including high-energy QCD effects, such as EKRT or IP-Glasma (34). Our results from the initial state models are shown in Fig. 5. The correlation is moderate and positive already at 0% centrality, therefore, we do not find agreement between the model calculation and experimental data. Therefore, typical models of initial conditions do not suggest a linear correlation between and , in contradiction, then, with the finding of hydrodynamic calculations (15). We remark, though, that this failure may be equally due to either the fact that the fourth-order initial-state anisotropy is not correctly quantified by , or that is not captured by the initial-state models, or both.
- We always neglect contributions proportional to , which is subleading.
- Specifically, in ALICE notation (4), , , , .
- In the centrality range where ATLAS uses 5% bins and ALICE uses 10% bins, we take the average event-plane correlation in two consecutive bins in ATLAS data.
- Since ALICE data on are available in the range 5-45%, and on in the range 0-35%, we are able to extract the matrix elements only in the range 5-35%.
- These analyses further make the approximation , whose validity is discussed in Appendix A.
- Qian et al (13) report a significant difference between MC Glauber and MC KLN initial conditions, but only for 2 out of 8 coefficients, namely, and . This difference is not seen by Zhao et al (10).
- For , Qian et. al (13) find a different result depending on initial conditions. Their MC KLN result is similar to ours, while the MC Glauber result is much larger, and incompatible with data.
- U. Heinz and R. Snellings, Ann. Rev. Nucl. Part. Sci. 63, 123 (2013) doi:10.1146/annurev-nucl-102212-170540 [arXiv:1301.2826 [nucl-th]].
- G. Aad et al. [ATLAS Collaboration], Phys. Rev. C 86, 014907 (2012) doi:10.1103/PhysRevC.86.014907 [arXiv:1203.3087 [hep-ex]].
- S. Chatrchyan et al. [CMS Collaboration], Phys. Rev. C 89, no. 4, 044906 (2014) doi:10.1103/PhysRevC.89.044906 [arXiv:1310.8651 [nucl-ex]].
- S. Acharya et al. [ALICE Collaboration], Phys. Lett. B 773, 68 (2017) doi:10.1016/j.physletb.2017.07.060 [arXiv:1705.04377 [nucl-ex]].
- S. Tuo [CMS Collaboration], Nucl. Phys. A 967, 381 (2017). doi:10.1016/j.nuclphysa.2017.05.064
- M. Luzum and J. Y. Ollitrault, Nucl. Phys. A 904-905, 377c (2013) doi:10.1016/j.nuclphysa.2013.02.028 [arXiv:1210.6010 [nucl-th]].
- S. Chatrchyan et al. [CMS Collaboration], JHEP 1402, 088 (2014) doi:10.1007/JHEP02(2014)088 [arXiv:1312.1845 [nucl-ex]].
- F. G. Gardim, F. Grassi, M. Luzum and J. Y. Ollitrault, Phys. Rev. C 85, 024908 (2012) doi:10.1103/PhysRevC.85.024908 [arXiv:1111.6538 [nucl-th]].
- L. Yan and J. Y. Ollitrault, Phys. Lett. B 744, 82 (2015) doi:10.1016/j.physletb.2015.03.040 [arXiv:1502.02502 [nucl-th]].
- W. Zhao, H. j. Xu and H. Song, Eur. Phys. J. C 77, no. 9, 645 (2017) doi:10.1140/epjc/s10052-017-5186-x [arXiv:1703.10792 [nucl-th]].
- N. Borghini and J. Y. Ollitrault, Phys. Lett. B 642, 227 (2006) doi:10.1016/j.physletb.2006.09.062 [nucl-th/0506045].
- D. Teaney and L. Yan, Phys. Rev. C 86, 044908 (2012) doi:10.1103/PhysRevC.86.044908 [arXiv:1206.1905 [nucl-th]].
- J. Qian, U. W. Heinz and J. Liu, Phys. Rev. C 93, no. 6, 064901 (2016) doi:10.1103/PhysRevC.93.064901 [arXiv:1602.02813 [nucl-th]].
- L. V. Bravina et al., Phys. Rev. C 89, no. 2, 024909 (2014) doi:10.1103/PhysRevC.89.024909 [arXiv:1311.0747 [hep-ph]].
- J. Qian, U. Heinz, R. He and L. Huo, Phys. Rev. C 95, no. 5, 054908 (2017) doi:10.1103/PhysRevC.95.054908 [arXiv:1703.04077 [nucl-th]].
- P. Liu and R. A. Lacey, arXiv:1802.06595 [nucl-ex].
- G. Aad et al. [ATLAS Collaboration], Phys. Rev. C 90, no. 2, 024905 (2014) doi:10.1103/PhysRevC.90.024905 [arXiv:1403.0489 [hep-ex]].
- R. S. Bhalerao, J. Y. Ollitrault and S. Pal, Phys. Lett. B 742, 94 (2015) doi:10.1016/j.physletb.2015.01.019 [arXiv:1411.5160 [nucl-th]].
- F. Cooper and G. Frye, Phys. Rev. D 10, 186 (1974). doi:10.1103/PhysRevD.10.186
- D. Teaney, Phys. Rev. C 68, 034913 (2003) doi:10.1103/PhysRevC.68.034913 [nucl-th/0301099].
- R. S. Bhalerao, M. Luzum and J. Y. Ollitrault, Phys. Rev. C 84, 034910 (2011) doi:10.1103/PhysRevC.84.034910 [arXiv:1104.4740 [nucl-th]].
- K. Fukushima, D. E. Kharzeev and H. J. Warringa, Phys. Rev. D 78, 074033 (2008) doi:10.1103/PhysRevD.78.074033 [arXiv:0808.3382 [hep-ph]].
- S. Voloshin and Y. Zhang, Z. Phys. C 70, 665 (1996) doi:10.1007/s002880050141 [hep-ph/9407282].
- R. S. Bhalerao, J. Y. Ollitrault and S. Pal, Phys. Rev. C 88, 024909 (2013) doi:10.1103/PhysRevC.88.024909 [arXiv:1307.0980 [nucl-th]].
- K. Aamodt et al. [ALICE Collaboration], Phys. Rev. Lett. 107, 032301 (2011) doi:10.1103/PhysRevLett.107.032301 [arXiv:1105.3865 [nucl-ex]].
- S. Acharya et al. [ALICE Collaboration], Phys. Rev. C 97, no. 2, 024906 (2018) doi:10.1103/PhysRevC.97.024906 [arXiv:1709.01127 [nucl-ex]].
- J. Adam et al. [ALICE Collaboration], Phys. Rev. Lett. 117, 182301 (2016) doi:10.1103/PhysRevLett.117.182301 [arXiv:1604.07663 [nucl-ex]].
- J. Qian and U. Heinz, Phys. Rev. C 94, no. 2, 024910 (2016) doi:10.1103/PhysRevC.94.024910 [arXiv:1607.01732 [nucl-th]].
- M. Luzum and P. Romatschke, Phys. Rev. C 78, 034915 (2008) Erratum: [Phys. Rev. C 79, 039903 (2009)] doi:10.1103/PhysRevC.78.034915, 10.1103/PhysRevC.79.039903 [arXiv:0804.4015 [nucl-th]].
- P. Kovtun, D. T. Son and A. O. Starinets, Phys. Rev. Lett. 94, 111601 (2005) doi:10.1103/PhysRevLett.94.111601 [hep-th/0405231].
- G. Giacalone, L. Yan, J. Noronha-Hostler and J. Y. Ollitrault, Phys. Rev. C 94, no. 1, 014906 (2016) doi:10.1103/PhysRevC.94.014906 [arXiv:1605.08303 [nucl-th]].
- A. Bilandzic, C. H. Christensen, K. Gulbrandsen, A. Hansen and Y. Zhou, Phys. Rev. C 89, no. 6, 064904 (2014) doi:10.1103/PhysRevC.89.064904 [arXiv:1312.3572 [nucl-ex]].
- D. Teaney and L. Yan, Phys. Rev. C 83, 064904 (2011) doi:10.1103/PhysRevC.83.064904 [arXiv:1010.1876 [nucl-th]].
- J. S. Moreland, J. E. Bernhard and S. A. Bass, Phys. Rev. C 92, no. 1, 011901 (2015) doi:10.1103/PhysRevC.92.011901 [arXiv:1412.4708 [nucl-th]].