Quartic Horndeski, planar black holes, holographic aspects and universal bounds

In this work, we consider a specific shift-invariant quartic Horndeski model where the function $G_4$ is chosen to be proportional to $\alpha \sqrt{(\partial \psi)^2/2}$, $\alpha$ denoting the coupling of the theory and $\psi$ being some massless scalar field. We take $(D-2)$-many copies of this model, and we assume a linear dependence of the scalar fields on the coordinates of the $(D-2)$-dimensional Euclidean submanifold, $\psi^{I}=p\delta^{I}_i x^{i}$. These choices allow us to construct planar black holes with a non-trivial axion profile which we explore in terms of their horizon structure and their thermodynamic properties. Since the particular scalar field profile dissipates momentum in the boundary theory, we are able to derive a sensible DC transport matrix describing the linear thermoelectric response of the holographic dual to an external electric field and a thermal gradient. We comment on the form of the conductivities and show that the heat conductivity-to-temperature ratio cannot have a universal lower bound at all scales due to the new coupling. We also verify the Kelvin formula motivated by the presence of an $\mathrm{AdS}_2\times\mathbb R^2$ horizon. Using the constants $D_{\pm}$ describing coupled diffusion at finite chemical potential, we show that these decouple in the incoherent limit of fast relaxation flowing to the usual expressions of charge/energy diffusitivities respectively, with $\alpha$ playing no qualitative role at all in this process. The bound of the refined $TD_{\pm}/v_B^2$ ratio is investigated, $v_B$ being the butterfly velocity. The coupling enters the latter only through the horizon radius, and it does not affect the $\mathcal{O}(1)$ form of the ratio in the strong dissipation regime. Finally, the viscosity-to-entropy ratio is computed by means of a (weaker) horizon formula, and the simple $(4\pi)^{-1}$-bound is found to be violated.


Introduction
Scalar-tensor theories of gravity have been well studied in the past, their pros and cons, as alternative theories of gravitation, elucidated in detail. They enrich the dynamical field content of General Relativity by the inclusion of scalar fields in the latter which constitute additional degrees of freedom. One of the most renowned scalar-tensor family is Horndeski gravity, the most general four-dimensional scalar-tensor theory with equations of motion containing up to second order derivatives of the dynamical fields. The field content of Horndeski gravity consists of the spacetime metric g µν and a scalar field φ. We focus on the subclass where the latter enjoys a global shift symmetry under which φ → φ + c, with c being some constant. In this scenario, the Horndeski action involves four arbitrary functions of the canonical kinetic term X := −(∂φ) 2 /2, denoted by G 2 , G 3 , G 4 and G 5 . It reads where Here, the following shorthand notation is used: φ µν... := ∇ µ ∇ ν ...φ and (φ µν ) p := φ λ 1 λ 2 φ λ 2 λ 3 ...φ λp λ 1 .
Additionally, R is the Ricci scalar, and G µν is the Einstein tensor. In the context of covariant Galileon theory, these non-minimal couplings play the role of counterterms which cancel out with the higher derivative terms arising from the variation of the action [2][3][4].
As said above, the justification of the shift symmetry enjoyed by the scalar field, lies in the restriction G i (φ, X) → G i (X) which assures that the first derivatives of the scalar field account for the lowest derivative order being present in (1.1). Due to this symmetry, it becomes possible to write the equation of motion for φ, i.e. the Klein-Gordon equation, in terms of the Noether current associated with global shift symmetry, Usually, when introducing additional degrees of freedom, one tries to see whether the latter allow for solutions with non-trivial profiles. No-hair theorems such as [5,6] state assumptions which, if met, forbid the departure from the GR solution spectrum. Perhaps more relevant to our case, are the no-hair arguments of [7] applying to static and spherically symmetric spacetimes in the framework of Horndeski gravity. By using the SO(3)-symmetry of the ansatz and the time-reversal invariance of the action, the authors of [7] show that the only non-vanishing component of the aforementioned current is the radial one. Then, assuming (i) asymptotic flatness and (ii) regularity of diffeo-invariant quantities at the horizon, like J µ J µ , together with (iii) vanishing boundary conditions at infinity -that is, taking φ ′ to vanish there-, it is finally proven that J r (r 0 ) = 0 where J r assumes the form J r = g rr φ ′ H(φ ′ , g, g ′ , g ′′ ). (1.5) This leads to J r = 0 everywhere via the conservation law, followed by the conclusion φ ′ = 0 at all radii, provided H asymptotes to a nonzero constant when φ ′ → 0. The assumptions have to be supplemented with (iv) functions G i such that their derivatives with respect to X do not introduce negative powers of X as the latter approaches the origin, and (v) the canonical kinetic term X must be present in the action. Then, the theorem guarantees that static, spherically symmetric black holes with non-trivial scalar field profiles cannot exist. One controversial beauty of no-hair theorems revolves around possible ways of circumventing their prohibitive results. Indeed, it has been shown that relaxing some of the hypotheses of [7] allows for non-trivial scalar hair. Giving up on (i), several (A)dS or Lifshitz black holes were reported [8][9][10][11][12][13]. In these cases, although H asymptotes to 0, the scalar field profile is non-trivial nevertheless. Another circumventing route goes through allowing the scalar field to linearly depend on time [14], providing several stealth solutions with asymptotic behaviors depending on whether (v) is broken or not. These solutions have also provided a natural scenario for the construction of neutron stars which avoid conflicts with Solar System tests [34,35]. Furthermore, failing (iv), static hairy solutions with asymptotic flatness can be constructed [15]. The lesson to be learned here is that by relaxing the hypotheses, one either strengthens no-hair theorems or avoids their confining results. This also constitutes part of the motivation behind this work. Mobilized by [14], instead of allowing for a linear time dependence, we will first start by introducing (D − 2)many scalar fields ψ I with I = 1, ..., D − 2 being some internal index labeling them, D being the spacetime dimension. The dynamics of these scalar fields will be dictated by the variational equations of motion derived from (D − 2)-many copies of the quartic Horndeski action with G 3 = G 5 = 0, G 2 ∝ X and G 4 ∝ √ −X. We will have the scalar fields reside in a Euclidean target space R D−2 , the latter being isomorphic to the spatial piece of the conformal boundary. This target space can be identified with the planar base sub-manifold of Euclidean signature and dimensions D − 2. We will focus on solutions breaking translation symmetry in all planar directions, i.e. massless Stückelberg fields with a linear bulk profile along the planar base sub-manifold directions, ψ I = pδ I i x i . They can be equally understood as magnetically charged 0-forms with their charge being proportional to the slope p of the profile; in this sense, they should classify as primary hair [16]. Having broken some of the assumptions of [7], we will work our way towards a new family of electrically charged hairy planar 1 black holes (PBH) characterized by a nonzero axion background. In general, the idea of looking for black hole solutions with various types of scalar fields or k-form fields homogeneously distributed along the planar directions, has been a fruitful practice as evidenced by some articles, e.g. see [17][18][19].
Additionally, the inclusion of such fields promotes another agenda: it leads to momentum dissipation in the dual field theory resolving the delta function multiplying the Drude weight in the conductivity formula; Dissipation happens because the bulk scalars -and hence the zeroth order terms in their Fefferman-Graham expansion which source the scalar operators of the dual field theory-have a spatial dependence. Simple boundary U 1 -and diffeomorphism-symmetry arguments suffice to derive the diffeomorphism Ward identity for the (non-)conservation of the boundary stress energy tensor, where one eventually sees that T ti fails to be conserved, exactly due to the spatial dependence of the ψ I 's. The presence of this relaxation mechanism will allow us to pursue the second objective of this work which is to apply holographic techniques in order to compute the DC transport coefficients [20][21][22] of the holographic dual in the broader gauge/gravity duality context of the renown AdS/CFT correspondence [25,26]. Many studies in GR and alternative theories of gravitation have been carried out, and their results, from a holographic point of view, have been elucidated [16,[27][28][29][30][31][32][33]36]. Furthermore, gauge/gravity duality has been also providing towards the study of fascinating phenomena in strongly correlated systems, indicative examples being the linear T -resistivity 2 and the universal Homes's law [37][38][39]64]. Insight from holography has also been given into various bounds and their possible universality, examples being the thermoelectric conductivity bounds [52,53], the Kovtun-Son-Starinets (KSS) bound [42,49] of the viscosity-to-entropy ratio, the universal bounds for the charge/energy diffusion constants [63,64,67] in the regime where diffusive physics dominate e.t.c. Especially interesting was the refinement [63] of the original Hartnoll proposal [64] -the latter supported by experimental data on dirty metals as well-by identifying the characteristic velocity of the system with the so called butterfly velocity v B [65,66,76], which measures how fast quantum information scrambles and proves to be a good candidate at strong coupling. As the refinement is an outcome of holographic methods, it always makes sense to probe it in various holographic models. In general, certain deformations of the bulk action can affect the universality of some of these bounds; Higher derivative corrections [40] can drive the charge diffusitivity bound all the way to zero, the inclusion of the Gauss-Bonnet term lowers the viscosity/entropy bound to a smaller O(1) number [41], e.t.c. Motivated by these works, and by the holographic treatment of scalar-tensor theories, e.g. [21,46,51], we test this particular quartic Horndeski model against the aforementioned conjectured bounds.
This paper is organized as follows: In Sec. 2 we formulate the action principle extracting its variational field equations. We derive electrically charged PBHs with axion hair in arbitrary spacetime dimensions D > 3, followed by a discussion about the horizon structure of the four-dimensional solution. We close the section by presenting a straightforward extension for a dyonic gauge field, while we also display the three-dimensional case which cannot be derived by simply taking the limit. Sec. 3 is dedicated to the study of the thermodynamic properties of the planar solution, where we derive the entropy via both the Wald formalism and the Euclidean path integral approach. Both methods agree on the result, and the first law is shown to hold true as well. The 1/4-area law for the entropy is modified, although it can be said to hold in the broader sense, in units of an effective gravitational coupling [43]. We move on to Sec. 4 where we probe holographic features of the bulk theory starting with DC transport coefficients. After we derive the thermoelectric response, we compute the Lorentz ratios, verifying the Kelvin formula and proceeding with a discussion about the bounding of heat conductivity κ. Due to the new coupling α, we find that there exists no fixed O(1) number bounding the κ/T ratio at all scales from below; the minimum value of the ratio rather depends on the relative strength of dissipation. Moreover, since the model possesses an AdS 2 × R 2 horizon, we check the Kelvin formula concluding that it holds. Next, we show that in the strong dissipation regime, regardless of α, the diffusion constants D ± decouple as the mixing term becomes negligible, and one can use the simpler formulas for the charge and energy diffusion, D c and D e , respectively. After determining the butterfly velocity v B , it is shown that α does not alter the physics in the incoherent limit, where the charge/energy diffusitivities-over-v 2 B ratios are bounded by the standard numbers from below. Neither does the new coupling play any role in the low temperature expansion. In general, it does not have any leading order contribution in these cases. Finally, we end the section with a very shallow investigation of the viscosity/entropy ratio where an explicit violation of the simple KSS bound is manifest, also due to the non-minimal coupling of the theory. In Sec. 5 we conclude.

The model: action, field equations and hairy solutions
Let us start by introducing the model we will work on. Let X I = −(∂ψ I ) 2 /2 be the canonical kinetic term of -what for the moment is-the I-th scalar field 3 . Then, the shorthand notation G I := G(X I ) will be convenient. Also, when writting ψ I µν... we mean ∇ µ ∇ ν ...ψ I . Having said that, we restrict to the quartic sector of (1.1) with meeting the no-hair assumption (v), i.e. G 3 = G 5 = 0 with G 2 =ηX − 2(16πG N ) −1 Λ/(D − 2),η ∈ R + . In order to fail (iv), we choose G 4 = α √ −X + (16πG N ) −1 /(D − 2). Then, where Ψ I := ( ψ I ) 2 − (ψ I µν ) 2 , and where we will set 16πG N to unity together with the AdS radius such that Λ = −(D − 1)(D − 2)/2. We consider the action functional with F 2 := F µν F µν and F µν = 2∂ [µ A ν] , A µ being the U 1 -field. The existence of a boundary action, necessary for a well-defined variational principle, is implied. In a more clear form, (2.2) can be rewritten as Observe that S ψ I contains the non-minimal coupling as well. Stationary variations with respect to g µν , A µ and ψ I yield Here, As mentioned in the introduction, we will focus on a particular class of solutions to the Klein-Gordon equation, namely ψ I ≡ ψ i = px i ; this breaks translation invariance in the planar directions but retains the little SO D−2 -symmetry. Such a solution will also simplify the calculations significantly.

Electrically charged PBH with a non-trivial axion profile
Let us consider a static spherically symmetric metric together with the bulk profile ψ i = px i , and let us start solving equations. The easiest equation to start with is the Maxwell one. We consider a purely electric field strength tensor F µν = −A ′ (r)δ tr µν for a Maxwell potential 1-form A = A(r)dt. Then, we need to This is solved by an electric field Let us proceed with the Klein-Gordon. First of all, X i = −p 2 /(2r 2 ). The only nonvanishing component of J Iµ is the i-th one, where we used the fact that It is evident that the Klein-Gordon is identically solved, because J Ix i is a function of the radius only. Moving on to the metric field equations, these are solved by Observe that as we approach α = 0, the limit where the non-minimal coupling vanishes and the theory reduces to Einstein-Maxwell-(linear) Axion (EMliA), the solution expands as i.e. the axionic solution in [20] is recovered as expected. For asymptotically flat spacetimes, the axion currents (2.15) are regular everywhere on the shell, except, of course, at the origin where a true curvature singularity resides. This can be seen by simply having a look at their asymptotic behavior; they go as O(1/r). The horizons are located at the real positive roots of Only for flat asymptotics, a change of variables z = r D−3 reveals a quadratic equation with roots , (2.20) in terms of the radial coordinate. Only r + is real in the physical domain of the solution, where α, p,η real with αp > 0,η > 0 (no ghost excitations) and D > 3. Hence, a planar asymptotically flat black hole is allowed to form, with its single horizon located at r 0 = r + . However, when the solution asymptotes AdS, a different causal structure emerges, as it will be shown for the upcoming four-dimensional case.

The four-dimensional solution
The four-dimensional solution is given by the D → 4 limit of (2.17). It reads and admits the asymptotic expansion where the effective mass and the effective charge read the former chosen to be positive. First, this is the most general solution; had we started with g tt = −U (r) in the ansatz, we would see that the system of field equations would force U to be a multiple of F times an integration constant which can always be fixed such that U = F . Second, we observe that if it wasn't for the presence of the r-term, our solution, taken at asymptotic infinity, would effectively behave as the solution in [20], without the need to kill the G 4 -copies. It is then easy to check that for flat asymptotics, the metric function does indeed have the last mentioned behavior, with the non-minimal coupling switched on as well, since the O(r) terms vanish and F approaches infinity as O(1).
The horizons of (2.21) are located at the real positive roots of the depressed quartic equation The multiplicity and reality of the roots depend on the sign of the discriminant. It is best if we study the extrema of the auxiliary function instead. These are located at the positive real solutions of We find that there exists only one positive real solution which corresponds to the global minimum which is independent of the new parameter α. For M < M * the singularity is naked, whereas when M > M * the black hole possesses two horizons, the outer one located at the largest root of the left hand side of (2.24). We remark that their location can be analytically determined since (2.24) is analytically solvable, but the explicit expressions are too lengthy to write down. When the inequality is saturated, an extremal black hole forms with its horizon located at r 0 = r * , the latter also being written as in terms of what will later be identified with the chemical potential of the holographic dual, given the expression A = µ(1 − r 0 /r). The extremal black hole asymptotes a unit-radius AdS 4 , while near the horizon, an AdS 2 × R 2 product structure appears with which of course agrees with the findings in [20] when the matter sector couples minimally via the kinetic term only. Finally, it is worth mentioning that the dyonic extension of (2.21), for a Maxwell potential 1-form A = A(r)dt + Q m x [1 dx 2] , is straightforward. The solution representing a dyonic PBH simply reads where the expected interchange duality Q e ↔ Q m is apparent. We would like to close this subsection by presenting the three-dimensional solution as well. The latter corresponds to a logarithmic branch and it cannot follow from the D → 3 limit of (2.17). In this separate case, and it is evident that this is simply a charged BTZ solution with the axion flux playing the role of the magnetic charge.

Black Hole thermodynamics
In this section, we will focus on the thermodynamic properties of the black hole solutions derived in Sec. 2.2. Indeed, even if (2.21) has the standard AdS asymptotic behavior, it nevertheless remains interesting to investigate it in terms of black hole thermodynamics. This study is further motivated by the presence of an unusual coupling between the scalar curvature and the square root of the kinetic term. Such a coupling is expected to modify the 1/4-area law of the entropy as we will see below. On the other hand, as it was pointed out in [23], the presence of a non-minimal coupling parameter generates some obscure facets when analyzing the thermodynamic properties of static Horndeski black holes. Indeed, in the last reference, an asymptotically AdS static black hole solution of a particular G 2 -and G 4 -Horndeski Lagrangian [9] was scrutinized from a thermodynamical point of view. It was observed that the Wald formalism [23], the regularized Euclidean method [9] and the quasilocal approach [45], all applied to this specific solution, give rise to distinct expressions of the thermodynamic quantities. This is somehow intriguing since these different approaches are usually consistent with each other 4 . It is clear that these discrepancies are essentially due to the non-minimal coupling between the geometry and the derivatives of the scalar field, but also due to the fact that the static scalar field and its radial derivative diverge at the horizon. As a first step, we will compute the so-called Wald entropy, S W , defined by for the solution (2.21). Here, the integral is taken over a slice of the horizon and L is the full Lagrangian. Also, ǫ µν denotes the unit bivector, normal to the horizon surface, whileǭ stands for the area of the slice. Finally,σ denotes the volume of the planar base sub-manifold. As previously anticipated, the non-minimal coupling between the scalar curvature and the square root of the kinetic term does indeed modify the standard 1/4area law of the entropy; if S o is the standard entropy 4σπr 2 0 in 16πG N = 1 units, then Wald's Noether charge entropy formula simply spits out In the sense of [43], the 1/4-area law still holds in 16πG EFF = 1 units where the effective running gravitational coupling takes the particular value an expression that could have been equally guessed from (2.2) as well. We will also confirm (3.1) by means of the Euclidean approach for which the thermal partition function is identified with the Euclidean path integral at the saddle point around the classical solution.
To do so, we consider the following Euclidean ansatz: together with an electric ansatz A µ dx µ = A(r)dτ . Here, τ is the Euclidean periodic time with period 0 ≤ τ < T −1 =: β where T stands for the Hawking temperature. The range of the radial coordinate r is given by r ≥ r 0 . In the mini superspace of the symmetry ansatz (3.4), the Euclidean action, I E , reads where Π = r 2 F rt corresponds to the canonical radial momentum conjugate to the gauge field. Moreover, B E is an appropriate boundary term codifying all the thermodynamic properties, while also ensuring that the solution corresponds to an extremum of the action, i.e. δI E = 0. Note that the volume element of the Euclidean action is not only radial, as it usually is; this is due to the fact that the scalar fields are assumed to depend on the planar coordinates x and y. A simple exercise shows that the Euler-Lagrange equations obtained from variation of the symmetry reduced action (3.5) with respect to the dynamical fields F, N, A and ψ i lead to the solution (2.21) with N = const., which we can set to unity without loss of generality. We now consider the formalism of the grand canonical ensemble, varying the Euclidean action while keeping fixed the temperature, the electric potential Φ = A(r 0 ) and the constant p controlling the strength of momentum dissipation in the dual theory. Under these considerations, the extremality condition δI E = 0 fixes the variation of the boundary term to read whereσ is the volume of the two-dimensional compact planar sub-manifold. At infinity, (3.6) evaluates to δB E (∞) =σβδM , while at the horizon, the lake of conical singularity -ensured by requiring that δF (r 0 ) = −4πT δr 0 -yields The latter can be re-written as As it is well-known, the Euclidean action is related to the Gibbs free energy G in the following manner: where Φ is the electrostatic potential and Q e the electric charge. Now, from the boundary term (3.8), it is easy to identify the various thermodynamic quantities. We read off that As a first observation, we remark that the expression of the entropy S perfectly coincides with the Wald entropy (3.1). Finally, in order to discuss the validity of the first law, we derive the Hawking temperature which reads It is a matter of fact to check that the first law, namely dM = T dS + A(r 0 )dQ e , holds if the slope p of the axion profile is considered as a fixed parameter. Moreover, the temperature is a monotonically increasing function of the horizon radius which vanishes at r 0 = r * , the latter defined in (2.28) as the horizon radius of the extremal black hole. It will be also useful to invert (3.11) so that

Holographic aspects
PBHs with a non-trivial axion backgrounds along the planar directions provide an ideal configuration for the computation of holographic DC responses. Their presence ensures that translation symmetry is broken which would otherwise lead to nonsensical infinite results. The breaking of the aforementioned symmetry allows for momentum dissipation in the dual field theory which in turn opens the door to finite associated DC conductivities.
In the boundary language, momentum relaxation simply means that ∇ t T ti = 0. To see that this holds, let us present a heuristic argument. Consider the four-dimensional renormalized/regularized version of the bulk action (2.2); dub it S REN . We saw that (2.21) is asymptotically AdS, hence it is conformally compact Einstein [44], and it can be brought to the Fefferman-Graham form the boundary being at ρ = 0. We take a, b, ... to be boundary indices, whereas i, j, ... are used for the spatial piece of the latter. The various fields admit the near boundary expansionsĝ ab =ĝ Since D = 4 there is no term in the series expansion related to the holographic conformal anomaly. For our purposes, determining the coefficients is irrelevant, since we particularly care about the sources. For a euclidean boundary signature, the on-shell variation of S REN yields where now summation is implied for all repeated indices. Symmetry under the boundary Here, ∇ a is the covariant derivative associated with theĝ (0) -compatible connection and It is clear that since ψ i (0) ∼ x i by assumption, T ti will not be conserved for O i = 0. Hence, boundary momentum gets dissipated in the spatial directions, whereas the energy is of course conserved.

Thermoelectric DC transport
It was shown in [22,47] that the electric, thermoelectric and thermal conductivities can be computed in terms of the black hole horizon data alone without the need to invoke direct calculations in the (boundary) field theory. This is achieved by properly manipulating the bulk field equations, revealing electric and heat currents which are manifestly independent of the holographic radial coordinate. These can be then evaluated at the horizon radius instead of the boundary. We start by considering the four-dimensional limit of (2.2), an action functional of the metric g µν , the gauge field A µ and the two axion fields ψ x and ψ y , where we take the bulk coordinates to be x µ = {t, r, x, y}. Studying the gauge field equations in the bulk, we observe that the only non-vanishing component is Defining the current density J t = r 2 F tr , this corresponds to the charge density of the dual field theory when the right hand side is evaluated at the boundary, i.e. Q ≡ J t , where Q is the charge of the black hole, what will be Q e in our case. Moreover, we assume the existence of a regular horizon at r 0 (in the case of two horizons, the outer one is chosen), about which we assume the Taylor expansions F ∼ 4πT (r − r 0 ) + ... and A ∼ A ′ (r 0 )(r−r 0 )+..., namely we take the electric potential to vanish at the horizon radius. We will use Eddington-Finkelstein coordinates (v, r) in order to make the regularity at the horizon apparent with v = t + (4πT ) −1 ln(r − r 0 ). We will also assume the asymptotic behavior A ∼ µ − Qr −1 + ... where µ is the chemical potential in the dual theory, defined as ∞ r 0 drF rt , while the dominant power in the asymptotic expansion of F will be ∼ r 2 .
It is time to proceed with the perturbations. For starters, we will turn on a constant electric field of magnitude E in the x-direction such that supplemented by the small perturbations g tx = ǫh tx (r), g rx = ǫr 2 h rx (r), ψ 1 = px + ǫX (r) (4.9) about the black hole background given by (2.21) and A = µ − Q e /r = µ(1 − r 0 /r). Here, ǫ is introduced as a small parameter helping us keep track of the perturbation order. We will now study the gauge field current density which possesses only one non-trivial component, the one in the x-direction, This can be evaluated at any r, and it is radially conserved since it is derived by integrating the equation ∂ r ( √ −gF xr ) = 0. This means that we are allowed to evaluate it at the horizon radius instead of the boundary. Next, we look at the metric field equations, We observe that G rx = 0 is an equation algebraic in h rx which is solved by where F is always on the background shell since we used the fact that G yy = 0 to arrive at this particular expression. The linearized axion field equations also follow from (4.12). In addition, we also have the second order inhomogeneous ODE: which corresponds to G tx = 0. To move on, we need to impose boundary conditions. We first need to check the gauge field perturbation and its regularity at the horizon. In Eddington-Finkelstein coordinates, the full gauge field perturbation reads (4.14) Taylor-expanding this about r 0 one sees that its regularity is ensured only if It is also evident that near the horizon because 4πT = F ′ (r 0 ) and F ∼ F ′ (r 0 )(r − r 0 ) + .... Now, we can also see that (4.12) diverges as r → r 0 because of the presence of F in the denominator. In order to save this, we let h tx expand as near the horizon. Then, one can immediately see that (4.13) vanishes when evaluated at the horizon. As for the axion field perturbation X we just assume a constant value at r 0 and sufficient falloff at infinity. The remaining boundary conditions at radial infinity are discussed in [22] in detail. Having established well posed perturbations of the bulk fields, we can easily extract the electric DC conductivity, by first evaluating (4.10) at r 0 and at leading order in (r − r 0 ), further dividing by the external electric field of magnitude E, i.e. σ = J x /E at the horizon. We find that This is in perfect agreement with [20] when α = 0. The a la quartic Horndeski non-minimal coupling of the axion fields modifies the electric conductivity which deviates from the results obtained considering EMliA theory. However, and most importantly, the behavior at both ends is the same. As we saw when we studied the horizon structure, zero temperature corresponds to r * which is independent of α, and thus matches the radius of the extremal EMliA solution. The electric conductivity at T = 0 is obtained by the replacement r 0 → r * in (4.18). It is finite and obviously agrees with the result in [20], whereas when T → ∞, σ goes to unity which is again the standard conducting behavior extracted from an EMliA bulk. Such a behavior has also been observed in the pertinent cases [23,51]. Noticeably, the result (4.18) satisfies the σ ≥ 1 bound proposed in [52] at all dissipation strengths. It is clear that since the new coupling does not enter into the leading order of the expansions about the two temperature extremes, one cannot expect deviations. To continue, we need to consider a time-dependent source for the heat current in our perturbation ansatz. This will allow us to compute the thermoelectric conductivities, α,ᾱ and the thermal conductivityκ at zero electric field, thus filling the remaining entries of the transport matrix.
We consider the ansätze (4.8) and (4.9), but now we switch on a time-dependent part in g tx , namely g tx = ǫ(tf 2 (r) + h tx (r)), while we make a more general ansatz for the gauge field; in particular, A x = ǫ(tf 1 (r) + A x (r)). The x-component of the gauge equations of motion is neatly written as a radial conservation law for the only non-vanishing component of the current density in the spatial directions, J x = r 2 F xr , Using the radial conservation of (4.19) together with the unperturbed field equations, we can manage to find a first r-integral of −2G tx ; namely the radially-constant quantity which we can identify with the x-component of the heat current of the boundary theory when evaluated at r → ∞. Again, since this is radially conserved, we can evaluate it at r 0 instead. Additionally, the rx-component of the metric field equations is an algebraic equation for the relevant perturbation which is solved by Indeed, as a consistency check, killing f 2 and setting f 1 = −E yields (4.12) as it should. We see that if we choose f 2 = −γF and f 1 = γA−E, all time dependence vanishes in J x and in the tx-component of the metric field equations, the former assuming the expression (4.10) while the latter becoming (4.13). In the Eddington-Finkelstein coordinate system, the regularity of the bulk perturbations and the satisfaction of the perturbed field equation G tx = 0 near the horizon radius, both boil down to the series expansion which in turn leads to the radially-constant quantities where S is defined in (3.10), S o := 4σπr 2 0 and σ is as in (4.18). Clearly, we have all the necessary transport coefficients of the strongly coupled theory, and we can now explicitly write down the generalized Ohm/Fourier law, from which we can read off the linear DC reponse of the system to an external electric field and a thermal gradient. Here, (4.26) are the thermoelectric conductivities and the thermal conductivity at zero electric field, respectively. First of all, the transport matrix (4.25) is symmetric which constitutes a successful consistency check against the Onsager relations [50] for theories invariant under time reversal, the latter relating the current densities of the background geometry to their counterparts obtained from a time-reversed solution. Secondly, when α → 0, the coefficients successfully reduce to those of EMliA theory obtained as a special example in [22].

Bounds of thermal conductivity and diffusion constants
In this subsection we wish to probe the theory against various relevant bounds in the holography-related literature. With the complete set of conductivities at hand we can work out some interesting relations from (4.25). First of all, let J ≡ J x and Q ≡ Q x . We have that which ultimately represents the conductivity in the absence of heat flows. In addition, the simple relation discussed in [22] is modified 5 ; From the transport matrix we can also define the thermal conductivity at zero electric current as Moreover, the Lorentz ratios will be of interest as well.
We observe that σ, α andκ blow up as p → 0, whereas κ goes to the finite value 4πS o T /µ 2 ,L → 4πS o /µ 2 and L → 0. Moreover, at zero temperature, we notice an electric conductor/thermal insulator behavior which is reminiscent of the findings in the much simpler linear axion model [20]. Furthermore, positivity of the temperature suggests that where G EFF is defined in (3.3). Thus, since S ≥ S o we find that We are interested in comparing the (κ/T )-bound we found with respect to the universal bound proposal in [53], κ/T ≥ 4π 2 /3. We reformulate (4.32) in terms of the dimensionless ratiosp := p/T andμ := µ/T , The relevant expansions about small and largeμ read respectively. We notice that for small µ/T , providedp lies in the coherent regime close to the order of α (such that α cannot be neglected), the bound has an α-and ap-scaledependence. For largep the ratio will eventually reach the 4π 2 /3 bound, but from below.  We collect this information in Fig. 1. In general, we find that the minimum of κ/T happens atp =μ and reads namely, it is scale-dependent; ergo, there does not exist a fixed O(1) number valid at all scales, and this is solely due to the non-minimal coupling and nothing else. However, the upper bound proposed in [22] still holds, i.e. (4.36) Another celebrated relation, the Kelvin formula, attributed to holographic models flowing towards an AdS 2 × R 2 fixed point in the IR [59,60], reads where Q ≡ Q e is the charge density in our case. We show that this holds true for the proposed model as well. First of all, the Seebeck coefficient α/σ at zero temperature is Then, we can use the chain rule in order to write (4.39) Using the inverse relation (3.12) together with Q e = µr 0 , and taking the T → 0 limit afterwards, we indeed arrive at (4.38) and the claim is proven. To close this subsection, we wish to investigate the thermoelectric response in the diffusion-dominated regime.
The incoherent limit is defined by p ≫ T, µ for a fixed finite ratio T /µ. For very strong dissipation, the transport coefficients expand as (4.40) with κ having the same leading order coefficient asκ. In this regime, where diffusion takes over, the horizon radius goes as ∝ p; in particular, r 0 = p η/6. One sees that the offdiagonal elements of (4.25) have an O(1/p)-falloff for large p, whereas the diagonal ones go to a finite value. The ratio of charged to neutral degrees of freedom measured by the Q e /S ratio goes as ∼ µ/(αp) → 0, and the charge/heat currents decouple [61]. A priori, we will not assume that charge and energy diffusion decouples; namely, we will not neglect the mixing term, M, defined below. The coupled diffusion is described [62] by the constants where χ, ζ and c Q are the charge susceptibility, thermoelectric susceptibility and specific heat (at fixed charge density Q), respectively. In whatever regime ζ, α = 0, the diffusitivities do, in fact, decouple with D + → D c := σ/χ (charge diffusion constant) and D − → D e := κ/c ρ (energy diffusion constant). Let us first compute the thermodynamic susceptibilities. We have , (4.44) whereC := 16π 2 T 2 + 24 √ 2παpT + 6ηp 2 + 3µ 2 ≡ 6r 0 − 4πT for safety of space. The specific heat at fixed Q is given by where is the specific heat at fixed chemical potential. Plugging everything back into (4.41), we can get an explicit expression for the diffusitivities at all p/T scales. The explicit expressions are too lengthy, hence we simply plot the results in Fig. 2. From the last mentioned Figure, subfigures (c) and (d) in particular, there are two observations to be made. Clearly, as µ/T increases, the mixing term has a decreasing     impact and the charge/energy diffusitivities do decouple. On the other hand, for small µ/T , subfigure (c) shows that the mixing term leads to an opposite identification compared to the one in the incoherent regime; D + is identified with D e , whereas D − is identified with D c . Also, it is apparent that regardless of the value of µ/T , the diffusion constants completely decouple in the incoherent limit, and it is safe to say that D + → D c and D − → D e when dissipation becomes strong. To see this clearly, subfigure (c) shows that the mixing term has no impact for p/T 10. This corresponds to p/µ 10 3 ≫ 1, i.e. this lies in the diffusion-dominated region. Indeed, the large p-expansions of D ± read exhibiting a leading order agreement with the expansions of D c and D e , respectively, in the incoherent limit. These results have also been derived in [62], in much greater detail. Apart from a mild curve shifting in the coherent regime, the new coupling has no effect in the strong dissipation limit; it does not contribute to the leading order in (4.47), and it has no qualitative effect, either in the maximal mixing scenario at small µ/T or in the incoherent decoupling limit. After all this song and dance, the ultimate aim is to see if the new coupling affects the diffusitivity bound proposal in [63,67]. According to the Hartnoll conjecture [64] for v being some characteristic velocity. Instead of the original idea to match the latter with the speed of light, a reasonably natural candidate for v at strong coupling has been the butterfly velocity, a measure of the spatial propagation speed of chaos through the dual quantum system. This has been derived in [63] for a general IR geometry with matter coupling to an Einstein bulk. The holographic derivation depends on the black hole horizon data and its geometric picture is that of a shock-wave propagating in the bulk; the butterfly effect is manifest through the exponential boosting of the energy of an in-falling particle near the horizon, at late times. For more details, please see [65,66]. Following [62], we are interested in the -expected to be-O(1) numbers B ± which act as a lower diffusion bound where instead of using the decoupled charge/energy diffusion constants, we will work with the D ± constants as well; this will allow us to see in what regions the universality of the bound holds, and moreover, if the new coupling has any impact at all. A lightning quick calculation with Mathematica using standard methods in [63,65,66], and the more relevant [62] in particular, reveals that the screening length m is not modified and the uu-component of the perturbed equations at linearized order reads with m 2 ∼ ∂V (uv)/∂(uv) /A(0)| uv=0 . Here, the expression is in Kruskal coordinates (u, v) with A, V functions of uv, f is some -irrelevant to the solution-function with α being part of its arguments, β = 1/T , t w is the past time the particle was released on the boundary of AdS 6 , and E 0 is the initial energy of that particle. Since the screening length is not modified, there will be no deviation, either in the value of the Lyapunov exponent λ L , or in the expression for the butterfly velocity v B . The solution to (4.51) has essentially the same form as if the bulk was pure Einstein gravity. Thus, comparing it with the exponential formula determining the growth in the commutators of generic Hermitian local operators (see relation (4) in [63] for example), we deduce that λ L = 2π/β and   (c) µ/T = 0.01 and α = 0.5 the latter expression in terms of the dimensionless variablesμ := µ/T andp := p/T . In the incoherent limit, the butterfly velocity goes to 0 as ∼ 1/p, and hence, since α has no contribution to the leading order of the strong dissipation expansion of the diffusion constants as well, we conclude that it will be impossible for α to affect the universal bounds B + = 2 + ... and B − = 1 + ..., neither does it have any impact in the low temperature expansion. These conclusions are also graphically manifest in Fig. 3.

Shear viscosity to entropy density ratio via a (weaker) horizon formula
To compute the shear viscosity to entropy density ratio, η/S, we employ the method devised in [24]. We consider the bulk metric perturbation g xy = ǫr 2 h(r)e −iωt of the eigenmode type about the black hole background (2.30). This leads to the second order ODE, G xy = 0, which can be written as The explicit expression of the mass function will be stated at a later point. The shear viscosity is then computed in terms of the correlator: where h o is the solution to (4.53) at zero frequency ω = 0, which (i) is regular at r 0 and (ii) goes like unity near radial infinity. Then, the shear viscosity to entropy density ratio is Obviously, when α = 0, the entropy reduces to S o , and G EFF = (16π) −1 ; the expression for the ratio assumes the standard form [24]. Moreover, since fluctuations of g xy are massive, and assuming a positive effective mass squared, we know that h o (r 0 ) < 1 which follows from a simple argument illustrated in [24]. Since G EFF < 1 strictly for non-trivial α, we already know that the simple (4π) −1 bound is definitely violated at finite temperatures. The mass squared is given by where again F and its derivatives are understood to be on the background shell. First of all, we observe that as α → 0 we recover the standard mass termη(p/r) 2 . Moreover, (4.56) will be strictly positive at r 0 . It will be also finite positive in the T → 0 limit where the black hole becomes extremal with r 0 = r * . There is no general argument why the effective mass term m 2 needs to be positive in general, but in our case it so happens that it is a strictly positive function of the radial coordinate in the physical domain of interest r 0 < r < ∞. Now, we define b = p/r 0 , and we notice that (4.53) with ω = 0 has already terms linear in b; indeed, it reads where Here, a convenient change of the radial coordinate, z = r 0 /r, was performed, such that the horizon and boundary are located at z = 1 and z = 0, respectively. Also,η was set to unity. We observe that due to the new coupling there are odd powers of b introduced in the differential equation. Treating b perturbatively, if we were to expand the solution as plugging it back into (4.57) and solving the ODE order by order, we would find inconsistencies already at order b; ergo, a general expansion of the form is necessary. All we need to do now is solve order by order. We remind the reader that we work at zero chemical potential. At zeroth order we need to obtain a solution to where a prime denotes differentiation with respect to z. The general solution to this reads Regularity at the horizon suggests that c 1 =0, while the boundary conditions at z = 0 imply that c 2 = 1. Hence, h o,0 = 1 and we need to make all other h o,n , n > 0, vanish at z = 0 so that the asymptotic behavior of h o meets condition (ii). Moving on to linear order in b, we find that the solution to compatible with the aforementioned conditions, reads At second b-order, the solution is already too lengthy to write down. It involves many dilogarithmic and arctangent functions. All contributions due to the non-minimal coupling go as α 2 and schematically, Higher orders k will go as the solutions in [24] plus α k corrections if k even, while if k is odd, the solution will have an overall α k factor, such that when we switch off the coupling constant we recover the h o of the linear axion model. Just as a minor example, let us try a very crude and inelegant approximation of the η/S ratio in the high temperature regime. We find that Unfortunately, this approximation, being so crude, is not very helpful; a numerical solution to (4.57) is certainly necessary in order to obtain a better insight. This can be seen in Fig 4. In subfigure (b), for the approximation plots corresponding to the cases with non-vanishing α, we have neglected the α 2 contributions in (4.69) in order to avoid using information that is not explicitly presented here. Since this contribution comes with an overall negative   sign in the end, it is of course expected that the curves fit even better to the numerical results in the very small p/T region (where the approximation applies). From subfigure (a), apart from a manifestly more brutal violation of the bound, it is apparent that the ratio exhibits a qualitatively similar behavior at both temperature extremes, always in comparison to the results obtained in the case of the linear axion model. It goes to unity as T → ∞, while it tends to zero when T → 0. In general, the violation itself is nothing really unexpected 7 , mainly because the shear mode mass is everywhere non-vanishing. Albeit that, it is nevertheless interesting to visualize the manner in which the Horndeski coupling enters the scene.

Concluding remarks
In this work, we have started by considering a specific model of Horndeski gravity with G 3 , G 5 = 0. We deliberately chose G 2 ∝ X = −(∂ψ) 2 /2 and G 4 ∝ α √ −X such that we meet condition (iv), but fail condition (v), both mentioned in the second paragraph of the introduction. We took (D − 2)-copies of the model such that the massless scalar fields ψ I = pδ I i x i "dress" the (D − 2)-many planar directions; the final action can be intuitively expressed as Einstein gravity with a running (with √ −X 1 + √ −X 2 + ...) effective gravitational coupling plus matter fields with higher derivatives accompanied by a Maxwell term. By doing so, we managed to circumvent the restrictive results of [7], finding charged planar black holes with a non-trivial axion profile. We studied the horizon structure of the four-dimensional solution with AdS asymptotics, revealing a mass region in which the PBHs have two horizons which coalesce into one in the extremal case. A near horizon AdS 2 × R 2 structure was observed, whereas these solutions asymptote standard unit-radius AdS 4 . A straightforward dyonic extension of these black holes was given, while we also exhibited the three-dimensional solution which does not flow from limiting the D-dimensional result and requires separate integration of the field equations. We proceeded by studying the thermodynamic properties of (2.21) where the entropy was derived via two routes: first, using Wald's Noether charge entropy formula, and second, via the conventional Euclidean path integral approach. Contrary to the discrepancies advertised in the beginning of Sec. 3, we found that both methods agree on the result; the entropy does not obey the 1/4-area law, although it does so in some units where (3.3) equals unity. In this sense and to some extent, hints to the 1/4-area law are still there. Expressions for the mass and the Hawking temperature were also provided, and the first law was shown to hold true.
Next, we used the method devised in [22,47] to compute the linear thermoelectric DC response of the holographic dual to some external electric field and some thermal gradient; this was done by means of black hole horizon data only, exploiting the radial conservation of the electric/heat currents. Analytic expressions were found for the electric and heat currents, along with a detailed extraction of the elements of the DC transport matrix. As a consistency check, we verified that the matrix was symmetric, a consequence of invariance under time reversal, while we found that theκ/α ratio stated in [22] was modified, being (4.28). The comparison is always with respect to the linear axion model [20]. The behavior of the electric conductivity at the two termperature extremes (or in the incoherent regime) was not altered by the presence of the Horndeski coupling, a fact suggesting that it is rather governed by the choice of electrodynamics instead; it would be interesting to consider different types of non-linear electrodynamics [54][55][56][57][58]70] coupled to this gravitational toy model, and investigate how the parameters mingle with each other. Then, knowing the linear responses of the system, we computed the thermal conductivity κ at zero electric current. We found that there is no scale-independent O(1) number bounding the κ/T ratio from below, although in the incoherent regime, the strong dissipation expansion was 4π 2 /3 + ..., where dots stand for subleading terms. However, we argued that the bound depends on the p/T value in general due to the non-minimal coupling of the axion fields to gravity, the minimum happening at p = µ for a given α, being (4.35), whereas for increasing µ/T , we saw that the relevant curve approaches the 4π 2 /3-line (see Fig. 1 and (4.34)). Nevertheless, the upper bound κ/T ≤ S 2 /Q 2 e was still found to hold good. The Kelvin formula for the Seebeck coefficient at zero temperature was also verified, in favor of the argument [59,60] that its validity is related to the flow of holographic models towards an AdS 2 × R 2 fixed point in the IR.
Next, we considered a generalized version of Einstein's relation where the diffusitivities are mixed by a term M [64] with the constants describing the coupled diffusion being D ± . We explicitly showed that for this model the diffusion constants D ± do indeed decouple to the charge and energy diffusitivities, respectively, only in the incoherent regime where the mixing term becomes negligible. As observed in [62] as well, the mixing becomes maximal for small µ/T , in the sense that there is an opposite identification between charge/energy diffusion constants and D ± , that is, opposite to the way the identification is realized in the strong dissipation limit. The new coupling had no influence on the decoupling process, neither did it have any contribution to the leading order terms in the incoherent expansion of the diffusion constants (or in the low temperature expansion, respectively). Considering the Blake refinement [63] of the T D/v 2 lower bound conjectured in [64], where v = v B , the butterfly velocity, we calculated the latter only to find that the new coupling did not alter the shift equation drastically, nor did it have any impact on the screening length m. The velocity obeyed the generic formula obtained for a pure Einstein bulk with minimallycoupled matter [63]; the effect of the new coupling entered the velocity only through the explicit expression of the horizon (3.12). Since the behavior of r 0 in the incoherent limit, or for low temperatures, did not depend on α, we concluded and graphically demonstrated that T D ± /v 2 B was eventually bounded by the standard numbers from below, for large p/T that is. The low-temperature behavior was also identical to the one exhibited in EMliA theory. Next, we employed the weaker horizon formula [24] to determine the shear viscosityto-entropy density ratio, η/S, at zero chemical potential. Since the g xy fluctuations were massive, with a positive effective mass squared given by (4.56), it was not a surprise that the model led to a violation of the simple bound (4π) −1 . Performing a very crude approximation, we noticed that, at least in the high temperature regime, the new coupling did not only contribute to the even-power subleading terms of the ratio expansion, but allowed for odd-power corrective terms, the latter entirely imputed to the presence of α. We also managed to solve (4.57) numerically, displaying η/S in Fig. 4. In comparison to the linear axion model, a more brutal violation of the KSS bound was observed, accredited to the Horndeski coupling α. However, the behavior at both temperature extremes was found to be similar with the ratio going to zero (unity) when T → 0 (T → ∞). As a closing remark, we note that it would be very interesting to further investigate this model in the relevant context of [74,75].