Critical behaviour of hydrodynamic series

We investigate the time-dependent perturbations of strongly coupledN = 4 SYM theory at finite temperature and finite chemical potential with a second order phase transition. This theory is modelled by a top-down Einstein-Maxwell-dilaton description which is a consistent truncation of the dimensional reduction of type IIB string theory on AdS5×S. We focus on spin-1 and spin-2 sectors of perturbations and compute the linearized hydrodynamic transport coefficients up to the third order in gradient expansion. We also determine the radius of convergence of the hydrodynamic mode in spin-1 sector and the lowest non-hydrodynamic modes in spin-2 sector. Analytically, we find that all the hydrodynamic quantities have the same critical exponent near the critical point θ = 1 2 . Moreover, we establish a relation between symmetry enhancement of the underlying theory and vanishing the only third order hydrodynamic transport coefficient θ1, which appears in the shear dispersion relation of a conformal theory on a flat background. ar X iv :2 10 2. 03 58 4v 2 [ he pth ] 9 F eb 2 02 1


Introduction
The collisions of heavy ion at relativistic energies produce a hot and dense nuclear matter composed of deconfined quarks and gluons known as the strongly coupled quark-gluon plasma (QGP). The studies of the QGP at the Relativistic Heavy Ion Collider (RHIC) and the Large Hadron Collider (LHC) have led to the important results that explore a wide variety of QGP-related phenomena [1][2][3]. To describe in and out of equilibrium properties of this new phase at extreme conditions, Relativistic Hydrodynamics (RH) is a powerful tool [4][5][6]. The RH approach has great triumphs to explain the experimental events of the collider labs. Furthermore, the condensed matter physics has also benefited from the RH applications [7]. All of these evidences show the "unreasonable effectiveness" of the RH irrespective of the energy scales [8].
On the other hand, experimental observations imply that the QGP is a strongly interacting matter and the perturbative calculations cease to apply [2]. Therefore, non-perturbative methods such as the gauge/gravity duality may shed some lights on various properties of the QGP. The well-known example of the AdS/CFT correspondence states that the type IIB supergravity on the AdS 5 ×S 5 is dual to the 4−dimensional super Yang-Mills (SYM) theory living on the boundary of the AdS 5 [9][10][11]. Indeed, this duality is a strong-weak duality which maps a strongly-coupled quantum gauge field theory to a weakly-coupled classical gravity in one higher dimension. One of the main consequences of the AdS/CFT application is the prediction of universal ratio of shear viscosity to entropy density, η s = 1 4π , which agrees very well with the experimental data [12,13].
The RH approach is an effective theory for long wave-length regime in which each conserved quantity can be expressed in terms of a gradient expansion [4]. The coefficients of this series are the so called transport coefficients which contain the information of the underlying microscopic theory. One of the main concern regarding an infinite series expansions is the convergence features. The divergences which exist in the RH series in the real space [14] can be handled by using the Borel-Padé techniques [15][16][17][18][19]. This provides good information about the origin of divergent points and even the region in which the RH series might be convergent. On the other hand, the convergence of the gradient expansion in real space can be related to the radius of convergence of RH in momentum space [20,21]. In the momentum space, the location of the singular points reflect the existence of non-hydro modes and the credit of the RH depends on the strength of hydro modes over non-hydro modes [22,23].
Recently, the radius of convergence of hydrodynamics series have been investigated in various holographic cases. The radius of convergence of the shear-mode series in 4-dimensional AdS-Reissner-Nordstrom (AdS-RN) black brane has been studied in [24]. It has been extended to full range of the charge all the way to the extremal value both in 4 and 5-dimensional AdS-RN black branes in [25] (see also [26]) and different types of pole collisions at different values of the charge has been found. Likewise, the analytic properties of dispersion relations in all sectors have been studied for N = 4 SYM theory at zero chemical potential dual to AdS 5 black branes in [27,28].
In this paper, we consider the 1-R charge black hole (1RCBH) 1 model which is an analytical top-down string theory construction [29][30][31][32][33][34]. It is obtained from 5dimensional maximally supersymmetric gauged supergravity and is holographically dual to a 4-dimensional strongly coupled N = 4 SYM theory at finite chemical potential under a U (1) subgroup of the global SU (4) symmetry of R-charges. The interesting feature of this model is the existence of a second order phase transition at a certain point of the parameter space called critical point. This phase transition belongs to the model B dynamical universality class [35] and due to the large N c approximation, this phase transition shares same critical exponents as a mean-field model [36,37].
Different aspects of this background have been investigated. Indeed, it was shown that various quantities such as R-charge conductivity [38], complexity [39], mutual information [40] and entanglement of purification [41,42] remain finite, while their slopes diverge at the critical point with the same critical exponent θ = 1 2 . The value of the dynamical critical exponent was also confirmed in [43] and [44] by studying the non-hydrodynamical quasi normal modes (QNMs) of the external fields and quantum quench in this background, respectively.
Here, we would like to derive the transport coefficients for the dual theory of this background and benefit from the standard recipe [45][46][47][48]. We perform this calculation up to the third order of gradient expansion. In spin-2 sector we compute the relaxation time τ π and κ in the second order, and λ  16 in the third order of expansion from the Kubo formula. Furthermore, in the spin-1 sector we obtain another third order transport coefficient, 4 ), by investigating the dispersion relation of shear hydro mode.
According to the holographic dictionary, the hydrodynamic excitations correspond to the QNMs associated with poles of retarded Green's function of the conserved currents. Additionally, there is an infinite series of non-hydrodynamic QNMs in each sector resembling of the Christmas trees [49]. We investigate the radius of convergence of hydrodynamic series in the spin-1 sector for the whole range of the parameter space by computing the corresponding QNMs. In particular, we want to explore whether there is a critical behaviour and if so what is the associated critical exponent. We show that, depend on the value of the ratio µ T , different forms of mode collision can happen. Near the critical point there is a level-crossing between the hydro-mode and lowest non-hydro mode which is originally from the gauge field perturbation. While near the zero chemical potential the radius of convergence is determined by a level-crossing between the hydro-mode and the lowest non-hydro mode from the gravity perturbation.
To study the linearized equations of motion, enormous simplification will occur if we use the master equations formalism [50][51][52] which has two great advantages. First, we can investigate the set of equations analytically to compute higher order transport coefficients and also to find the radius of convergence of the shear hydrodynamic series. Second, it speeds up significantly the numerical computations of QNM in the spin-1 sector. The organization of this paper is as follows. In section 2 we review the building blocks of hydrodynamics and gradient expansion and introduce the transport coefficients up to the third order in gradient expansion for the conformal hydrodynamics. In section 3 we give the preliminary ingredients about the Green's function and how to derive the transport coefficients from the two-point functions in the AdS/CFT formalism via the Kubo formula. Section 4 is devoted to review the thermodynamics of the 1RCBH model and fix our notation. In section 5 we study the hydrodynamics of the 1RCBH model and illustrate the details for perturbations in spin-2 and spin-1 sectors. The second and third order transport coefficients as well as the dispersion relation of the shear mode are derived. In section 6 we compute the QNM frequencies in spin-2 and spin-1 sectors for complex momenta. We find the radius of convergence of the hydrodynamic series in shear channel using both analytical and numerical approaches. Likewise, in spin-2 sector the radius of convergence of lowest non-hydro modes is determined. In section 7 we investigate the behaviour of the transport coefficients and the hydrodynamic radius of convergence near the critical point of the phase transition and we show that all of them exhibit the same critical exponent θ = 1 2 . We conclude with a summary and an outlook to further directions in section 8.

Hydrodynamics as derivative expansion
In this section we will review the basic principles of the relativistic hydrodynamics of the boundary theory. The existence of the hydrodynamic equations is due to the conservation laws which are related to the continuous symmetries of the underlying theory. These symmetries yield the conserved quantities whose fluctuations are longlived and long distance propagating modes, i.e. ω → 0 as k → 0. If one considers a relativistic hydrodynamic including a local U (1) symmetry, then the hydrodynamic equations can be read as 2 where T µν is the energy-momentum tensor corresponding to the space-time symmetries and J µ is a conserved current corresponding to a local U (1) symmetry of the underlying theory. In order to solve the Equations (2.1), one may use the fact that T µν and J µ can be expressed as functions of the local temperature T (x ν ), the local four velocity u µ (x ν ) and the local chemical potential µ(x ν ) which are known as the hydrodynamic variables. Usually, the constitutive relations for the energy-momentum tensor and the current density are written in the Landau frame as where ( , p, n) are the equilibrium energy density, pressure and charge density, respectively. In flat space-time the operator ∆ µν ≡ η µν + u µ u ν is a projector perpendicular to the fluid velocity u µ and the flat metric η µν = diag(−1, 1, 1, 1). The Π µν and j µ tensors are dissipative contributions which can be expressed in terms of the derivatives of the hydrodynamic variables [4][5][6]. If the underlying theory possess the conformal symmetry, then the gradient expansion of the dissipative parts can be written as [53,54] Π µν = −ησ µν + ητ π Dσ µν + 1 3 3) The structure of the operators are chosen such that they transform homogeneously under the Weyl transformations [53]. Among the transport coefficients appeared in (2.3), some of them are more familiar such as the shear viscosity, η, the relaxation time, τ π and the conductivity, σ, while other ones are less recognized, like the non-linearized second order transport coefficients λ (2) i and the third order transport coefficients λ i . Due to the Landau matching condition, dissipative terms are transverse to the velocity profile u µ Π µν = u µ j µ = 0 and because of the conformal symmetry we have g µν Π µν = 0. The shear stress tensor is defined as σ µν ≡ ∇ µ u ν and for a given second rank tensor we use the following notation In addition, the vorticity field Ω µν and the convective derivative D are expressed as The third order operators O in Equation (2.3) are complicated Weyl-covariant tensors built out of the third order gradients whose detailed forms can be found in Ref. [54]. The corresponding third order transport coefficients λ (3) i are not uniquely defined but for a given underlying theory there exist at most 20 of them.
The transport coefficients (η, σ) of the 1RCBH model are calculated in the [38], while the other ones are still undetermined. It is shown that the shear viscosity and the conductivity exhibit a critical behavior close to the critical point of the second order phase transition [38]. A natural question is whether other hydrodynamic properties (e.g. the transport coefficients, the radius of convergence) of the underlying theory would exhibit similar behaviour. In the following section we present the basic ingredients and methods to derive the transport coefficients for N = 4 SYM theory at finite temperature and finite chemical potential by using the gauge/gravity duality. In section 5 we will benefit from our results to derive the corresponding second and third order transport coefficients.

Green's functions and transport coefficients
For a strongly coupled field theory which has dual gravity interpretation, the AdS/CFT correspondence has a recipe to derive the transport coefficients [45][46][47][48]. In this section, our main focus is to review the building blocks of this recipe. This includes the formalism to obtain the retarded Green's functions, the variational approach and the Kubo formula.

Retarded Green's functions and the holographic principle
To compute the n-point functions one needs to know the generating functional Z(J) of the theory where J is the source term. The n-point functions are nothing but the n'th functional derivatives of Z[J] with respect to J To our purpose, the two-point functions are of great importance. Consider a quantum field theory which has some gauge invariant operators, sayÔ i . Due to the causal structure of the relativistic hydrodynamics, we will focus only on the retarded two-point functions which are defined in momentum space as where k = (ω, − → q ). From the spectral representations, if one knows the retarded Green's function, then the two-point functions (advanced, symmetric and Feynman) will be easy to find [45]. In a holographic setup, there is a straightforward method to compute the retarded two-point functions of the boundary theory. Suppose that J is a source for an operator O in the boundary theory with an interaction Lagrangian density term JÔ. According to the Euclidean description of the AdS/CFT correspondence, we have the following relation where the left-hand side is the expectation value of the generating functional for the operatorÔ living on the boundary ∂M , while the right-hand side is the exponent of the classical on-shell action subject to the boundary condition Φ ∂M = J. Note that the Equation (3.1) is written in Euclidean signature while we are interested in Minkowski two-point functions. The procedure to make this connection is explained neatly in section 4 of Ref. [46] and we follow their recipe to compute G R ij (k).

Linear response theory
The transport coefficients can be determined via the variational approach and Kubo formula [4]. To quantify the hydrodynamics response, we use the conserved currents, T µν (x) and J µ (x), since they have microscopic definitions. One can use the variational approach to derive a specific two-point functions in the following steps [4]: 1. Write down the constitutive relations and expand them up to the first order in the hydrodynamic fluctuations (δT, δµ, δu µ ) and the source fluctuations (δA µ , δg µν ).
3. Plug the solutions into the constitutive relations (2.2) and make them on-shell.
4. Then, we use the on-shell currents to define the retarded Green's functions as Now we are in a position to obtain the transport coefficients. Each Green's function corresponds to a specific set of transport coefficients [4]. We are interested in a 3+1 dimensional hydrodynamic theory with a SO(3) symmetry in the spacial directions. Using this symmetry we can set the momentum along the z coordinate. To our purpose, we will focus on a so-called shear retarded Green's function G R xy,xy ≡ G R T xy T xy which has the following gradient expansion up to the third order [53,54] From this relation it is transparent that the transport coefficients can be found by taking special limits of G R xy,xy which is the well-known Kubo formula [4,53,54]. Here is the explicit form of the Kubo formula for each transport coefficient in this channel, These relations are universal which can be applied to any conformal field theory, including one with dual gravity description. For example, in the context of AdS/CFT correspondence and by using Equation (3.3) one can find the first and second order transport coefficients of the strongly interacting N = 4 SYM [45][46][47][48]53]. Although, some difficulties may arise such as implementing the holographic renormalization to remove the divergent terms [55][56][57], but the calculation is straightforward. In section 5, we will show how to utilize the holographic machinery as well as the Equation (3.3) to obtain the transport coefficients of a 4-dimensional CFT with second order phase transition which is dual to the 1RCBH backgrounds.

Holographic model: 1RCBH
In this section, we present a short review on the gravity setup which mimics a second order phase transition of the boundary theory. The bulk theory is a top-down string theory construction which is a consistent truncation of the super-gravity on AdS 5 ×S 5 geometry keeping only one scalar field and one gauge field coupled to the Einstein gravity.
The thermal solution is an asymptotically AdS black brane geometry with nontrivial profile of the scalar and gauge fields, so-called 1RCBH [29][30][31][32][33][34]. This geometry is dual to a four dimensional strongly coupled gauge theory N = 4 SYM at finite temperature and chemical potential.

The background geometry
The 1RCBH model is described by the following Einstein-Maxwel-dilaton (EMD) action where G 5 is the five dimensional Newton's constant. The self-interacting dilaton potential V (φ) and the Maxwell-dilaton coupling f (φ) are given by where L is the AdS 5 radius and without loose of generality we set L = 1 now-on. The corresponding equations of motion are The charged static stationary black brane solutions are given by wherer = 0 is the boundary of asymptotically AdS 5 background and M and Q are the mass and charge of the black brane, respectively. The black hole horizon is given by the largest root of the blackening function, h(r)|r =r h = 0, , and the gauge field A has only a time component which vanishes at the horizon due to regularity conditions. We use the following reparametrization to introduce the dimensionless parameters. In this manner the horizon radius can be fixed to one which leads to a simple relation between the mass and charge, namely Using the above dimensionless parameters simplifies our calculation presented in the following sections. We will use this parameterization now-on.

Thermodynamic
The Hawking temperature T of the black brane is given by which according to the gauge-gravity correspondence is equal to the temperature of the boundary theory. In addition, the chemical potential of the dual theory reads as It is straightforward to see that once the charge parameter vanishes, we will recover the geometry of the AdS 5 -Schwarzschild background as one may expect. One can characterize the 1RCBH model either by two non-negative parameters (M, Q) from gravity point of view or (µ, T ) from the boundary point of view. Nevertheless, the reparameterization introduced in Equation (4.4) can manifest an extra scaling symmetry in the system. Using Equations (4.6) and (4.7) one may make it more transparent to find the charge as a ratio of the boundary theory parameters, Since Q is real the above equation indicates that µ T ∈ [0, π/ √ 2]. Moreover, Equation (4.8) shows that for a given value of µ T ∈ 0, π/ √ 2 there are two distinct solutions corresponding to values of Q, while Q = √ 2 (µ = πT / √ 2) is the merging point of the two branches. We will show that thermodynamic quantities of the 1RCBH background diverge at the merging point declaring that this is the critical point of a second order phase transition point. In Ref. [38] it was shown that the solutions with −/+ sign in Equation. (4.8) are thermodynamically stable/unstable black branes. In this paper we are only interested in stable geometries and choose the branch corresponding to the − sign in Equation. (4.8). It turns out that introducing a new variable y as both simplifies our equations and makes the investigations close to the critical point more intelligible. Accordingly, one may express the charge Q, the Hawking temperature T and the chemical potential µ of the black brane in terms of this new dimensionless parameter It is easy to see that the critical point (Q = √ 2) and the AdS 5 -Schwarzschild background (Q = 0) correspond to y = 0 and y = 1, respectively.
In the context of the AdS/CFT correspondence, we have L 3 where κ 5 is the five dimensional gravitational constant given by κ 2 5 = 8πG 5 . By using the Equation (4.9) the Bekenstein entropy density for the black brane geometry (4.3) can be computed as As stated by the holography principles, the entropy density of the boundary theory equals to the entropy density of the black brane. Likewise, the R-charge density of the boundary model ρ = limr →0 δS δΦ is given by Having obtained the entropy and charge density with the aid of the thermodynamic equations s = ( ∂p ∂T ) µ and ρ = ( ∂p ∂µ ) T , one can derive the pressure of the dual strongly coupled theory p N 2 The conformal symmetry of the underlying boundary theory leads to a relation between the energy density ε and pressure, namely ε = 3p which is nothing but the tracelessness of the stress tensor. On the other hand, the Hessian matrix of the boundary thermodynamic quantities is , whose Jacobian J , the determinant of the above matrix, is given by Clearly, the Jacobian is divergent at y = 0 ( µ T = π √ 2 ) signaling that there is a second order phase transition. The thermodynamics of the 1RCBH model is reviewed comprehensively near the critical point, y = 0, in the Ref. [43]. Moreover, it is easy to show that at y = 1 we recover the thermodynamics of the N = 4 SYM at zero chemical potential.

Hydrodynamics dual to the 1RCBH model
We will extend the previous studies of the hydrodynamics dual to the 1RCBH model [38] to the higher orders of gradient expansion. There are two distinct computational ways to study the hydrodynamics of the current model which, in principle, at the end of the day should lead to the same results. The first method, is to use the gauge fixing process, following KS [49], and classify the perturbations according to the diffeomorphism and gauge transformations. On top of that, the little SO(D − 2) group can be used to specify various sectors of the perturbations. The second method, is to adopt the master equations formalism [50][51][52]. 3 Unlike the KS method, in master formalism one should write the gauge invariant combination of all perturbations before using any gauge fixing which leads to gauge invariant linearized equations. There exist particular combinations of fields such that the equations of motion reduce to the Schrodinger type equations. While for spin-2 sector there is no difference between those methods, the latter method has great advantages and simplifies vastly the equations in other sectors. In this paper, we study the hydrodynamics in spin-2 and spin-1 sectors, utilizing the master formalism [50][51][52], and leave the spin-0 sector to the future works.
We close this part by addressing the importance of total action. As stated in section 3, the transport coefficients are derived from two-point functions which are nothing but the derivatives of on-shell total action with respect to the sources. The total action is S tot = S bulk + S GH + S ct , 3 We thank Andrzej Rostworowski for constructive discussion on master equations formalism.
where S bulk is the bulk action given in Equation (4.1), S GH is the Gibbons-Hawking boundary term and S ct is the counter term [58], where γ ij is the boundary induced metric, K is the corresponding extrinsic curvature, K = 1 2 γ ij ∂ r γ ij . We would like to emphasize that the S ct contributes to the second order transport coefficients κ and τ Π , while it has no impact on the first and third order ones.

Spin-2 sector
We start by investigating the spin-2 sector which includes only one, namely perturbation h xy (t, z, r) = g xx (r)H xy (r)e −iωt+iqz in which we use the SO(3) symmetry of the background to fix the momentum along z direction. The linearized Einstein equation for this perturbation in Schwarzschild coordinate is given by where we use the dimensionless frequency and momentum defined as The Equation (5.2) is a subtle linear second order differential equation which is very hard to find the exact analytical solutions. However, we are interested in the hydrodynamic limit, i.e. w 1 and q 1. The solution representing the ingoing wave at the horizon r = 1 can be written as (for more details see Appendix A), where H(r) is regular at the horizon. In hydrodynamic limit we use the following ansatz to solve Equation (5.2) order by order in . The explicit form of the solutions up to third order is given in Equation (A.11). The next step is to compute the on-shell action and keep all the terms quadratic in H xy following the recipe explained in Ref. [46] which leads to We would like to highlight a practical point at this stage in our computation. The contribution of the counter term in (5.4) starts from the second order which is given in the second line of Equation (5.4). On the other hand, each term in the counter term (5.1) is composed of even number of derivatives. Therefore, the general form of S tot is the same in 2k and 2k + 1 order of expansion. Since we are interested in computing the transport coefficients up to the third order, we only need to plug the solution of H xy to the corresponding order in (5.4). The retarded Green's function can be obtained by taking the second derivative of this total action with respect to the source [46], 17 ) .
Now we are in the position to apply the Kubo formula (3.3) and compute the transport coefficients up to the third order of gradient expansion, It is easy to check that the shear viscosity satisfies the universal relation which is in complete agreement with [13,38]. Moreover the y = 1 limit of those transport coefficients coincide with the counterparts of the N = 4 SYM theory dual to AdS 5 black brane background [53,54]. One of the remarkable features of these transport coefficients is their behaviour near the critical point of the phase transition. In the section 7 we will show that all the transport coefficients reach to their critical value with the same critical exponent θ = 1 2 .

Spin-1 sector
In this section we study the hydrodynamic limit corresponding to the spin-1 sector. In our setup, by employing the master equations formalism [52] the coupled equations of motion of the gauge invariant perturbations can be written in a decoupled form. Indeed the chief advantage of using this approach is twofold. First, it authorizes the analytical investigations of this sector, as we will discuss in this section. Second, the numerical computation of the QNM frequencies will be less costly. We will come back to this point in section 6. Let us summarize the key steps of the master equations formalism [50][51][52] in the following. At the first step, we should find the gauge-invariant combinations of the perturbation in spin-1 sector with plane wave ansatz Φ(r)e −iωt+iqz for all the fields The next step is to write the linearized equations of motion for the gauge invariant perturbations Z i (i = 1, 2, 3) which are coupled equations. At the final step we should rewrite the gauge invariant combinations as linear combinations of the master scalars and their derivatives such that they satisfy master equations. The master scalars corresponding to the spin-1 sector are defined by Z 1 (r) ≡ r 2 h(r)e 3A(r)−B(r) q 2 (3A (r)Ψ 2 (r) + Ψ 2 (r)) , and the master equations are given by where Ψ i ≡ Ψ i (r)e −iωt+iqz , the d'Alembert operator is the wave operator on the background and the potential matrix W is given in terms of the background functions.
While in general the final equations are still coupled they have a simple form in terms of the master scalars. Nevertheless, in our case, since the following condition holds The master equations (5.6) can be further simplified to decoupled equations as This form of master equations make our analytical studies more intelligible. It turns out that the equation with " + (−)" sign describes the shear channel perturbation (transverse perturbation of the gauge field) . In the rest of this section we solve Equation (5.6) with " + " sign perturbatively to find the hydrodynamic dispersion relation of the shear mode following the prescription given in Ref. [46]. We solve the relevant equation in Eddington-Finklestein (EF) coordinate by using the following ansatz Ψ + (r) = r 1/2 n=0 n ψ n (r), and we scale the frequency and momentum as (w, q) → ( 2 w, q). Again we utilize the variation of parameters method reviewed in Appendix A to find the ψ n 's. At the zeroth order in one can solve a second order differential equation for ψ 0 with regularity condition at the horizon ψ 0 (r) = M r where the second integration constant is fixed by imposing ψ 0 (1) = 1. In higher orders, without loose of generality, we impose ψ i (1) = 0 for i > 0. We find the solutions analytically up to the second order and the source-less boundary condition leads to the following form of the spectral curve [28] in shear channel By solving this equation perturbatively in q, one can find the dispersion relation of the shear hydrodynamic mode where we use Q = 2(1−y) 1+y . By comparing our results with general form of the dispersion relation in the shear channel for conformal theories (see for example [54]), one can find a relevant third order transport coefficient θ 1 as Note that at the critical point of the phase transition the θ 1 vanishes. This could be a hint to the symmetry enhancement of the underlying theory at the critical point. One may compare this phenomenon with vanishing bulk viscosity for theories with conformal symmetry. We will elaborate on this point in the discussion section. Let us close this section by emphasizing that without using the master equations formalism [52] the above computation could not be achieved in this straightforward manner (as pointed out in a simpler case such as AdS-RN black holes [26]). In principle, one can continue and find the higher order terms in the dispersion relation by solving ordinary second order differential equations for higher ψ n using the variation of parameters method.

QNMs and convergence of the hydrodynamic series
In this section we compute the QNM frequencies in spin-2 and spin-1 sectors associated with the poles of the corresponding retarded Green's function of the boundary theory. We consider the complex momentum square q 2 = |q 2 |e iϕ and compute the QNM frequencies.
That is simply because of the symmetries of the static background which leads to the fact that the linearized equations are functions of even powers of momentum.
To compute the QNMs in our background we used the pseudo spectral Chebyshev discretization along the radial coordinate in EF parameterization. The boundary conditions we should impose are ingoing wave at the horizon which translates to the regularity condition in EF coordinate, and sourceless Dirichlet boundary condition at the asymptotic region which comes from the holographic dictionary. In all the cases we consider in this work, the 30 number of grid points along the radial coordinate are Figure 1. The trajectory of the lowest QNMs in spin-2 sector at the critical point (y = 0). In the first row we show the trajectories before and after the first collision between the lowest modes at |q 2 * | = 1.505 and w * = −1.637i. In the second row we show the trajectories before and after the second collision at |q 2 * | = 2.621 and w * = ±1.551 − 3.081i.
enough to find ∼ #10 lowest eigenmodes and we compare the results with 60 number of grid points to select the reliable modes with 10 −6 relative accuracy. To cross-cheek our numerical findings we compute the dispersion relation of the hydro mode in small real momenta in shear channel and compare it with the analytical result presented in Equation (5.7) which shows perfect agreement. Finding the radius of convergence for series is an interesting topic in complex analysis [59]. Suppose we have an analytic (spectral) curve T (u, v) = 0 of complex variables (u, v) in the complex plane and we want to find v = v(u). These solutions are classified as regular points and critical points. Regular points are the zeros of curve which ∂ i T ∂v i | v=v (u) = 0 for i ≥ 1, while critical points of order j are the zeroes which ∂ i T ∂v i | v=v (u) = 0 for 1 ≤ i ≤ j. By definition, the radius of convergence for v is determined by the location of the nearest critical point to the origin. Therefore, at the critical point there is a degeneracy of solutions which specifies the radius of convergence. As explained in details in [23], if the spectral curve is non-analytic at some points then there may be other sources of singularities. In all the cases studied in this paper we did not find any sign of non-analyticity in the spectral curves for momenta smaller or equal to the critical momentum. 4

Spin-2 sector
In this section we consider the perturbations in the spin-2 sector and compute the corresponding QNMs numerically. Since there is no hydro-mode in this sector we can only find the radius of convergence of the non-hydro modes and we will focus only on the lowest modes. In small values of |q 2 | each mode has a closed trajectory for ϕ from 0 to 2π while for larger values they may collide and share their trajectories. Due to the symmetries the collision between the pair modes is always on the imaginary axes in the complex frequency plane with a purely imaginary momentum. The trajectories of the QNMs are qualitatively the same for the whole range of y. In Figure 1, as an example, we illustrate the trajectories of the modes before and after the first collision for y = 0. The general features of the lowest collision in this sector does not change in the whole range of y. In Figure 2 we show the radius of convergence of the lowest mode in the spin-2 sector and also the corresponding frequency as functions of y. Close to the critical point the radius of convergence and the corresponding frequency are linear functions in y which can be fitted by In section 7 we will discuss on these relation. 4 We thank Michal Heller for pointing this out.

Spin-1 sector
In this section we consider the perturbations in the spin-1 sector, including the shear and the transverse gauge field channels, and compute the lowest QNMs using numerical techniques.
In small values of |q 2 | each mode has a closed trajectory for ϕ from 0 to 2π. In other words, in this regime each mode can be found uniquely in complex momentum square plane. On the other hand for larger values of |q 2 | this may change due to level-crossing or collision of the modes [27]. In fact, depends on how far the medium is from the critical point either of phenomena may happen. In the rest of this section we present various plots of the trajectory of the QNM frequencies for the complex momentum square and discuss the main features at different regimes. One of the purposes is to investigate the convergence radius of the hydrodynamic series in the whole regime of our model. In particular its behaviour near the second order phase transition will be explored in the next section.
To compute the QNMs of the system in complex momentum square, we employ the master equations formalism [52] explained in section 5.2. The decoupled linearized equations are given in Equation (5.6) and as pointed out in [24,25] there is a collision of the modes due to the appearance of a radical in the last term of the potentials W ± . If the lowest collision of the hydrodynamic mode occurs at this |q 2 |, then the radius of convergence can be related to this phenomena, Figure 4. The trajectory of the lowest QNMs in shear channel (solid lines) and in the transverse channel of the gauge field (dashed lines) for y = 1 (AdS-Schwarzschild black brane) in the first row and for y = 0.8 in the second rows. Top-left (right) panel: the trajectory of the lowest modes before (after) the joining of the hydro mode to the lowest non-hydro shear mode. The joining happens at |q 2 * | = 2.224. Bottom-left (right) panel: the trajectory modes before (after) the collision occurs at |q 2 * | = 2.553. The collision of hydro and the lowest non-hydro shear mode happen leads to exchange of their trajectories and not joining to a common trajectory. This is so called level-crossing [28].
It is easy to show that the critical momentum given in Equation (6.2) has a simple expression in terms of the thermodynamic variables This branch point is of square-root type and at this value of the momentum the shear and transverse gauge modes satisfy the same equations of motion. Therefore, for each mode in the former there is a cousin in the latter which they meet at the critical momentum. In other words, the convergence radius of the hydrodynamics is constrained by a crossing between the hydrodynamic mode and one of the modes in the transverse gauge channel. Figure 5. The trajectory of the lowest QNMs in shear channel for various values of the amplitude of the complex momenta square |q 2 | for y = 0, 0.2. We use the same colors as Figure 4 to keep track of the shear modes and transverse gauge field modes in the coupled system. The joining happens at w = 0.244i and w = 0.370i respectively.
Our main results can be summarized in Figure 3 in which the y dependency of the convergence radius is demonstrated, where it indicates that there exist at least two different structures for the first mode collision of the hydrodynamic QNM. In the left panel, the red and blue points are what we found by studying the QNMs of the system in q 2 numerically while the solid green line is the analytical formula for the convergence radius (6.2) computed by employing the master equations formalism [52]. In the right panel we show the corresponding frequency for 0 ≤ y ≤ 0.596 (solid green line) and the linear fit near critical point given by −iw = 0.244 + 0.531 y. (6.4) There are couple of intriguing regimes of the y parameter that we will investigate separately in the following. Let's start by recomputing the modes for AdS-Schwarzschild black brane corresponding to y = 1. In this case we can compare our results with the one which is already known [28]. We also introduce our conventions for the presentations of the trajectories of the modes. In the first row of Figure 4, we show the trajectories of the lowest modes in spin-1 sector as functions of complex momenta for y = 1. In this case, for real momenta the QNM frequencies are two Christmas trees including one for the shear modes and one for the transverse gauge modes [49]. Interestingly, the first level-crossing occurs between the latter modes at |q 2 | = 0.479, while the level-crossing between the hydro and first shear non-hydro mode occurs at |q 2 | = 2.224 which are in agreement with [28]. Therefore, by definition the radius of convergence of the hydrodynamic is |q 2 * | = 2.224. Also note that for real momenta the non-hydro modes from two channels are always in pair and in the regime we are interested in, they are such that the one from shear perturbation is less damped compare to the one from transverse gauge field perturbation. This will be the case all the way to y = 0. In all plots presented in this section, we use the same colors for the modes and their trajectories. The corresponding modes, for the real momenta are shown by dots and the positions of mode collisions are shown by crossed-circles.
While close to y = 1 the convergence radius of the hydrodynamic series is due to a level-crossing between the hydro mode and the lowest non-hydro shear mode, for 0.596 < y < 0.89 it's not the case anymore. The level-crossing will be replaced by a collision and the hydro-mode will still have a closed trajectory. 5 This is demonstrated in the second row of the Figure 4 for y = 0.8. For this background the radius of convergence of the hydrodynamic series is |q 2 * | = 2.553 which is associated with a mode collision at q 2 * = 2.10475 ± 1.44493i and w = ±1.48974 − 1.06122i. As we show in Figure 3 the radius of convergence has a maximum value at y = 0.596 and therefore one may expect that around this point something interesting could happen. At y = 0.596 the convergence radius of the hydrodynamic series is |q 2 * | = 2.854 where two lines are crossed. Associated to this point, there is a collision between the hydro mode and non-hydro shear modes at q 2 * = 2.28767 ± 1.70642i and w = ±1.5247 − 1.0434i as well as a level-crossing with transverse gauge field non-hydro mode at q 2 * = −2.854 and w = 0.86i. The most important point here is that the convergence radius of hydrodynamic series is found due to multi-phenomena: a level-crossing and a collision between the QNMs simultaneously.
Finally, in Figure 5 we show the trajectory of the least damped modes in spin-1 sector as functions of complex momenta at the critical point of the second order phase transition y = 0 in the first row and y = 0.2 in the second row. Again, the first level-crossing happens between the lowest non-hydro transverse gauge field modes. This is always true for any value of 0 ≤ y ≤ 1. The radius of convergence of the hydrodynamic series is |q 2 * | = 9 8 for y = 0 and |q 2 * | = 1.47 for y = 0.2. Note that in this regime there are always infinite number of level-crossing between the shear and gauge field modes at the |q 2 * |, since the linearized equations are the same at this critical value of the momenta, see the Equation (5.6). The radius of convergence of the hydrodynamic series at the critical point, given by y = 0, is almost half of its counterpart N = 4 SYM without chemical potential, corresponding to y = 1.
Last but not least, in Figure 5 there is a second damping mode which can be seen in small y cases, as one of the lowest QNMs. This observation may give us a hint about the QNM structure in the whole range of parameter space, in particular for y = 1. As we already discussed at the critical value of momentum given in Equation (6.2) two equations of motion in spin-1 sector coincide with each other. In other words, for each shear mode there should be a cousin in the transverse gauge field mode. On the other hand, due to the symmetries of the background and QNM equations in each sector, either the modes are purely imaginary or they are in pair as ω ± = ±ω R − iω I . Having said that, we can conclude that there should be odd numbers of purely imaginary modes in the transverse gauge field channel to accomplish the fact that there are odd numbers of modes in the shear channel (one hydro and infinite number of pair non-hydro modes). Based on this argument we infer that in AdS-Schwarzschild black brane geometry the transverse perturbation of a gauge field should have at least one purely imaginary QNM which has not been addressed in the literature to our best knowledge.
To summarize our results, for the convergence radius of the hydrodynamic series in shear channel we would like to point out that there is a competition between the lowest non-hydro modes in the shear channel and the ones in the transverse channel to join to the hydro-mode. In the range 0 ≤ y < 0.596 which ends to the phase transition point the latter is responsible, while in 0.596 < y ≤ 1 the former takes the main role to specify the radius of convergence of hydrodynamic series.

Critical exponent
As mentioned before, the theory we study in this paper enjoys a critical point at µ = πT / √ 2. The behavior of different observables near this critical point has been studied in [38][39][40][41][42][43][44]. In all cases it was shown that the critical exponent is θ = 1 2 . All quantities remain finite (except the diffusion constant which vanishes), while their slopes diverge at the critical point. By using the chain rule, it is easy to show that if a quantity P has the following expression near the critical point in terms of the dimensionless parameter y = 1 − 2µ 2 then, the behaviour of the associated quantity near critical point at either fixed temperature or fixed chemical potential will be given by, Let us note that, although the background (4.3) has linear expansion close to any value of y 0 , namely the physical quantities associated to the boundary theory do not inherit this behaviour necessarily.

Conclusion
After reviewing general aspects of hydrodynamics as a gradient expansion in section 2, Green's function and linear response theory in section 3, we focused on a specific example from section 4 onwards, namely the 4-dimensional N = 4 SYM theory at finite temperature and finite chemical potential which poses a second order phase transition at µ = πT / √ 2. This theory is dual to a Einstein-Maxwell-dilaton theory which is constructed by a consistent reduction and truncation of 10-dimensional type IIB supergravity on AdS 5 ×S 5 , keeping only the metric, one scalar field and one gauge field [29][30][31][32][33][34]. The black hole solutions with planar horizon and with sourceless scalar field boundary condition in this theory are so-called 1RCBHs which are asymptotically AdS. In sections 5-7 we computed the linear hydrodynamic transport coefficients up to the third order in spin-1 and spin-2 sectors as well as the radius of convergence of the hydrodynamic series in shear channel for arbitrary temperature and chemical potential.
Surprisingly we found that close to the critical point of the phase transition all the quantities share the same critical exponent θ = 1 2 which might be related to the large N c limit of our setup [36,37]. 7 We speculate that the other hydrodynamic transport coefficients, including the non-linear ones, exhibit the same critical behaviour. To support our proposal let us recall a universal relation among the second order transport coefficients [61,62], which we expect to be hold in our setup too.
As one of our main results we found that the only third order transport coefficient which appears in shear hydrodynamic dispersion relation dies out at the critical point, 7 Although the Einstein-dilaton theory with second order phase transition exhibit a universal critical exponent θ = 2 3 in thermodynamic quantities in various dimensions [60]. θ 1 = 0. This peculiar observation might be related to the symmetry enhancement of the underlying theory at the critical point, in a same spirit that the bulk viscosity vanishes due to conformal symmetry. Let us emphasize that both analytical and numerical studies in the spin-1 sector have been carried out by employing the master equations approach [50][51][52] to solve the linearized equations of motion. To our best knowledge this is the first non-trivial example of Einstein-Maxwell-dilaton theory which benefits enormously from master equations formalism. We analyzed the QNM structure in the spin-1 sector with complex momentum to compute the convergence radius of the hydrodynamic series. We found that in different regime different mode collisions and/or level-crossing is responsible for radius of convergence. In a range of parameters which ends to the critical point (0 ≤ y ≤ 0.596) we found a perfect agreement between our numerical results and the analytical findings which is a consequence of having the same linearized equations for both shear channel and transverse channel at special momentum. Therefore, in this regime the lowest level-crossing that fixes the hydrodynamics radius of convergence is between the shear hydro mode and the lowest transverse non-hydro modes. On top of that we found that the minimum radius of convergence occurs at the critical point. On the other hand, in the second regime which ends to the zero chemical potential (0.596 ≤ y ≤ 1) we found that a mode collision or a level-crossing between the shear hydro mode and lowest shear non-hydro modes lead to the radius of convergence of the hydrodynamic series. For completeness, we have also studied the spectrum of QNMs in spin-2 sector. Due to the fact that there is no hydro mode in this sector, the radius of convergence is set by the collision of non-hydro modes and we have analyzed it in the entire range of the parameter space. Now we would like to address some generalizations of our results. As we already mentioned, the only third order transport coefficient which shows up in shear dispersion relation, θ 1 , vanishes at the critical point. This suggest that the full hydrodynamic expansion may contain less number of transport coefficients at the critical points, even at the higher order of expansion. In this manner, further investigations for other cases deserve more attention in future studies. To complete our results, the next step is to study the spin-0 sector of the same model which includes three coupled linearized equations. With some efforts we found numerically that the other third order transport coefficient θ 2 ≡ −(λ 6 ) which appears in the sound hydro mode dispersion relation [54] does not vanish at the critical point of the phase transition. It turns out that a complete investigation of this sector with complex momentum is much more involved and has some subtleties that we leave this study for future work.
It could also be of interest to venture other models with second order phase transition, e.g. bottom-up holographic cases [60,[63][64][65][66], to study the radius of convergence of the hydro-mode and to check our proposal of a relation between the symmetry enhancement of the underlying theory and vanishing higher order transport coefficients at the critical point.
Finally, it will be interesting to consider the large but finite coupling impacts to our findings with an ultimate goal of interpolating between weak and strong coupling results. This can be achieved for example by adding the Gauss-Bonnet term to the action. The QNMs and some of the transport coefficients of such a theory at zero chemical potential have been already computed in [67].