Anomalous transport from holography: Part I

We revisit the transport properties induced by the chiral anomaly in a charged plasma holographically dual to anomalous $U(1)_V\times U(1)_A$ Maxwell theory in Schwarzschild-$AdS_5$. Off-shell constitutive relations for vector and axial currents are derived using various approximations generalising most of known in the literature anomaly-induced phenomena and revealing some new ones. In a weak external field approximation, the constitutive relations have all-order derivatives resummed into six momenta-dependent transport coefficient functions: the diffusion, the electric/magnetic conductivity, and three anomaly induced functions. The latter generalise the chiral magnetic and chiral separation effects. Nonlinear transport is studied assuming presence of constant background external fields. The chiral magnetic effect, including all order nonlinearity in magnetic field, is proven to be exact when the magnetic field is the only external field that is turned on. Non-linear corrections to the constitutive relations due to electric and axial external fields are computed.


Introduction and summary
Hydrodynamics is an effective long-distance description of most QFTs at nonzero temperature. Within the hydrodynamic approximation, the entire dynamics of a microscopic theory is reduced to that of macroscopic currents, such as of charge current operators computed in a locally near equilibrium thermal state. An essential element of any hydrodynamics is a constitutive relation which relates the macroscopic currents to fluid-dynamic variables, such as charge densities, and to external forces. The most simple example of constitutive relation is the diffusion approximation for the electric current J J = −D 0 ∇ρ (1.1) Chiral anomalies emerge and play an important role in relativistic QFTs with massless fermions. The chiral anomaly is reflected in three-point functions of currents associated with global symmetries. When the global U (1) currents are coupled to external electromagnetic fields, the triangle anomaly renders the axial current into non-conserved, where J µ /J µ 5 are vector/axial currents, and κ is the anomaly coefficient. In SU (N c ) gauge theory with a massless Dirac fermion in fundamental representation, κ = eN c /(24π 2 ) and e is electric charge. Here E, B ( E a , B a ) are vector (axial) external electromagnetic fields 1 .
The chiral separation effect (CSE) [36,37] is another interesting phenomenon induced by triangle anomalies. It is reflected in separation of chiral charges along external magnetic field at finite density of vector charges. Chiral charges can be also separated along external electric field, when both vector and axial charge densities are nonzero, the so-called chiral electric separation effect (CESE) [38,39].
Without external fields, triangle anomalies affect transport properties through hydrodynamic flows. Particularly, there exists an anomaly induced chiral vortical effect [40][41][42], which relates the current to fluid's vorticity. In the fluid's local rest frame, the chiral vortical effect is J = 1 2 ξ ∇ × u, where u is the fluid velocity. The transport coefficient ξ was first calculated in [40,41] using the fluid/gravity correspondence [43][44][45]. Later, in [42] it was shown that the chiral vortical term is required by existence of a positive-definite entropy current associated with the hydrodynamic system and furthermore that ξ is uniquely determined by the anomaly coefficient κ.
In [54] we derived the most general off-shell constitutive relation for a globally conserved U (1) current driven by non-dynamical external electromagnetic fields. The derivation involved a resummation of all-order gradient terms in U (1) current. The gradient resummation was implemented via the technique of [55][56][57][58]. The latter was devised to resum, in linear approximation, all-order velocity derivatives in the energy-momentum tensor of a 1 A possibility of experimentally creating axial electromagnetic fields in a laboratory was recently discussed in [1]. 2 See also [7][8][9] for earlier related works.
holographic conformal fluid. For a holographically defined theory involving a probe Maxwell field in the Schwarzschild-AdS 5 geometry, the constitutive relation for the boundary current was found to be parameterised by three momenta-dependent transport coefficient functions: diffusion and two conductivities. The key element in the derivation was "off-shellness" of our method. That is the transport coefficients were uniquely determined via solution of dynamical components of the Maxwell equations in the bulk. The constraint component was shown to be equivalent to continuity equation of thus derived current.
In the present work we extend the study of [54] and account for the effects induced by the triangle anomaly: when the triangle anomaly is present for both left/right-handed chiralities, we derive off-shell constitutive relations for vector/axial currents. As mentioned above, anomalies contribute to the stress-energy tensor [40,41]. However, as in [54] we chose to work in the probe limit in which the currents and stress-energy tensor decouple. In the dual gravity, the probe limit ignores the backreaction of the gauge dynamics on the bulk geometry. The holographic model in study consists of two Maxwell fields in the Schwarzschild-AdS 5 black brane geometry. The triangle anomaly is holographically modeled via the gauge Chern-Simons action for both Maxwell fields (with opposite signs). This holographic setup can be realised via a top-down brane construction of D4/D8/D8 [59].
Below we will consider the charge densities (chemical potentials) as constant with small inhomogeneous fluctuations on top: whereρ,ρ 5 ,μ,μ 5 are constant backgrounds and δρ, δρ 5 , δµ, δµ 5 are the inhomogeneous fluctuations. The parameter is formally introduced as a small parameter. It will be used to set up a perturbative procedure. We will be particularly interested in linearisation in inhomogeneous fluctuations and most of our results will be accurate up to first order in . The charge densities ρ, ρ 5 are the hydrodynamic variables. They can be related to corresponding chemical potentials via (3.7). For constant parts,μ =ρ/2 andμ 5 =ρ 5 /2. Our study is divided into two largely independent parts. In the first part, the external fields are assumed to be weak and scale linearly with To first order in , we are able to derive the most general all order off-shell constitutive relations for the vector/axial currents.
In the second part of our work, the external fields are assumed to have constant background values E, B, E a , B a plus small inhomogeneous fluctuations δ E, δ B, δ E a , δ B a , (1.5) We will see that the constant backgrounds induce interesting nonlinear anomaly-induced structures in both currents.
Our study goes in several directions beyond the results reported in the literature.
• For the linearised setup, we present a rigorous derivation of off-shell constitutive relations for the vector/axial currents, coupled to non-dynamical external vector/axial electromagnetic fields. Apart from the well-known chiral magnetic/separation effects, we obtain two additional anomaly-related transport coefficients, magnetic conductivity and its axial analogue. Furthermore, all the transport coefficients are generalised to momenta-dependent functions as a result of exact all order gradient resummation.
• Beyond the linear regime, we calculate some nonlinear effects induced by constant background external fields. While the phenomena that we discover are largely not new, most of them have not been reported for the present holographic model. In the absence of external electric and axial fields, the chiral magnetic/separation effects are proven to be exact for arbitrary strong constant magnetic field.
• In a follow-on publication [60], we will extend the study of non-linear CME to spatially-varying magnetic field and will evaluate derivative corrections to it in the constitutive relation. Furthermore, we will consider a case when the axial chemical potential is dynamically generated through E · B term for the case of constant magnetic field and a weak time-dependent electric fields. Such setup is experimentally realisable in condensed matter systems 3 . Dependence of AC conductivity on magnetic field is in the focus of this study.
In the first, "linear", part of our study, all-order derivatives are resummed into the following constitutive relations for the vector/axial currents 4 where D, σ e/m , σ χ/κ and σ a are scalar functionals of spacetime derivative operators Transforming to the Fourier space via the replacement (∂ t , ∂) → (−iω, i q), these scalar functionals become functions of the frequency ω and momentum squared q 2 . We refer to these functions as transport coefficient functions (TCFs). Here ω and q are dimensionless, while the dimensionfull momenta are πT ω and πT q with T being the temperature. The axial current is TCFs contain information on infinitely many derivatives (and transport coefficients) in conventional hydrodynamic expansion. The latter is recovered via small momenta expansion. While most of the results on transport coefficients reported in the literature are obtained order-by-order in the expansion, our results are exact to all orders. Furthermore, the TCFs account for collective effects of non-hydrodynamic modes, which never emerge in the strict low frequency/momentum expansion. The diffusion D, electric conductivity σ e and magnetic conductivity σ m were studied in [54]. The additional three TCFs are induced by the anomaly: σ χ is the momenta-dependent chiral magnetic conductivity [10]; σ κ generalises the chiral separation effect [36,37]; σ a is an axial analogue of the magnetic conductivity σ m .
(1.6,1.7) can be equivalently presented using the chemical potentials µ, µ 5 instead of the charge densities When the triangle anomaly is switched off (κ = 0), (1.8) for J µ coincides with the constitutive relation of [61] derived from an effective action.
In the present holographic model, we succeeded to uniquely determine all the TCFs in (1.6,1.7). As for (1.8,1.9), the coefficients α 1 , α 2 , D , σ e/m are in fact frame-dependent. We postpone further discussion until section 4 where these TCFs are presented in a certain frame (see (4.20,4.21)).
In the hydrodynamic limit ω, q 1, the TCFs are series expandable: σ a = 144κ 2μμ 5 (2 log 2 − 1) + · · · , (1.14) We were unable to obtain an analytical result for the anomalous correction to iω-term in σ m . In section 4.2.2 we will reveal that anomalous correction to iω-term in σ m is linear in (μ 2 +μ 2 5 ). Turning the anomaly off, the magnetic conductivity σ m coincides with that of [54]. Interestingly, σ m is being corrected by the anomaly 5 . Appearance of anomalous corrections in σ m could be explained as the effect of two triangular anomaly-generating Feynman diagrams inserted in the current-current correlator 6 [63].
In [54], we made a full comparison of D, σ e and anomaly-free part of σ m with known results in the literature, particularly with the ones that could be extracted from the currentcurrent correlators. The underlined terms in D, σ e/m cannot be fixed from the correlators, 5 This fact was previously noticed in [62] where the authors went beyond the probe limit and included backreaction of the bulk gauge fields onto the geometry. However, taking the probe limit in [62] does not seem to coincide with our results. The reasons behind this discrepancy remain unclear to us. 6 We thank Ho-Ung Yee for bringing this argument to our attention.
while the constant pieces of σ m , σ a can be. In [64], the constant piece σ 0 m was evaluated for a pure QED plasma with one Dirac fermion at one loop level: the result was found to be positive and was interpreted as an anti-screening of electric currents in the plasma medium. In contrast, in [65] σ 0 m was argued to be zero based on Boltzmann equations. In strongly coupled regime, to our best knowledge, σ m , σ a have not been reported in the literature.
In the same holographic model as considered here, the constant terms in σ χ and σ κ were originally presented in [29]. We find full agreement with those results. The constitutive relations (1.6,1.7) can be used to derive Kubo formula (4.28) for σ χ/κ , which comes out to be consistent with [10]. Beyond the constant terms, the analytical results in σ χ/κ are new as far as we can tell.
Away from the hydrodynamic limit, the TCFs are known numerically only: D, σ e and the anomaly-free part of σ m were already reported in [54]; numerical results for σ m , σ a , σ χ , σ κ will be presented below, in section 4.2.2. Comparison of σ χ with the results of [24,66,67] will be discussed as well.
In the second part of our study, new nonlinear structures emerge in the constitutive relations for both currents. The complete set of results will be displayed in section 5. Meanwhile we will only flash the results with the axial background fields turned off. To zeroth order in fluctuations, the vector/axial currents are where the subscript (0) denotes zeroth order in fluctuations (O( 0 )). V i (1) are functions of E, B, µ, µ 5 , and are perturbatively computed in section 5.1. When E = 0, (1.16) confirms exactness of the CME [2][3][4][5][6] for arbitrary constant magnetic field, in agreement with [68,69]. Both µ, µ 5 depend on E, B and these dependences account for nonlinearity of the CME with respect to external fields. The electric field E introduces corrections to the original form of the CME and is a source of new structures.
In the weak field limit, V j (1) are expandable in the amplitudes of E, B. We quote the results up to third order: (1.18) where µ, µ 5 are also expandable in E and B, where · · · are terms of higher powers in E and B. The second order term E × B is the Hall current induced by the anomaly, referred to as chiral Hall effect in [70]. In classical electromagnetism, the Hall current is generated by the Lorentz force and the imposition of steady state condition. However, the usual Hall current is generated non-linearly and cannot be generated in a holographic model with the Maxwell action only, which does not induce any nonlinearity. Beyond the probe limit, the Hall current does emerge [62,71].
The last term in (1.18) is induced by the anomaly and was derived in [72] within the chiral kinetic theory. Quite naturally, the transport coefficient associated with the last term in (1.18) calculated in [72] is different from our result. Alternatively, using the identity (1.18) can be split into two pieces: one represents E 2 -correction to the chiral magnetic conductivity; the other contributes to the chiral electric effect [73]. Analogous to (1.18), the last term in (1.19) is of interest too. On the one hand, it makes E 2 -correction to the chiral separation effect; on the other hand, it contributes to the chiral separation effect induced by the electric field. Our study implies that, via higher order corrections, the chiral electric separation effect exists even when there is no axial chemical potential.
That is, the axial charge density will linearly grow with time leading to instability. It will manifest itself in a breakdown of the constitutive relations (1.18) and (1.19), which would have to be amended by derivative terms.
Beyond the constant background field approximation, numerous new structures emerge, which involve inhomogeneous fluctuations of the external fields and charge densities. Particularly, we notice the emergence of δρ 5 B (δρ B) in J ( J 5 ), see (5.38,5.39). The interplay between these two terms predicts the chiral magnetic wave [68]. We list the new structures in section 5.2, but leave computation of corresponding transport coefficients to future work. This paper is structured as follows. In section 2 we present the holographic model. In section 3 we outline the strategy of deriving the boundary currents through solving the anomalous Maxwell equations in the bulk. Section 4 presents the first part of our study. In section 4.1, we derive the constitutive relations (1.6,1.7) by solving the dynamical components of the bulk anomalous Maxwell equations near the conformal boundary. The boundary external fields and charge densities appear as source terms in these equations. The main technique is based on decomposition of the bulk gauge fields in terms of basis constructed from the external fields and charge densities. Dynamics of the Maxwell equations is translated into ordinary differential equations (ODEs) for the decomposition coefficients, which are nothing else but the components of the inverse Green function matrix. In section 4.2 we determine the TCFs by solving these ODEs. In section 5 we turn to the second part of this work, corresponding to the scheme (1.5). Finally in section 6 we again outline the main results of this work and make some discussions. Some technical details are deposited into several appendices. In appendix A, we write down all the ODEs for the decomposition coefficients and derive constraint relations among them. In appendix B we summarise analytic perturbative solutions for these coefficients. In appendix C we prove that the TCFs in (1.6,1.7) are frame-independent.

The holographic model: from
The holographic model is the U (1) L × U (1) R theory in the Schwarzschild-AdS 5 black brane spacetime. The triangle anomaly of the boundary field theory is introduced via the Chern-Simons terms (of opposite signs for left/right fields) in the bulk action where the Lagrangian density L 1 is (2.2) In the ingoing Eddington-Finkelstein coordinates, the spacetime metric is where f (r) = 1 − 1/r 4 , so that the Hawking temperature (identified as temperature of the boundary theory) is normalised to πT = 1. On the constant r hypersurface Σ, the induced metric γ µν is The bulk theory can be reformulated as U (1) V × U (1) A via the combination where the gauge coupling e (e ) is associated with the vector (axial) field V M (A M ). In terms of V and A fields, the Lagrangian density L 1 becomes which can equivalently be written as where F V,a are field strengths of V, A, respectively. κ and V M are defined as (2.9) While the underlined total derivative term in (2.8) does not affect the equations of motion, it results in non-conservation of the vector current (dual to V M ). Following [24] we are to cancel this total derivative by adding the Bardeen counter-term so that the vector current becomes conserved, as in real electromagnetic theory. In terms of V and A fields, the counter-term action (2.5) is From now on we will work with a new action S (2.12) Equations of motion for V and A fields are derived via standard variational procedure. Under the variation 14) and the variation of where ∇ µ is covariant derivative compatible with the induced metric γ µν . Then, equations of motion for V and A fields are Imposing the dynamical equations (2.16), the action variation δS reduces to where n M is the outpointing unit normal vector of the slice Σ. The last line of (2.20) vanishes either using the constraint equations (2.17) or through a radial gauge choice δV r = δA r = 0.
The boundary currents are defined as In terms of the bulk fields, the boundary currents are It is important to stress that the currents in (2.21) are defined independently of the constraint equations (2.17). Throughout this work, the radial gauge V r = A r = 0 will be assumed. Thus, in order to completely determine the boundary currents (2.22) it is sufficient to solve the dynamical equations (2.16) for the bulk gauge fields V µ , A µ only, leaving the constraints aside. The constraint equations (2.17) give rise to the continuity equations (1.2). In this way, the currents' constitutive relations to be derived below are offshell. In subsequent presentation, the couplings e and e will be absorbed into redefinition of V and A fields, while the notations for V and A will remain unchanged for convenience.
It is useful to reexpress the currents (2.22) in terms of coefficients of near boundary asymptotic expansion of the bulk gauge fields. Near r = ∞,

23)
The holographic dictionary implies that V µ , A µ are gauge potentials of the external fields E, B, E a and B a , (2.27) As mentioned above, only the dynamical equations (2.16) were utilized in order to get (2.23-2.26). The near-boundary data V (2.28) Note explicit dependence of the currents J µ and J µ 5 on the axial gauge potential A µ . The last term in J µ of (2.28) is crucial in guaranteeing conservation of J µ , that is the gauge invariance under the vector gauge transformation V µ → V µ + ∂ µ φ. Clearly, explicit dependence of physical quantities on the axial potential A µ is because that the transformation In presence of anomaly one distinguishes between consistent current and covariant current [76]. Consistent current is defined as a functional derivative of effective action with respect to external gauge field. Covariant current is obtained by subtracting a suitably chosen Chern-Simons current from the consistent one, so that the current becomes invariant under both vector and axial gauge transformations. The currents defined in (2.21,2.28) are consistent. The associated covariant currents are Obviously, when the axial field A µ = 0, both consistent and covariant currents coincide.

Anomalous Maxwell equations in the bulk
To derive constitutive relations for the currents J µ and J µ 5 , we consider finite vector/axial charge densities exposed to external vector and axial electromagnetic fields. Holographically, the charge densities and external fields are encoded in asymptotic behaviors of the bulk gauge fields. In the bulk, we will solve the dynamical equations (2.16) assuming some charge densities and external fields, but without specifying them explicitly. In this section we outline the strategy for deriving of currents' constitutive relations.
Following [54], we start with the most general static and homogeneous profiles for the bulk gauge fields which solve the dynamical equations (2.16), where V µ , A µ , ρ, ρ 5 are all constants for the moment. Regularity requirement at r = 1 fixes one integration constant for each V i and A i . Through (2.28), the boundary currents are Hence, ρ and ρ 5 are identified as the vector/axial charge densities. Next, following the idea of fluid/gravity correspondence [43], we promote V µ , A µ , ρ, ρ 5 into arbitrary functions of the boundary coordinates (3.1) ceases to be a solution of the dynamical equations (2.16). To have them satisfied, suitable corrections in V µ and A µ have to be introduced: where V µ , A µ will be determined from solving (2.16). Appropriate boundary conditions have to be specified. First, V µ and A µ have to be regular over the whole integration interval of r, from one to infinity. Second, at the conformal boundary r = ∞, we require which amounts to fixing external gauge potentials to be V µ and A µ . Additional integration constants will be fixed by a frame choice. In this work we adopt the Landau frame convention for covariant currents, The Landau frame choice can be identified as a residual gauge fixing for the bulk fields. Most of our results, however, would be independent of this choice. Appendix C is entirely devoted to this discussion. The vector/axial chemical potentials are defined as For the homogeneous case, the definition (3.7) results in µ = ρ/2, µ 5 = ρ 5 /2. Beyond the homogeneous case, µ, µ 5 are nonlinear functions of densities and external fields. For generic configurations of external fields and charge densities, (2.16,2.17) become rather involved. In terms of V µ and A µ , the dynamical equations (2.16) are Triangle anomaly is a source of nonlinearity in all these equations. In the context of fluid/gravity correspondence [43], external fields V µ , A µ and charge densities ρ, ρ 5 are assumed to vary slowly from point to point. Consequently, the corrections V µ and A µ can be constructed through order by order expansion in derivatives of the external fields and charge densities. Contrary to the approach adopted below, the method of [43] is implemented using "on-shell" relations. That is, the bulk solutions are constructed with the help of the constraint equations.
To extract the TCFs to all order in derivative expansion, we do linearisation in external fields and charge densities. As announced in section 1, we will solve (3.8-3.11) under two different linearization schemes (1.4) and (1.5).

Study I: linear transport
In this section we study linear TCFs corresponding to the linearisation scheme (1.4), whereρ,ρ 5 are constants. All calculations below are accurate to linear order in . Obviously, The presentation is split into two subsections: one is devoted to derivation of the constitutive relations (1.6,1.7) while the other one focuses on determination of transport coefficients.

Derivation of constitutive relations from the dynamical equations
Under the scheme (4.1), the dynamical equations (3.8-3.11) are (4.6) At linear level, V µ and A µ are still coupled together through the anomaly-induced terms.
To order O( ), the constraint equations (2.17) are which do not feel the effect of triangle anomaly at this order in . The corrections V µ and A µ are decomposed as where S i , V i ,S i ,V i are elements of the inverse Green function matrix. They are scalar functionals of the boundary derivative operators and functions of radial coordinate r. In momentum space, the derivative operators turn into scalar functions of frequency ω and momentum squared q 2 : which satisfy partially decoupled ODEs (A.1-A. 22). Dynamics of the bulk theory is now reflected by these ODEs. Accordingly, the boundary conditions for the decomposition coefficients in (4.9,4.10) are Additional integration constants will be fixed by the frame convention (3.6). Pre-asymptotic expansions (2.23-2.26) translate into pre-asymptotic behaviour of the decomposition coefficients in (4.9,4.10). Near r = ∞, The frame convention (3.6) leads to the following relations s 4 = s 6 =s 1 =s 3 = 0.
Following the definition (3.7), the chemical potentials µ, µ 5 are These relations can be used to replace the charge densities ρ, ρ 5 in (1.6,1.7) in favour of µ, µ 5 . The results are presented in (1.8,1.9) with the coefficients given by where S 2 , S 3 are to be determined in section 4.2.
Putting the currents on-shell, the constitutive relations (1.6,1.7) bear a standard form of linear response theory, from which current-current correlators read While the dispersion relation iω − q 2 D(ω, q 2 ) = 0 is not affected by the anomaly, residues of the correlators (4.26,4.27) get modified. Kubo formulas are usually used to relate transport coefficients to thermal correlators. However, evaluated on-shell, the correlators partially lose information about dynamics of off-shell one-point currents. As a consequence, they are insufficient to determine all order transport coefficients. For example, beyond their constant values, D, σ e/m cannot be fully extracted from the current-current correlators [54]. Yet, there are exact relations between the correlators and σ χ , σ κ , which are valid for arbitrary ω and q. The relation (4.28) for σ χ was first derived in [10] by promoting constant magnetic field in the original CME into an inhomogeneous perturbation. Our constitutive relation translates into rigorous derivation of (4.28). Contrary to σ m , σ a can be determined from the correlators, particularly from the mixed correlator J i J j 5 . This became possible thanks to the absence of the E a ( E) term in the constitutive relations for J µ (J µ 5 ). We suspect this is accidental and specific to the model in study.

Results: solving the bulk equations
To determine all the TCFs in (1.6,1.7), we merely need to solve the following ODEs: We first solve the ODEs analytically in the hydrodynamic limit and then numerically for arbitrary ω and q.

Hydrodynamic expansion: analytical results
In the hydrodynamic limit ω, q 1, the ODEs (A.3, A.5-A.10) can be solved perturbatively. Let introduce a formal expansion parameter λ ω → λω, q → λ q. (4.30) At each order in λ, the solutions are expressed as double integrals over r, see appendix B. Below we collect the series expansions of The series expansions of v 2 , v 4 were worked out in [54], Once substituted into (4.18), these perturbative results generate the hydrodynamic expansion of all the TCFs as quoted in (1.10-1.15).

Beyond the hydro limit: numerical results
To proceed with the all-order derivative resummation, we have to go beyond the conventional hydrodynamic limit and solve the ODEs (A.3, A.5-A.10) for generic values of momenta. We were able to do it numerically only. We deal with a boundary value problem for a system of second order linear ODEs. We performed our numerical calculations within two different numerical methods, a shooting technique and a spectral method. Within numerical accuracy, both approaches give the same results. Within the shooting technique our numerical procedure is much like that of [54]. One starts with a trial initial value for the functions to be solved for at the horizon r = 1 and integrates the ODEs up to the conformal boundary r = ∞. The solutions generated in this step have to fulfil the boundary conditions at r = ∞. If not, the trial initial data have to be adjusted and the procedure is repeated until the requirements at the boundary are satisfied with a satisfactory numerical accuracy. The fine-tuning process of finding the correct initial data is reduced to a root-finding routine, which can be implemented by the Newton's method.
The spectral method converts the continuous boundary value problem of linear ODEs into that of discrete linear algebra. We distribute a number of points on the integration domain. These points are collectively referred to as collation grid. The functions to be solved for are then represented by their values on the grid. For given values of the functions on the grid, their derivatives at the grid are approximated by differentiating interpolation functions (normally based on polynomial or trigonometric interpolation). Thus, a differential operation is mapped into a matrix. Eventually, this collation procedure allows to discretise the original continuous problem and turns it into a system of algebraic equations involving values of the functions on the grid. The boundary conditions are mapped into algebraic relations among the values of the functions at the outermost grid points. They will be imposed by replacing suitable equations in above-mentioned algebraic equations. For more details on spectral method, we recommend the references [77][78][79][80]. As for other non-periodic problems, we choose a Chebyshev grid and use polynomial interpolation to calculate differentiation matrices.
Since D and σ e are not affected by the anomaly, they are the same as those presented in [54], and we would not display them here. As σ κ can be obtained from σ χ via the crossing rule as pointed out below (4.18), we will focus on numerical plots for σ m , σ a , σ χ only. Figure 1 is a reproduction of a 3D plot for the magnetic conductivity σ m from [54]: compared to [54], we extend the plot domain to larger momenta so that asymptotic regime is more clearly seen. Figures 2,3 show anomaly-modified σ m for sample values of κμ and κμ 5 . In Figure 4 we show momenta-dependent σ a . Note that σ a is non-vanishing only when κμμ 5 = 0.  In Figures 5,6 we show 2D slices of Figures 1,2,3,4 when either ω = 0 or q = 0. It is demonstrated that asymptotically σ m approaches a nonzero value, while σ a goes to zero after damped oscillation. The asymptotic regime is achieved around ω 5 for both σ m and σ a . So, σ m encodes some UV physics while σ a decouples asymptotically. In contrast with the viscosity function [55,56] and the diffusion function [54], Figure 6 illustrates that the dependence of σ m , σ a on q is more pronounced. This feature is shared by the chiral magnetic conductivity σ χ (see Figure 9). As an axial analogue of σ m , σ a shows qualitatively similar dependence on ω as can be seen from Figure 5. However, q dependence of σ a differs from that of σ m as shown in Figure 6. When κμ and/or κμ 5 get increased, ω-dependence of σ m becomes more enhanced, which signifies a stronger response to time-dependent external fields.     We now turn to the chiral magnetic conductivity σ χ . For weakly coupled theories, momenta-dependence of σ χ was studied in [10,67]: in high temperature regimes Re(σ χ ) drops from its DC limit σ 0 χ at ω = 0 to σ 0 χ /3 just away from ω = 0. For strongly coupled theories with dual gravity description, frequency-dependence of σ χ was initially considered in [24] for the case q = 0. Two different holographic models were considered in [24]: RN-AdS 5 geometry and a finite-temperature Sakai-Sugimoto model [59]. Numerical plots of [24] look rather similar, suggesting a certain universality of σ χ . Within the RN-AdS 5 model, ref. [67] also explored the momenta-dependence of anomalous TCFS, particularly original results on σ κ (ω, q) were presented there.

Study II: nonlinear transport induced by constant external fields
In this section we turn on constant backgrounds for external fields, corresponding to the linearisation scheme (1.5), whereρ andρ 5 are treated as constants.V µ andĀ µ depend linearly on x α so that their field strengthsF V µν andF a µν are constant backgrounds E, B, E a , B a . The corrections V µ , A µ of (3.4) are expanded to linear order in .

Solutions for V
So, (5.3-5.6) can be rewritten in integral forms, where using the definition (3.7), the chemical potentials µ, µ 5 are From (2.28), the boundary currents are which reduce to (1.16,1.17) whenĀ µ = 0. Despite the fact that the equations for V i are linear, the solutions involve complex inverse propagators, which are non-linear functions of the background fields. To proceed, we resort to a weak field approximation by introducing yet another formal expansion µ are formally expanded in powers of α, We analytically solved (5.3-5.6) up to order O( 0 α 2 ). The results are summarised below.

Solutions for V
(1) The corrections V At the lowest order O 1 α 0 , the dynamical equations (3.8-3.11) are exactly (4.3-4.6) with V µ , A µ in (4.3-4.6) replaced by δV µ , δA µ . Therefore, solutions for V where S 3 , V 2 , V 4 were studied in [54] while S 2 , V 3 , V 5 ,V 3 ,V 5 were adressed in section 4.2. As a result, at order O 1 α 0 , the boundary currents J µ , J µ 5 are As α → 0, the above results coincide with those of section 4. To the order O 1 α 1 , the dynamical equations (3.8-3.11) are , which then get propagated into the constitutive relations for the currents J µ and J µ 5 . The full solutions can be constructed in parallel with section 4. We have decided not to solve (5.34-5.37) in this publication and to leave a comprehensive study of V (1)(1) µ and A (1)(1) µ and corresponding transport coefficients for future work. Here, we merely list all the basic structures that would emerge in the constitutive relations, as dictated by the source terms in (5.34-5.37). In J we anticipate to have the following terms each multiplied by its own TCF, (5.38) In the axial current J 5 , (5.39) The terms δρ 5 B, δρ B in (5.38,5.39) would lead to the chiral magnetic wave [68], which reflects density fluctuations δρ, δρ 5 at constant external magnetic field. According to [68], the speed of the chiral magnetic wave can depend on B nonlinearly, the property which is inherited from a nonlinear B-dependence of µ, µ 5 . In [68], this nonlinearity of µ, µ 5 was realised in a top-down holographic QCD model based on a DBI action for the bulk gauge fields. In contrast, working with the canonical Maxwell action, our study demonstrates that similar non-linear phenomena can emerge solely from the Chen-Simons term.
The terms ∇δρ 5 × E and ∇δρ × E in (5.38,5.39) were studied in [81] within the chiral kinetic theory. It would be interesting to compare the results once the corresponding transport coefficients are computed in the holographic model.

Conclusions
In this paper, we have revised anomaly induced transport in a holographic model containing two U (1) fields interacting via Chern-Simons terms. For a finite temperature system, we have computed off-shell constitutive relations for the vector and axial currents responding to external vector and axial electromagnetic fields.
Within the linear approximation, the TCFs D(ω, q 2 ) and σ e (ω, q 2 ) are left unaffected by the anomaly and were computed previously in [54]. While σ m (ω, q 2 ) gets an anomaly induced correction, the remaining TCFs, σ χ/κ/a , are induced by the anomaly. In the hydrodynamic regime, we have analytically reproduced all the known results in the literature and succeeded to extended the gradient expansion to third order, see (1.10-1.15). Beyond the hydrodynamic regime, these transport coefficient functions were numerically calculated up to large values of momenta so that the asymptotic regime is reached. The results are displayed by the plots in section 4.2. The electric/magnetic conductivities σ e/m are the only TCFs that survive at asymptotically large ω (the asymptotics is reached around ω 5).
Nonlinear transport has been studied in a specific setting of the external fields having constant backgrounds. When the only non-vanishing external field is a constant magnetic field, the CME has been shown to be exact, relating the induced vector current and the magnetic field, see (1.16). This exact relation is nonlinear in the magnetic field and the entire nonlinearity is absorbed into the axial chemical potential µ 5 , see (1.20). Electric fields lead to new nonlinear effects, see (1.18). Small time-dependent/non-homogenous perturbations introduce many more interesting anomaly-induced structures in the constitutive relations. We have merely listed these structures, leaving determination of associated new transport coefficients for a future study.
Additional non-linear anomaly-induced effects are explored in our forthcoming publication [60]. Particularly, CME in presence of a space-varying magnetic field B( x) is found to be modified by derivative corrections. An interplay between constant magnetic and time-dependent electric fields is another focus of [60].
At the lowest order in derivative expansion, the anomalous transport coefficients σ 0 χ/κ are known to be dissipationless [42,71]. Particularly, they do not contribute to entropy production in a hydrodynamic system. It would be interesting to classify all the higher derivatives/nonlinear terms in the constitutive relations in accord with their dissipative nature [82]. We do expect that the higher order gradient terms would introduce dissipation and this would potentially affect various phenomena such as the chiral drag force [83,84].
Our study has been carried out in the probe limit, in which the currents get decoupled from the dynamics of the energy-momentum tensor. Beyond this limit, new phenomena emerge such as the normal Hall and chiral vortical effects [62,71]. Interplay between the vorticity and strong magnetic field [85] is an interesting direction worth further study in a holographic setup beyond the probe limit.
A ODEs and the constraints for the decomposition coefficients in (4.9,4.10) In this appendix, we collect the ODEs satisfied by the decomposition coefficients in (4.9,4.10), and derive the constraint relations relating these coefficients. Substituting (4.9,4.10) into (4.3-4.6) and making Fourier transform ∂ µ → (−iω, i q), we arrive at the ODEs, which are grouped into several partially decoupled sub-sectors. (A.10) Furthermore, This comes from the fact that these functions satisfy identical ODEs (A.1,A.2) with identical boundary conditions. Additionally, as in [54] consider the combinations which satisfy homogeneous equations. Therefore, under (4. 16,4.17) and boundary conditions (4.11,4.12) we have The relations (A.23, A.24, A.27, A.28) could be further translated into constraints among The relations (A.23, A.24, A.27, A.28) reduce the number of equations that needs to be solved. The remaining independent sub-sectors are Note that under the interchange The ODEs satisfied by the sub-sectors V 1 ,V 1 , V 5 ,V 5 and V 6 ,V 6 , V 10 ,V 10 get exchanged in the following way, Given that V 1 ,V 1 , V 5 ,V 5 and V 6 ,V 6 , V 10 ,V 10 obey the same boundary conditions (4.11,4.12), the following "symmetric" relations hold (A.32) Based on (A.27,A.28), (A.32) implies Eventually, we only need to solve and obtain all the other functions through the relations revealed above.
Expanding these functions in powers of λ, we solve (A.3,A.5-A.10) perturbatively, order by order in λ. Final results are quoted below.