Including off-diagonal anisotropies in anisotropic hydrodynamics

TL;DR

This paper introduces an efficient method to incorporate off-diagonal local rest frame momentum anisotropies into anisotropic hydrodynamics by diagonalizing the anisotropy tensor, simplifying calculations of relevant moments.

Contribution

The authors develop a novel diagonalization-based approach to include off-diagonal anisotropies in leading-order anisotropic hydrodynamics, reducing computational complexity.

Findings

Method efficiently includes off-diagonal anisotropies.

Framework applicable to energy-momentum tensor calculations.

Simplifies viscous update equations in anisotropic hydrodynamics.

Abstract

In this paper we present a method for efficiently including the effects of off-diagonal local rest frame momentum anisotropies in leading-order anisotropic hydrodynamics. The method relies on diagonalization of the space-like block of the anisotropy tensor and allows one to reduce the necessary moments of the distribution function in the off-diagonal case to a linear combination of diagonal-anisotropy integrals. Once reduced to diagonal-anisotropy integrals, the results can be computed efficiently using techniques described previously in the literature. We present a general framework for how to accomplish this and provide examples for off-diagonal anisotropy moments entering into the energy-momentum tensor and viscous update equations which emerge when performing anisotropic pressure matching.

Equations

$$\displaystyle f_{a}(x,p)=f_{\rm iso}\bigg{(}\frac{1}{\lambda}\sqrt{p_{\mu}\Xi^{\mu\nu}p_{\nu}}\bigg{)}\,,$$

$${\Xi }^{\mu \nu }=u^{\mu }{u}^{\nu }+{\xi }^{\mu \nu }-\Phi {\Delta }^{\mu \nu },$$

$$u^{\mu }{u}_{\mu }$$

$${{\xi }^{\mu }}_{\mu }$$

$$u_{\mu }{\xi }^{\mu \nu }$$

$${\xi }^{00}={\xi }^{0i}={\xi }^{i0}=0.$$

$$p\cdot \Xi \cdot p=\mathbf{p}\cdot \kappa \cdot \mathbf{p}+m^2,$$

$$f_a(x,p)=f_{\mathrm{iso}}\left(\frac{1}{\lambda }\sqrt{\mathbf{p}\cdot \kappa \cdot \mathbf{p}+m^2}\right),$$

$$\kappa \equiv \mathbf{I}(1+\Phi )+\mathbf{\xi },$$

$$\mathbf{\xi }=\mathrm{diag}({\xi }_x,{\xi }_y,{\xi }_z),$$

$$\kappa =\left(\begin{array}{ccc}1/{\alpha }_x^2& {\xi }_{xy}& {\xi }_{xz}\\ {\xi }_{xy}& 1/{\alpha }_y^2& {\xi }_{yz}\\ {\xi }_{xz}& {\xi }_{yz}& 1/{\alpha }_z^2\end{array}\right).$$

$$\kappa =A{\kappa }_D{A}^{†},$$

$$\displaystyle{\bf p}\cdot\kappa\cdot{\bf p}={\bf p}^{T}\kappa\,{\bf p}=\Big{[}{\bf p}^{T}A\Big{]}\Big{[}A^{\dagger}\kappa A\Big{]}\Big{[}A^{\dagger}{\bf p}\Big{]}=\tilde{\bf p}^{T}\,\kappa_{D}\,\tilde{{\bf p}}=\tilde{\bf p}\cdot\kappa_{D}\cdot\tilde{\bf p}\,,$$

$$\mathbf{p}=A\tilde{\mathbf{p}}\Rightarrow p_i=\sum _j{A}_{ij}{\tilde{p}}_j.$$

$$p_i=\sum _j{A}_{ij}{\tilde{p}}_j=\sum _j{v}_i^{(j)}{\tilde{p}}_j,$$

$$E=\sqrt{{\mathbf{p}}^2+m^2}=\sqrt{{\tilde{\mathbf{p}}}^2+m^2}=\tilde{E}.$$

$$f_a(x,p)=f_{\mathrm{iso}}\left(\frac{1}{\lambda }\sqrt{\mathbf{p}\cdot \kappa \cdot \mathbf{p}+m^2}\right)=f_{\mathrm{iso}}\left(\frac{1}{\lambda }\sqrt{\tilde{\mathbf{p}}\cdot {\kappa }_D\cdot \tilde{\mathbf{p}}+m^2}\right)\equiv f_a^D(x,\tilde{p}).$$

$$T^{00}$$

$$T^{ij}$$

$$\mathcal{E}$$

$$T^{ij}$$

$${\mathcal{P}}_i=T^{ii}=\tilde{\alpha }{\lambda }^4\sum _{k=1}^3{\left[v_i^{(k)}\right]}^2{\tilde{\alpha }}_k^2{Q}_3^k({\tilde{\alpha }}_x^2,{\tilde{\alpha }}_y^2,{\tilde{\alpha }}_z^2,\hat{m}).$$

$${\partial }_{\mu }{T}^{\mu \nu }=0,$$

$$D_u{\pi }^{⟨\mu \nu ⟩}+\frac{1}{{\tau }_{\mathrm{eq}}}{\pi }^{\mu \nu }$$

$$D_u\mathcal{P}+\frac{1}{{\tau }_{\mathrm{eq}}}\left(\mathcal{P}-{\mathcal{P}}_{\mathrm{eq}}\right)$$

$$\begin{array}{rl}{D}_u& =u^{\mu }{\partial }_{\mu },\\ D_i& =X_i^{\mu }{\partial }_{\mu },\\ \theta & ={\nabla }_{\mu }{u}^{\mu },\end{array}\begin{array}{rl}{\nabla }^{\mu }& ={\Delta }^{\mu \nu }{\partial }_{\nu },\\ {\omega }^{\mu \nu }& =({\nabla }^{\mu }{u}^{\nu }-{\nabla }^{\nu }{u}^{\mu })/2,\\ {\sigma }^{\mu \nu }& ={\Delta }_{\alpha \beta }^{\mu \nu }{\partial }^{\alpha }{u}^{\beta }.\end{array}$$

$$T^{\mu \nu }=T_{\mathrm{eq}}^{\mu \nu }+{\pi }^{\mu \nu }+\Pi {\Delta }^{\mu \nu },$$

$$\left({\sigma }_{\rho \sigma }+\frac{1}{3}\theta {\Delta }_{\rho \sigma }\right)p^{\rho }{p}^{\sigma }=p^{\rho }{p}^{\sigma }{\nabla }_{\rho }{u}_{\sigma }=p^{\sigma }(\mathbf{p}\cdot \mathbf{D})u_{\sigma },$$

$$p^{⟨\mu }{p}^{\nu ⟩}=p^{\alpha }{p}^{\beta }{\Delta }_{\alpha }^{\mu }{\Delta }_{\beta }^{\nu }+\frac{1}{3}{\Delta }^{\mu \nu }{\mathbf{p}}^2=p^i{p}^j{X}_i^{\mu }{X}_j^{\nu }+\frac{1}{3}{\Delta }^{\mu \nu }{\mathbf{p}}^2,$$

$$p\cdot \Delta \cdot p=-{\mathbf{p}}^2.$$

$$\displaystyle-\bigg{(}\sigma_{\rho\sigma}+\frac{1}{3}\theta\Delta_{\rho\sigma}\bigg{)}\int\frac{dP}{E^{2}}p^{\rho}p^{\sigma}p^{\langle\mu}p^{\nu\rangle}f_{a}$$

$$=-\int \frac{dP}{E^2}{f}_a{p}^i{p}^j{X}_i^{\mu }{X}_j^{\nu }{p}^{\sigma }(\mathbf{p}\cdot \mathbf{D})u_{\sigma }-\frac{{\Delta }^{\mu \nu }}{3}\int \frac{dP}{E^2}{f}_a{\mathbf{p}}^2{p}^{\sigma }(\mathbf{p}\cdot \mathbf{D})u_{\sigma }$$

$$=\int \frac{dP}{E^2}{f}_a{p}^i{p}^j{X}_i^{\mu }{X}_j^{\nu }{p}^l(\mathbf{p}\cdot \mathbf{D})u_l+\frac{{\Delta }^{\mu \nu }}{3}\int \frac{dP}{E^2}{f}_a{\mathbf{p}}^2{p}^l(\mathbf{p}\cdot \mathbf{D})u_l$$

$$=\left[{\mathcal{F}}^{ijkl}{X}_i^{\mu }{X}_j^{\nu }+\frac{{\Delta }^{\mu \nu }}{3}{\mathcal{F}}^{iikl}\right]D_k{u}_l,$$

Full text

Quark-gluon plasma, Relativistic heavy-ion collisions, Anisotropic hydrodynamics, Equation of state, Boltzmann equation, Off-diagonal anisotropy

I Introduction

Ultra-relativistic heavy ion collision (URHIC) experiments, e.g. RHIC at BNL and LHC at CERN, aim to study the dynamics and properties of matter at extremely high-energy density. In these experiments, matter is heated to temperatures exceeding the QCD pseudo-critical temperature, $T_{\rm pc}\simeq 155$ MeV, using ultra-relativistic collisions among heavy nuclei, protons, deuterons, etc. The strongly interacting droplet of matter produced during high-energy and high-multiplicity URHICs is called the quark-gluon plasma (QGP). In high-multiplicity events, the QGP demonstrates strong collective behavior during evolution from hydrodynamization ($\tau\sim 1$ fm/c) to hadronic freeze out ($\tau\sim 10$ fm/c). During this time period it has been found that relativistic fluid dynamics formalisms can effectively describe the evolution of the system and one finds that information about initial state geometry of the target (average eccentricity and fluctuations) is reflected in final state observables, e.g. the azimuthal dependence of hadron production. In other words, one can track the correlations between the eccentricity of the initial state’s geometry and the flow harmonics observed in the final state hadron spectra using dissipative hydrodynamics. The success of relativistic dissipative hydrodynamics Heinz and Snellings (2013); Jeon and Heinz (2016); Romatschke and Romatschke (2017); Alqahtani et al. (2018) has inspired theoreticians to make the underlying formalisms more complete and robust with respect to large deviations from isotropic thermal equilibrium using standard fixed-order viscous hydrodynamics (vHydro) treatments Muller (1967); Israel (1976); Israel and Stewart (1979); Muronga (2002, 2004); Muronga and Rischke (2004); Heinz et al. (2006); Baier et al. (2006); Romatschke and Romatschke (2007); Baier et al. (2008); Dusling and Teaney (2008); Luzum and Romatschke (2008); Song and Heinz (2009); Heinz (2010); Schenke et al. (2011a, b); Bozek (2011); Niemi et al. (2011); Denicol et al. (2011); Niemi et al. (2012); Bożek and Wyskiel-Piekarska (2012); Denicol et al. (2012a, b); Jaiswal (2013a, b); Denicol et al. (2014a, b); Jaiswal et al. (2014) and resummed anisotropic hydrodynamics (aHydro) treatments Florkowski and Ryblewski (2011); Martinez and Strickland (2010); Ryblewski and Florkowski (2011); Florkowski and Ryblewski (2012); Martinez et al. (2012); Ryblewski and Florkowski (2012); Bazow et al. (2014); Tinti and Florkowski (2014); Nopoush et al. (2014); Florkowski et al. (2014a); Tinti (2016); Bazow et al. (2015, 2016); Nopoush et al. (2015a); Alqahtani et al. (2015); Molnar et al. (2016a, b); Bluhm and Schaefer (2015, 2016); Alqahtani et al. (2017a, b, 2018); Almaalol et al. (2018).

The introduction of the aHydro formalism was driven by the fact that, due to the strong early-stage longitudinal expansion of the QGP, one finds large momentum-space anisotropy in the local rest frame (LRF) of the QGP which persists for many fm/c. The magnitude of the momentum-space anisotropy has cast some doubt on the quantitative accuracy of standard vHydro which assumes that one can linearize around isotropic equilibrium. aHydro is a non-equilibrium hydrodynamics model which takes into account the strong momentum-space anisotropy of the QGP at leading order and in doing so resums an infinite number of terms in inverse Reynolds number Strickland et al. (2018). In contrast to standard vHydro, aHydro is based on Taylor expansion about an anisotropic distribution function instead of an isotropic one. This allows one to capture the dominant anisotropic contributions to the distribution function in the leading order term, thereby guaranteeing positivity of the one-particle distribution at all space-time points at leading-order. aHydro and vHydro have been tested against exact solutions of the Boltzmann equation for systems subject to Bjorken Florkowski et al. (2013a, b, 2014b); Strickland et al. (2018); Florkowski et al. (2018); Strickland (2018) and Gubser flows Denicol et al. (2014c); Nopoush et al. (2015b); Denicol et al. (2014c); Martinez et al. (2017); Behtash et al. (2018). In all cases, it was found that aHydro provided the best approximation to the exact solutions for both hydrodynamic and non-hydrodynamic moments of the distribution function Strickland (2018).

This provided motivation to compare the aHydro framework with experimental results. Despite the success of these early comparisons, in all phenomenological applications of aHydro to date, leading-order aHydro codes have been implemented using an anisotropy tensor which possesses only diagonal (elliptical) anisotropies (see Ref. Alqahtani et al. (2018) for a recent review). This was done mainly because of the difficulty of efficiently evaluating the necessary moment integrals in the presence of off-diagonal anisotropies $\xi^{ij}$ with $i\neq j$. However, to be complete, one must also include the possibility of off-diagonal leading-order anisotropies. Near equilibrium, this is equivalent to including off-diagonal components in the LRF shear viscous tensor $\pi^{ij}$.

In this paper, we present a technique that can be used to efficiently include non-vanishing $\xi^{ij}$. This is done by a change of variables in the generic moment integrals which diagonalizes the anisotropy tensor. Once cast into diagonal form, a previously developed technique for the efficient application of diagonal moment integrals can be used to compute the necessary off-diagonal moment integrals (see Appendix B of Ref. Alqahtani et al. (2017b)). We present the general method of diagonalization and provide some concrete examples for the application to aHydro frameworks which use the so-called anisotropic-pressure- or Tinti-matching Tinti (2016); Molnar et al. (2016a).

Conventions and notation

The Minkowski metric tensor is taken to be “mostly minus”, i.e. $g^{\mu\nu}={\rm diag}(+,-,-,-)$. The vector $u^{\mu}$ is the flow velocity which satisfies normalization condition $u_{\mu}u^{\mu}=1$. The transverse projection operator $\Delta^{\mu\nu}\equiv g^{\mu\nu}{-}u^{\mu}u^{\nu}$ is used to project four-vectors and/or tensors into the space orthogonal to $u^{\mu}$. Parentheses and square brackets on indices denote symmetrization and anti-symmetrization, respectively, i.e. $A^{(\mu\nu)}\equiv\frac{1}{2}\left(A^{\mu\nu}{+}A^{\nu\mu}\right)$ and $A^{[\mu\nu]}\equiv\frac{1}{2}\left(A^{\mu\nu}{-}A^{\nu\mu}\right)$. Angle brackets on indices indicate projection with a four-index transverse projector, $A^{\langle\mu\nu\rangle}\equiv\Delta^{\mu\nu}_{\alpha\beta}A^{\alpha\beta}$, where $\Delta^{\mu\nu}_{\alpha\beta}\equiv\Delta^{(\mu}_{\alpha}\Delta^{\nu)}_{\beta}-\Delta^{\mu\nu}\Delta_{\alpha\beta}/3$ projects out the traceless and $u^{\mu}$-transverse components of a rank-two tensor. The Lorentz-invariant momentum-space integration measure is indicated as $dP=\tilde{N}d^{3}{\bf p}/(p\cdot u)$, with $\tilde{N}=N_{\rm dof}/(2\pi)^{3}$ where $N_{\rm dof}$ is the number of degrees of freedom.

In order to write the equations of motion in a manifestly Lorentz-covariant manner it is useful to introduce the LRF basis vectors, as $u^{\mu}_{\rm LRF}=(1,{\bf 0})$ and $X^{\mu}_{i,{\rm LRF}}=(0,{\boldsymbol{\delta}}^{\mu}_{i})$ with $i\in\{1,2,3\}$. By applying a sequence of Lorentz transformations, one can construct the lab frame basis vectors, i.e. $u^{\mu}$ and $X^{\mu}_{i}$ with $i\in\{1,2,3\}$, where the dynamical equations are solved and particle spectra are computed Ryblewski and Florkowski (2011); Martinez et al. (2012). It is also useful to define the transverse projection operator in terms of the space-like basis vectors, i.e. $\Delta^{\mu\nu}=-\sum_{i}X^{\mu}_{i}X^{\nu}_{i}$. Finally, note that the Latin indices sum over space-like indices (components of three-vectors) and Greek indices sum over components of four-vectors.

II Leading-order anisotropic hydrodynamics

In leading-order aHydro, the one-particle distribution function is parametrized by an anisotropy tensor which results in the deformation of the argument of an isotropic distribution function into an anisotropic one Martinez et al. (2012); Nopoush et al. (2014)

[TABLE]

where $\lambda$ has dimensions of energy and can be identified with temperature only in the isotropic equilibrium limit. In practice, $f_{\rm iso}$ can be a Bose-Einstein, Fermi-Dirac, or Maxwell-Boltzmann distribution depending on particle statistics and/or energy. In the non-conformal (massive) case, the rank-2 tensor $\Xi^{\mu\nu}$ specifying the shape of the distribution in momentum space is defined as Martinez et al. (2012); Nopoush et al. (2014)

[TABLE]

where $\xi^{\mu\nu}$ denotes a symmetric traceless anisotropy tensor, i.e. $\xi_{x}+\xi_{y}+\xi_{z}=0$ in the LRF.

The quantities $\lambda$, $u^{\mu}$, and $\xi^{\mu\nu}$ are spacetime fields which satisfy the following identities

[TABLE]

The third condition above, indicating orthogonality of $\xi^{\mu\nu}$ to $u^{\mu}$ which implies that, in the LRF, $\xi^{\mu\nu}$ obeys the following conditions

[TABLE]

Working in the LRF, this allows us to focus on the non-trivial space-like components of $\xi^{\mu\nu}$ as $\bm{\xi}$, which is a $3\times 3$ matrix. The argument of distribution function subject to the mass-shell condition can be simplified as

[TABLE]

which gives

[TABLE]

where

[TABLE]

with $\bm{I}$ being a $3\times 3$ identity matrix.

If $\bm{\xi}$ is diagonal, i.e.

[TABLE]

which implies the ellipsoidal distribution, the $\kappa$ matrix is automatically diagonal, i.e. $\kappa={\rm diag}(1/\alpha_{x}^{2},1/\alpha_{y}^{2},1/\alpha_{z}^{2})$ with $\alpha_{i}=(1+\xi_{i}+\Phi)^{-1/2}$ Nopoush et al. (2014). For a non-ellipsoidal distribution function, generalizing $\bm{\xi}$ to include off-diagonal components, one has

[TABLE]

Note that, in a general frame one has $\xi^{\mu\nu}=\kappa_{ij}X^{\mu}_{i}X^{\nu}_{j}$ where the summation over $i$ and $j$ is implied.

III Diagonalization

Calculating the bulk variables in aHydro requires computing momentum-space moments of the distribution function. However, the distribution function in Eq. (8) is a complicated function of momentum and there is no way to perform the integrals analytically except in some special cases. In this section, we introduce an algebraic method to diagonalize the $\kappa$ matrix so that we can reduce the computation of moment-integrals including off-diagonal anisotropies to a linear combination of diagonal momentum-space moment integrals.

For any $N\times N$ real and symmetric matrix $\kappa$ there exists a unitary matrix $A$ such that

[TABLE]

where $A$ is constructed such that its columns are the eigenvectors of $\kappa$. The combination ${\bf p}\cdot\kappa\cdot{\bf p}$ can be written as

[TABLE]

with $\tilde{\bf p}\equiv A^{\dagger}{\bf p}$. By definition we have

[TABLE]

For example

[TABLE]

where the vector $v^{(i)}=(v_{x}^{(i)},v_{y}^{(i)},v_{z}^{(i)})$ is the $i^{\rm th}$ eigenvector of $\kappa$. Therefore, we have two frames, i.e. the original frame and the rotated frame, where the components of the momentum vector are ${\bf p}_{i}$ and $\tilde{\bf p}_{i}$, respectively. The $\kappa$ matrix in the original frame is defined in Eq. (11) and, in the rotated frame, is defined as $\kappa_{D}\equiv{\rm diag}(1/\tilde{\alpha}_{i}^{2})$. These two frames are connected by rotations through a set of Euler angles. Note that the Jacobian for transforming between two frames is unity.

It is obvious that the length of $\bf p$ is invariant under this coordinate transformation. Accordingly, as expected, $E$ is the same in both coordinate systems

[TABLE]

Using Eq. (15), one can simplify the general anisotropic distribution function to the anisotropic distribution function (8) with diagonal anisotropy tensor (11) in the rotated frame

[TABLE]

IV The energy-momentum tensor

We begin by demonstrating how this method can be used to efficiently evaluate the components of the energy-momentum tensor including off-diagonal anisotropies. In the general case, we have six independent anisotropy parameters ($\alpha_{x}$, $\alpha_{y}$, $\alpha_{z}$, $\xi_{xy}$, $\xi_{xz}$, and $\xi_{yz}$), one momentum-scale parameter ($\lambda$), and the three independent components of the fluid four-velocity ($u^{i}$), resulting in ten space-time fields for which we must obtain equations of motion. In the LRF, the non-vanishing components of the energy-momentum tensor are

[TABLE]

Using the techniques introduced in the previous section, one finds

[TABLE]

and

[TABLE]

The $Q$-functions appearing above only depend on the diagonal anisotropies $\tilde{\boldsymbol{\alpha}}$ and are defined in Appendix A. The scaled mass variable is defined as $\hat{m}\equiv m/\lambda$ and we have introduced a compact notation as $\tilde{\alpha}\equiv\tilde{\alpha}_{x}\tilde{\alpha}_{y}\tilde{\alpha}_{z}$. Based on the symmetry of $T^{ij}$ under exchanging the indices, out of 9 possible values there are only 6 unique terms that must be calculated. Note that for the diagonal terms (pressures) one obtains

[TABLE]

In all cases above, we have reduced the problem to computing $Q$-functions with only diagonal anisotropies. The diagonal anisotropy tensor integrals can be well-approximated by Taylor expanding to high-order around an isotropic point, e.g. $\tilde{\boldsymbol{\alpha}}_{\rm iso}=(\alpha_{0},\alpha_{0},\alpha_{0})$. At each order in this expansion the required integrals can be performed analytically. In order to cover the space using truncated Taylor expansions, one can utilize multiple expansion points which are then pieced together to accurately span the range of diagonal anisotropies which are generated in typical simulations. Using modern computerized algebra systems one can extend the Taylor expansion expressions described above to high order. In practice, phenomenological codes have used $12^{\rm th}$ order truncations in $\tilde{\boldsymbol{\delta}}=\tilde{\boldsymbol{\alpha}}-\tilde{\boldsymbol{\alpha}}_{\rm iso}$ (see Appendix B of Ref. Alqahtani et al. (2017b)).

V Dynamical equations - Anisotropic pressure matching

To further demonstrate the utility of this method, we now consider equations for the viscous tensor obtained by anisotropic pressure matching Tinti (2016). In relaxation-time approximation (RTA) the dynamical equations for the shear and bulk viscous corrections based on anisotropic matching are

[TABLE]

In the above relations, $f_{a}$ is the general distribution function defined at Eq. (1). The tensor $\pi^{\mu\nu}$ is the shear tensor, which is traceless and orthogonal to flow velocity $u^{\mu}$. In the relations above one has ${\cal P}={\cal P}_{\rm eq}+\Pi$ with ${\cal P}_{\rm eq}$ being the LRF equilibrium pressure which can be obtained by evaluating any component of Eq. (22) with $\kappa$ equal to an identity matrix and $\lambda$ set to the local effective temperature $T$. The other symbols appearing in Eqs. (24) and (25) above are defined as

[TABLE]

The equations (25) represent a set of ten dynamical equations for the ten independent macroscopic variables of the system. Microscopically, one has three components of flow velocity $u^{i}$, six independent anisotropy parameters, and the temperature-like scale $\lambda$, resulting in ten dynamical microscopic variables. Correspondingly, when coding up these equations, one can choose between using macroscopic or microscopic variables. In addition, if using the macroscopic variables, one can evolve the ten independent components of (symmetric) energy-momentum tensor $T^{\mu\nu}$ or one can use the standard decomposition McNelis et al. (2018); Tinti (2016)

[TABLE]

which has as dynamical variables ${\cal E}$, three components of flow velocity $u^{i}$, five independent components of shear tensor $\pi^{\mu\nu}$, and the bulk viscous correction $\Pi$, again added up to ten.

In practice, it is preferable to evolve the macroscopic (thermodynamics) variables, since modern flux-conserving algorithms are better suited to these equations than the microscopic ones. However, this procedure is non-trivial because, although the above equations evolve macroscopic variables, they explicitly contain microscopic ones as well, e.g. the distribution function $f_{a}$ appearing in Eq. (8). Therefore, in order to close the system of equations one must update the microscopic variables in parallel to the macroscopic ones during the evolution. Roughly speaking, the procedure is as follows: The equation $\partial_{\mu}T^{\mu\nu}=0$ provides the evolution of ${\cal E}$ and $u^{i}$. The other equations evolve the components of the shear tensor. Using these, one can construct the full $T^{\mu\nu}$ using (27). Once the lab frame $T^{\mu\nu}$ is evolved forward one time step, the updated microscopic variables can be obtained by boosting to the LRF and solving a set of seven coupled matching equations which match $T^{00}_{\rm LRF}$, and six components of upper diagonal space-like block of $T_{\rm LRF}^{\mu\nu}$ to their microscopic definitions as a function of $\alpha_{i}$, $\xi_{ij}$, and $\lambda$, i.e. Eqs. (20) and (21).

In order to further develop the necessary formalism, one must expand and simplify the dynamical equations (25) for the case of a non-ellipsoidal anisotropic distribution function. Note that we will expand the equations in the lab frame, where the dynamical equations are solved. However, whenever a scalar quantity is obtained, we have the freedom to choose a covariant Lorentz frame, e.g. local reference frame, where the calculation is simpler.

There are two terms in Eqs. (24) and (25) needing detailed expansion. The first one is

[TABLE]

where ${\bf D}$ is defined in (26). We also have

[TABLE]

where the Einstein summation convention for repeated spatial indices is applied. The very last step is performed in order to make the dependence of components of momentum 3-vector explicit, which is useful in evaluating the integrals necessary.

The other term is

[TABLE]

Using the above relations, one can expand the following integrals

[TABLE]

where the four-index function introduced above is defined as

[TABLE]

Note that for Eq. (31) to be non-vanishing one must have an even number of spatial momenta with matching indices, appearing in ${\cal F}^{ijkl}$. To see this, consider the integral above containing an odd number of spatial momenta. Using the map (15) it will contain an odd number of $\tilde{p_{i}}$ even in the rotated frame and the rest of the integrand will be an even function of the momenta. Therefore, the integral will vanish by symmetry in this case. This suggests that in the third line of the equation (31) defined above, $p^{\sigma}u_{\sigma}\rightarrow-p^{i}u_{i}$.

Similarly, the non-trivial term appearing in the bulk viscous equation of motion (25) is

[TABLE]

To complete the simplification of the non-trivial terms in (31) and (33), we now consider the ${\cal F}$ function. Using similar techniques as used for the $T^{ij}$, one obtains

[TABLE]

The operator ${\cal P}_{mn}$ introduced above is the permutation operator which sums over all possible permutations of $m$ and $n$ in the operand (including repeated ones). Based on the symmetry of ${\cal F}^{ijkl}$ under exchanging the indices, out of 81 possible values there are only 15 unique terms that must be calculated. The function $Q^{mn}$ introduced above is defined in Appendix A.

VI Discussion and summary

As we demonstrated in the previous two sections, one can reduce the problem of evaluating complicated off-diagonal anisotropy moment integrals to a sum of diagonal anisotropy integrals. In practice, one can use Eqs. (23), (24), and (25) to evolve the energy-momentum tensor, shear viscous tensor, and the bulk viscous correction, respectively. Given an initial condition specified in terms of all anisotropies and the momentum scale, $\lambda$, one can construct the full energy-momentum tensor at the initial time. One can then evolve the coupled partial differential equations (23), (24), and (25) forward in time by one infinitesimal step making use of the methods explained in the previous section to evaluate the non-trivial integrals involving $f_{a}$ in Eqs. (24) and (25).

Once the update is complete, one can solve a set of seven non-linear equations to extract the updated LRF anisotropies and scale parameter. These can then be used to compute the non-trivial integrals involving $f_{a}$ in the next time step. Repeating this procedure, one can evolve all dynamical fields using Eqs. (23), (24), and (25). Critical to accomplishing this is the efficient evaluation of the integrals involving $f_{a}$ in Eqs. (24) and (25) and the subsequent extraction of the local anisotropy tensor from the full energy-momentum tensor. The diagonalization method described in the previous two sections solves this problem by removing the bottleneck of evaluating complicated three dimensional integrals on demand.

VII Conclusions

In this paper we presented a method for efficiently including the effects of off-diagonal local rest frame momentum anisotropies in leading-order anisotropic hydrodynamics. The method relies on diagonalization of the space-like block of the anisotropy tensor and allows one to reduce the necessary moments of the distribution function in the off-diagonal case to a linear combination of diagonal-anisotropy integrals. Once reduced to diagonal-anisotropy integrals, the results can be computed efficiently using techniques described previously in the literature Alqahtani et al. (2017b). We presented a general framework for how to accomplish this and provided examples for off-diagonal anisotropy moments entering into the energy-momentum tensor and viscous update equations which emerge when performing anisotropic pressure matching Tinti (2016). With this method in hand one can implement a leading-order anisotropic hydrodynamics code that takes into account off-diagonal anisotropies non-perturbatively. Additionally, since the equations are formulated at the level of the energy-momentum tensor and shear viscous tensor, this more easily allows for the use of advanced numerical techniques for solving the necessary partial differential equations (see e.g. Bazow et al. (2018)).

Acknowledgements.

We thank L. Tinti for early conversations about how to handle off-diagonal anisotropies. M. Strickland was supported by the U.S. Department of Energy, Office of Science, Office of Nuclear Physics under Award No. DE-SC0013470.

Appendix A Q-functions

The Q-functions used in expanding the equations are defined as follows

[TABLE]

We note that the functions above functions are related, e.g.

[TABLE]

This fact allows us to reduce the number of underlying $Q$-functions that have to be computed to the “master function” $Q_{3}$.