Black holes in vector-tensor theories

TL;DR

This paper investigates black hole solutions in generalized Proca theories with nonminimal vector couplings, finding exact and numerical solutions that include vector hair, which could be tested through gravitational wave observations.

Contribution

It provides new exact and numerical black hole solutions in vector-tensor theories with nonminimal derivative couplings, highlighting the existence of vector hair and potential observational signatures.

Findings

Exact solutions include RN, stealth Schwarzschild, and extremal RN types.

Numerical solutions show black holes with non-trivial metric component differences.

Vector hair can be secondary and regular at the horizon.

Abstract

We study static and spherically symmetric black hole (BH) solutions in second-order generalized Proca theories with nonminimal vector field derivative couplings to the Ricci scalar, the Einstein tensor, and the double dual Riemann tensor. We find concrete Lagrangians which give rise to exact BH solutions by imposing two conditions of the two identical metric components and the constant norm of the vector field. These exact solutions are described by either Reissner-Nordstr\"{o}m (RN), stealth Schwarzschild, or extremal RN solutions with a non-trivial longitudinal mode of the vector field. We then numerically construct BH solutions without imposing these conditions. For cubic and quartic Lagrangians with power-law couplings which encompass vector Galileons as the specific cases, we show the existence of BH solutions with the difference between two non-trivial metric components. The quintic-order power-law couplings do not give rise to non-trivial BH solutions regular throughout the horizon exterior. The sixth-order and intrinsic vector-mode couplings can lead to BH solutions with a secondary hair. For all the solutions, the vector field is regular at least at the future or past horizon. The deviation from General Relativity induced by the Proca hair can be potentially tested by future measurements of gravitational waves in the nonlinear regime of gravity.

Equations

$$S=\int d^{4}x\sqrt{-g}\left(F+\sum _{i=2}^6{\mathcal{L}}_i\right),$$

$${\mathcal{L}}_2$$

$${\mathcal{L}}_3$$

$${\mathcal{L}}_4$$

$${\mathcal{L}}_5$$

$${\mathcal{L}}_6$$

$$L^{\mu \nu \alpha \beta }=\frac{1}{4}{\mathcal{E}}^{\mu \nu \rho \sigma }{\mathcal{E}}^{\alpha \beta \gamma \delta }{R}_{\rho \sigma \gamma \delta },$$

$$d{s}^2=-f(r)d{t}^2+h^{-1}(r)d{r}^2+r^2(d{\theta }^2+{\sin }^2\theta d{\phi }^2),$$

$$A_{\mu }=\left(A_0(r),A_1(r),0,0\right),$$

$${\mathcal{L}}_2=G_2(X,F),$$

$$X_0=\frac{A_0^2}{2f},X_1=-\frac{h{A}_1^2}{2}.$$

$$rf\left[2fh(r{A}_0^{\prime \prime }+2A_0^{\prime })+r(f{h}^{\prime }-f^{\prime }h)A_0^{\prime }\right](1+G_{2,F})+r^{2}h{A}_0^{\prime 2}\left[2fh{A}_0^{\prime \prime }-(f^{\prime }h-f{h}^{\prime })A_0^{\prime }\right]G_{2,FF}-2r^2{f}^2{A}_0{G}_{2,X}$$

$$-2r^{2}f{A}_0^{\prime }\left(f{h}^2{A}_1{A}_1^{\prime }-h{A}_0{A}_0^{\prime }+f^{\prime }h{X}_0-f{h}^{\prime }{X}_1\right)G_{2,XF}-rf{A}_0\left[2rfh{A}_1^{\prime }+(r{f}^{\prime }h+rf{h}^{\prime }+4fh)A_1\right]G_{3,X}$$

$$+4f^2{A}_0(r{h}^{\prime }+h-1)G_{4,X}-8f{A}_0\left[rf{h}^2{A}_1{A}_1^{\prime }-(r{f}^{\prime }h+rf{h}^{\prime }+fh)X_1\right]G_{4,XX}$$

$$-f{A}_0\left[f(3h-1)h^{\prime }{A}_1+h(h-1)(f^{\prime }{A}_1+2f{A}_1^{\prime })\right]G_{5,X}-2fh{A}_0{X}_1\left[2fh{A}_1^{\prime }+(f^{\prime }h+f{h}^{\prime })A_1\right]G_{5,XX}$$

$$-2f\left[f(3h-1)h^{\prime }{A}_0^{\prime }+h(h-1)(2f{A}_0^{\prime \prime }-f^{\prime }{A}_0^{\prime })\right]G_6-4fh{A}_0^{\prime }{X}_1\left(h{A}_0{A}_0^{\prime }-2f{h}^2{A}_1{A}_1^{\prime }-2f^{\prime }h{X}_0+2f{h}^{\prime }{X}_1\right)G_{6,XX}$$

$$-2f\left[4f{h}^2{X}_1{A}_0^{\prime \prime }-2h(hX-X_0)f^{\prime }{A}_0^{\prime }+2f(6h-1)h^{\prime }{X}_1{A}_0^{\prime }+h(h-1)A_0{A}_0^{\prime 2}-2f{h}^2(3h-1)A_0^{\prime }{A}_1{A}_1^{\prime }\right]G_{6,X}$$

$$-4fh\left[2rfh{A}_1{A}_0^{\prime \prime }-\{(r{f}^{\prime }h-3rf{h}^{\prime }-2fh)A_1-2rfh{A}_1^{\prime }\}{A}_0^{\prime }\right]g_5$$

$$-4rfh{A}_0^{\prime }\left[h{A}_0{A}_0^{\prime }{A}_1+4fh{X}_1{A}_1^{\prime }-2A_1(f^{\prime }h{X}_0-f{h}^{\prime }{X}_1)\right]g_{5,X}=0.$$

$$A_1\left[r^{2}f{G}_{2,X}-2(r{f}^{\prime }h+fh-f)G_{4,X}+4h(r{A}_0{A}_0^{\prime }-r{f}^{\prime }X-f{X}_1)G_{4,XX}-h{A}_0^{\prime 2}(3h-1)G_{6,X}-2h^2{X}_1{A}_0^{\prime 2}{G}_{6,XX}\right]$$

$$=r[r(f^{\prime }X-A_0{A}_0^{\prime })+4f{X}_1]G_{3,X}+2f^{\prime }h{X}_1{G}_{5,X}+(A_0{A}_0^{\prime }-f^{\prime }X)\left[(1-h)G_{5,X}-2h{X}_1{G}_{5,XX}\right]$$

$$-2rh{A}_0^{\prime 2}(g_5+2X_1{g}_{5,X}).$$

$$\left(c_1+\frac{c_2}{r}+\frac{c_3}{r^2}\right)h^{\prime }+c_4+\frac{c_5}{r}+\frac{c_6}{r^2}=0,$$

$$-\frac{h}{f}\left(c_1+\frac{c_2}{r}+\frac{c_3}{r^2}\right)f^{\prime }+c_7+\frac{c_8}{r}+\frac{c_9}{r^2}=0,$$

$$\left(c_{10}+\frac{c_{11}}{r}\right)f^{\prime \prime }+\left(c_{12}+\frac{c_{13}}{r}\right)f^{\prime 2}+\left(\frac{c_2}{2f}+\frac{c_{14}}{r}\right)f^{\prime }{h}^{\prime }+\left(c_{15}+\frac{c_{16}}{r}\right)f^{\prime }+\left(-\frac{c_8}{2h}+\frac{c_{17}}{r}\right)h^{\prime }+c_{18}+\frac{c_{19}}{r}=0,$$

$$G_4=\frac{M_{\mathrm{pl}}^2}{2},G_2=G_3=G_5=G_6=0,g_5=0,$$

$$2fh\left(r{A}_0^{\prime \prime }+2A_0^{\prime }\right)+r\left(f{h}^{\prime }-h{f}^{\prime }\right)A_0^{\prime }=0,$$

$$2f(r{h}^{\prime }+h-1)M_{\mathrm{pl}}^2+r^{2}h{A}_0^{\prime 2}=0,$$

$$2\left[rh{f}^{\prime }+(h-1)f\right]M_{\mathrm{pl}}^2+r^{2}h{A}_0^{\prime 2}=0,$$

$$r{A}_0^{\prime \prime }+2A_0^{\prime }=0,$$

$$f^{\prime }=-\frac{r^2{A}_0^{\prime 2}+2M_{\mathrm{pl}}^2(f-1)}{2M_{\mathrm{pl}}^{2}r},$$

$$f=h=1-\frac{2M}{r}+\frac{Q^2}{2M_{\mathrm{pl}}^2{r}^2},$$

$$A_0=P+\frac{Q}{r},$$

$$G_4=\frac{M_{\mathrm{pl}}^2}{2},G_2=m^{2}X,G_3=G_5=G_6=0,g_5=0,$$

$$A_1=0.$$

$$2fh\left(r{A}_0^{\prime \prime }+2A_0^{\prime }\right)+r\left(f{h}^{\prime }-h{f}^{\prime }\right)A_0^{\prime }-2m^{2}fr{A}_0=0,$$

$$2f(r{h}^{\prime }+h-1)M_{\mathrm{pl}}^2+r^{2}h{A}_0^{\prime 2}+m^2{r}^2{A}_0^2=0,$$

$$2\left[rh{f}^{\prime }+(h-1)f\right]M_{\mathrm{pl}}^2+r^{2}h{A}_0^{\prime 2}-m^2{r}^2{A}_0^2=0.$$

$${\left(\frac{f}{h}\right)}^{\prime }=\frac{m^2{A}_0^{2}r}{M_{\mathrm{pl}}^2{h}^2}.$$

$$A_0=0,$$

$$f=h=1-\frac{2M}{r},$$

Full text

I Introduction

The direct detection of gravitational waves by Advanced LIGO LIGO opened up a new opportunity for probing the physics of black holes (BHs) and their surroundings. The Solar System experiments of gravity have shown that General Relativity (GR) holds in high precision in the weak gravitational regime of the local Universe Will . The future gravitational wave experiments will allow us to test the validity of GR and possible deviations from GR in the nonlinear regime of gravity Vitor ; Yagi .

From the cosmological side, the observational evidence of dark energy and dark matter SNIa ; CMB implies that there may be some extra propagating degrees of freedom (DOFs) beyond the realm of GR. In particular, the infrared modifications of gravity with new DOFs have been extensively studied as sources for the late-time cosmic acceleration moreview . In the local Universe with weak gravitational backgrounds, the propagation of fifth forces mediated by new DOFs can be suppressed under the operation of screening mechanisms, e.g., Vainshtein Vain or chameleon chame mechanisms. In the strong-gravity regime associated with BHs and neutron stars (NSs), the behavior of new DOFs is generally more involved due to tensor nonlinearities.

In the Einstein-Maxwell system of GR without matter, there is a uniqueness theorem stating that the asymptotically flat and stationary BH solutions are described only by three parameters, i.e., mass, electric charge, and angular momentum Israel ; Carter ; Wheeler ; Hawking . The “no-hair” BH theorem is valid for a canonical scalar field minimally coupled to gravity Chase ; BekenPRL . The same no-hair property also persists for standard scalar-tensor theories in which the field has a direct coupling to the Ricci scalar Hawking72 ; Beken95 ; Soti12 . However, the no-hair theorem of Ref. Beken95 loses its validity in modified gravitational theories with nonminimal derivative couplings to gravity. The typical examples of such derivative couplings are Galileons Galileon1 ; Galileon2 , whose equations of motion respect the Galilean symmetry in the Minkowski limit. The extension of Galileons to more general couplings led to the rediscovery of Horndeski theories redis –most general scalar-tensor theories with second-order equations of motion Horndeski .

In the case of shift-symmetric Horndeski theories including Galileons, Hui and Nicolis Hui argued conditions for the no-hair properties of BHs by utilizing the properties of the conserved Noether current $J^{\mu}$. The existence of the shift symmetry gives rise to the field equation $\nabla_{\mu}J^{\mu}=0$, where $\nabla_{\mu}$ denotes the covariant derivative. Under the assumptions that i) the spacetime is static, spherically-symmetric and asymptotically flat, ii) the scalar field $\phi$ respects the symmetry of spacetime, i.e., $\phi=\phi(r)$, where $r$ is the distance from the center of symmetry, and iii) the scalar product $J_{\mu}J^{\mu}$ is regular everywhere, it can be shown that the radial current $J^{r}$ needs to vanish at all the distance $r$ due to the regularity on the horizon, $J^{r}=0$. They further employed the fact that the current takes the form $J^{r}=\phi^{\prime}g^{rr}{\cal F}(\phi^{\prime};g,g^{\prime},g^{\prime\prime})$, where $\phi^{\prime}=d\phi/dr$, $g^{rr}$ is the radial component of the metric $g^{\mu\nu}$, and ${\cal F}$ is a function containing $\phi^{\prime}$ and derivatives of $g^{\mu\nu}$. Provided ${\cal F}$ does not contain negative powers of $\phi^{\prime}$ so that the canonical kinetic term dominates in the asymptotically flat region, $g^{rr}$ and ${\cal F}$ approach non-vanishing constant values and hence $\phi^{\prime}=0$ at the infinity. Moving inward from infinity, $g^{rr}$ and ${\cal F}$ vary continuously taking non-zero values, so only the allowed profile consistent with $J^{r}=0$ is the no-hair solution satisfying $\phi^{\prime}=0$. In other words, in order to have a non-trivial hairy BH solution, one has to break at least one of the assumptions made by Hui and Nicolis.

If the constancy of ${\cal F}$ is not imposed in the limit $\phi^{\prime}\to 0$, there exist some hairy BH solutions in shift-symmetric Horndeski theories with non-vanishing values of $\phi^{\prime}$ Rinaldi ; Anabalon ; Minami13 ; Soti1 ; Soti2 ; Babi17 . In the case where the scalar field is linearly coupled to a Gauss-Bonnet term, the function ${\cal F}$ contains a negative power of $\phi^{\prime}$, which allows the existence of solutions with $\phi^{\prime}\neq 0$ Soti1 ; Soti2 . Another approach to construct BH solutions in Horndeski theories is to assume a linearly time-dependent scalar field $\phi=qt+\psi(r)$ Babi14 , where ${\cal F}=0$ is ensured as the consequence of the equations of motion, leaving $\phi^{\prime}$ unfixed by the condition $J^{r}=0$. This leads to a family of BH solutions with a static metric Koba14 ; Babi16 ; Lefteris , especially the stealth Schwarzschild BH solution Babi14 .

If we consider a massless vector field $A_{\mu}$ with the Lagrangian $F=-F_{\mu\nu}F^{\mu\nu}/4$, where $F_{\mu\nu}=\nabla_{\mu}A_{\nu}-\nabla_{\nu}A_{\mu}$ is the field strength tensor, the resulting static and spherically symmetric solution in GR is given by the Reissner-Nordström (RN) metric with mass $M$ and electric charge $Q$. For a massive vector field with the Lagrangian $-m^{2}A^{\mu}A_{\mu}/2$, the $U(1)$ gauge symmetry is explicitly broken, so there is a longitudinal propagation besides two transverse polarizations. In the framework of GR, Bekenstein Beken72 showed that this massive Proca field $A_{\mu}$ needs to vanish due to the regularity of a physical scalar constructed from $A_{\mu}$ on the horizon. Hence the resulting solution is given by the Schwarzschild solution without the vector hair.

This no-hair theorem for the massive Proca field cannot be applied to theories with vector derivative couplings. The action of generalized Proca theories with nonminimal derivative couplings to gravity was first constructed in Refs. Heisenberg ; Tasinato from the demand of keeping the three propagating DOFs besides two tensor polarizations. The theories can be further extended Allys to include intrinsic vector-mode couplings with the double dual Riemann tensor $L^{\mu\nu\alpha\beta}$ Jimenez2016 , such that the $U(1)$-invariant interactions derived by Horndeski Horndeski76 can be accommodated as a specific case. The equations of motion in these theories remain of second order, but one can build more general vector-tensor interactions beyond the second-order domain without introducing extra DOFs associated with the Ostrogradski instability HKT ; Kimura . In second-order generalized Proca theories and their extensions the derivative interactions can drive the late-time cosmic acceleration DeFelicecosmo with some distinct observational signatures obsig1 ; obsig2 , while satisfying local gravity constraints in the Solar System DeFeliceVain . See Ref. Heisenberg:2017mzp for a short summary.

The study of hairy BH solutions in generalized Proca theories with an Abelian vector field $A_{\mu}$ has recently received attention for probing physics in the nonlinear regime of gravity Chagoya ; Fan ; Minami ; Cisterna ; Chagoya2 ; Babichev17 ; HKMT , see also Refs. colored for early works of BH solutions in the presence of non-Abelian Yangs-Mills fields. In theories whose Lagrangians contain the coupling $\beta_{4}G^{\mu\nu}A_{\mu}A_{\nu}$ as well as $-F^{\mu\nu}F_{\mu\nu}/4$ and the Einstein-Hilbert term, Chagoya et al. Chagoya derived an exact spherically symmetric and static BH solution for the specific coupling $\beta_{4}=1/4$. This exact BH solution was further extended to asymptotically non-flat solutions Minami ; Babichev17 , non-exact solutions for $\beta_{4}\neq 1/4$ Babichev17 ; Chagoya2 , rotating solutions Minami , and NSs Chagoya2 .

On a static and spherically symmetric background with the radial coordinate $r$, the vector field is characterized by the temporal component $A_{0}(r)$ and the longitudinal mode $A_{1}(r)$. In general, the equation of motion for $A_{1}$ can be written in the form ${\cal F}(A_{1},A_{0},A_{0}^{\prime};g,g^{\prime})=0$. Unlike scalar-tensor theories, the presence of two vector components naturally allows the solution $A_{1}\neq 0$ without restricting the functional form of ${\cal F}$. Hence it is not difficult to find hairy BH solutions even for simple power-law couplings like those of vector Galileons HKMT . For some derivative interactions the equation for $A_{1}$ reduces to the form $A_{1}\tilde{\cal F}(A_{1},A_{0},A_{0}^{\prime};g,g^{\prime})=0$, so there is the branch $A_{1}=0$ besides $\tilde{\cal F}=0$. Even with the branch $A_{1}=0$, the deviation from GR can arise in the metric components due to a modification of the temporal component $A_{0}$ induced by derivative couplings, e.g., the BH solution arising from the $U(1)$-invariant interaction HorndeskiBH .

In this paper, we will present a detailed study of BH solutions in second-order generalized Proca theories by extending the analysis of Ref. HKMT . We will consider the full set of Lagrangians ${\cal L}_{2,3,4,5,6}$ of the generalized Proca theories Heisenberg ; Jimenez2016 . Our analysis also covers the generalized quadratic-order Lagrangian $G_{2}(X,F,Y)$, where $X=-A_{\mu}A^{\mu}/2$ and $Y=A^{\mu}A^{\nu}{F_{\mu}}^{\alpha}F_{\nu\alpha}$. We also explain in details how to construct non-exact BH solutions in power-law coupling models containing the dependence of $X^{n}$ in each Lagrangian, where $n$ is a positive integer.

We organize our paper as follows. In Sec. II we present the full equations of motion on a static and spherically symmetric background and revisit the Bekenstein’s no-hair BH solution for a massive Proca field. In Sec. III we review the exact BH solution present for the quartic derivative coupling $G_{4}(X)$. In Sec. IV we construct a family of exact BH solutions in the presence of other couplings under the conditions that the two metric components are identical and that the norm $X$ of the vector field is constant. In Sec. V we study the BH solutions for cubic-order power-law coupling models $G_{3}(X)\propto X^{n}$ including vector Galileons ($n=1$) and numerically confirm the existence of regular hairy BH solutions outside the horizon. Similarly, in Secs. VI, VII, VIII, IX we clarify the cases in which the BH solutions with primary or secondary Proca hairs are present for the power-law models containing the $X^{n}$ dependence in the couplings $G_{4}(X),G_{5}(X)$ and the intrinsic vector-mode couplings $G_{6}(X),g_{4,5}(X)$, respectively. The last Sec. X is devoted to conclusions.

II Generalized Proca theories

In the presence of a vector field $A_{\mu}$ with the field strength tensor $F_{\mu\nu}=\nabla_{\mu}A_{\nu}-\nabla_{\nu}A_{\mu}$, we consider second-order generalized Proca theories given by the action

[TABLE]

where $F=-F_{\mu\nu}F^{\mu\nu}/4$ is the standard Maxwell term, $g$ is a determinant of the metric tensor $g_{\mu\nu}$, and Heisenberg ; Jimenez2016

[TABLE]

with $X=-A_{\mu}A^{\mu}/2$, $Y=A^{\mu}A^{\nu}{F_{\mu}}^{\alpha}F_{\nu\alpha}$ and $G_{i,X}=\partial G_{i}/\partial X$. The functions $G_{3,4,5,6}$ and $g_{5}$ depend on $X$ alone, whereas the function $G_{2}$ is dependent on $F$ as well as $X$ and $Y$. The quantity $\tilde{F}^{\mu\nu}$ is the dual strength tensor given by $\tilde{F}^{\mu\nu}=\mathcal{E}^{\mu\nu\alpha\beta}F_{\alpha\beta}/2$, where $\mathcal{E}^{\mu\nu\alpha\beta}$ is the Levi-Civita tensor satisfying the normalization $\mathcal{E}^{\mu\nu\alpha\beta}\mathcal{E}_{\mu\nu\alpha\beta}=-4!$. The vector field has derivative couplings to the Ricci scalar $R$, the Einstein tensor $G_{\mu\nu}$, and the double dual Riemann tensor $L^{\mu\nu\alpha\beta}$ defined by

[TABLE]

where $R_{\rho\sigma\gamma\delta}$ is the Riemann tensor.

The original Proca theory with the mass term $m$ corresponds to the Lagrangian $G_{2}=m^{2}X$, in which case the longitudinal propagation arises besides two transverse polarizations. The action (1) has been constructed to keep the propagating DOFs of scalar and vector modes unchanged, i.e., three DOFs. Taking the scalar limit $A_{\mu}\to\nabla_{\mu}\pi$, the quantity $F$ as well as the Lagrangians $-g_{5}(X)\tilde{F}^{\alpha\mu}\tilde{F}^{\beta}_{~{}\mu}\nabla_{\alpha}A_{\beta}$ and ${\cal L}_{6}$ vanish, so they correspond to intrinsic vector modes. In Ref. HKMT there is the term $-2g_{4}(X)F$ in ${\cal L}_{4}$, but such a term is now absorbed into the Lagrangian ${\cal L}_{2}=G_{2}(X,F,Y)$.

To study BH solutions on a static and spherically symmetric background, we take the line element

[TABLE]

where $t$, $r$ and $(\theta,\varphi)$ represent the time, radial, and angular coordinates, respectively, $f(r)$ and $h(r)$ are functions of $r$ such that $f(r)>0$ and $h(r)>0$ outside the event horizon $r>r_{h}$, $r_{h}$ is the position of the horizon at which $f(r_{h})=h(r_{h})=0$. Expressing the vector field in the form $A_{\mu}=(A_{0},A_{i})$, the spatial vector $A_{i}$ can be decomposed into the transverse and longitudinal components, as $A_{i}=A_{i}^{(T)}+\nabla_{i}\chi$, where $A_{i}^{(T)}$ obeys the transverse condition $\nabla^{i}A_{i}^{(T)}=0$ and $\chi$ is the longitudinal scalar. From the regularity at the origin of the static and spherically symmetric background the transverse mode $A_{i}^{(T)}$ needs to vanish DeFeliceVain , so we are left with the longitudinal scalar $\chi$ in $A_{i}$. Hence the vector-field profile compatible with the background (8) is given by

[TABLE]

where $A_{1}(r)=\chi^{\prime}(r)$, and a prime represents the derivative with respect to $r$.

At this stage, we would like to make a comment about the additional intrinsic vector-mode contribution $Y$ in ${\cal L}_{2}=G_{2}(X,F,Y)$. On the static and spherically symmetric background (8) with the vector components (9) we have that $Y=4FX$, so the additional dependence of $Y$ in Eq. (2) can be removed. Thus, we will work on the quadratic-order Lagrangian

[TABLE]

in the rest of the paper.

The term $X=-A_{\mu}A^{\mu}/2$ can be expressed as $X=X_{0}+X_{1}$, where

[TABLE]

Varying the action (1) with respect to $A_{0}$ and $A_{1}$ respectively gives rise to the vector-field equations. The equation of motion for $A_{0}$ results in

[TABLE]

And similarly, the equation of motion for $A_{1}$ reads

[TABLE]

For the theories in which only the couplings $G_{i}(X)$ with even index $i$ are present, Eq. (13) admits the branch of the solution $A_{1}=0$. This is not the case for the theories containing the couplings $G_{3},G_{5}$ and $g_{5}$. To derive the gravitational equations of motion, we write the metric (8) in a more general form $ds^{2}=-f(r)dt^{2}+h^{-1}(r)dr^{2}+r^{2}e^{2\zeta(r)}(d\theta^{2}+\sin^{2}\theta\,d\varphi^{2})$ and express the action (1) in terms of $f,h,\zeta,A_{0},A_{1}$. Varying the resulting action with respect to $f,h,\zeta$ and setting $\zeta=0$ in the end, we obtain the following equations of motion

[TABLE]

where the coefficients $c_{1,2,\cdots,19}$ are given in Appendix. Among Eqs. (12)-(16) four of them are independent, so we will mostly employ Eqs. (12), (13), (14), and (15) in the discussions below.

II.1 RN solutions in GR

As a warm up, we first review the solutions in GR characterized by the functions

[TABLE]

where $M_{\rm pl}$ is the reduced Planck mass. Then, Eqs. (12), (14), and (15) reduce, respectively, to

[TABLE]

where Eq. (13) is trivially satisfied. In this case, the longitudinal vector component corresponds to an unphysical gauge mode in that the value of $A_{1}$ is undetermined from Eqs. (18)-(20). From Eqs. (19) and (20) we have that $f^{\prime}/f=h^{\prime}/h$. The solution recovering the asymptotically flat geometry at spatial infinity $r\to\infty$ corresponds to $f=h$ after the proper rescaling of time. Then, it follows that

[TABLE]

which are integrated to give

[TABLE]

where $P,Q,M$ are integration constants. This corresponds to the RN solution with mass $M$ and charge $Q$ of the BH. Note that $P$ is an arbitrary constant without having any physical meaning.

II.2 Massive Proca field

Let us also revisit the massive Proca field in GR given by the functions

[TABLE]

where $m$ is a non-vanishing constant. Since Eq. (13) reduces to $A_{1}r^{2}\,f\,m^{2}=0$, the longitudinal mode is constrained to be

[TABLE]

From Eqs. (12), (14), and (15) it follows that

[TABLE]

Combining Eq. (28) with Eq. (29), we obtain

[TABLE]

On the horizon characterized by the distance $r_{h}$, we have that $f=h=0$. Since the metrics can be expanded as $f=\sum_{i=1}f_{i}(r-r_{h})^{i}$ and $h=\sum_{i=1}h_{i}(r-r_{h})^{i}$ around the horizon, the l.h.s. of Eq. (30) is finite at $r=r_{h}$. This means that $A_{0}$ needs to approach 0 as $r\to r_{h}$ for the consistency with the r.h.s. of Eq. (30). Imposing the asymptotic flatness at spatial infinity, $f\to 1$ and $h\to 1$ as $r\to\infty$, it follows that $A_{0}\to 0$ as $r\to\infty$. If we choose the boundary condition $A_{0}^{\prime}\neq 0$ at $r=r_{h}$, then Eq. (27) gives rise to the growing-mode solution $A_{0}\propto e^{mr}/r$ at spatial infinity. Provided that $A_{0}$ starts to deviate from 0 at some distance, this growing mode manifests for $r\gtrsim 1/m$. Hence the solution consistent with the regularity in two asymptotic regimes is given by

[TABLE]

throughout the horizon exterior Beken72 . Substituting Eq. (31) into Eqs. (28) and (29), the integrated solutions read

[TABLE]

which corresponds to the Schwarzschild geometry.

The reason why we obtained the solution $A_{0}=0$ is attributed to the appearance of terms containing $A_{0}$ in Eqs. (27)-(29). Let us consider the more general case in which $G_{2}$ depends on both $X$ and $F$. From Eq. (13) we have $A_{1}r^{2}fG_{2,X}=0$, so there exists the branch $A_{1}=0$ again for $G_{2,X}\neq 0$. For the massive Proca field discussed above, the terms $-2r^{2}f^{2}A_{0}G_{2,X}$ in Eq. (12), $G_{2}-2X_{0}G_{2,X}$ in Eq. (14), and $-G_{2}+2X_{1}G_{2,X}$ in Eq. (15) give rise to those containing $A_{0}$ in Eqs. (27)-(29). For the theories with $G_{2}=g_{4}(X)h(F)$, where $g_{4}(X)$ and $h(F)$ are functions of $X$ and $F$ respectively, they are factored out by $h(F)$. Since $F=hA_{0}^{\prime 2}/(2f)$, the terms containing $A_{0}$ can be multiplied by the power of the derivative $A_{0}^{\prime}$. Then, the equations corresponding to (27) and (30) admit the solution where $A_{0}$ approaches a non-vanishing constant $P$ at spatial infinity [like the solution (24)].

The above argument shows that the property (31) does not generally hold for the theories with $G_{2}=g_{4}(X)h(F)$. As we will discuss later in Sec. IX, this is actually the case for the coupling $G_{2}=-2g_{4}(X)F$. If we consider the theories with $G_{2}=g(X)+h(F)$, the terms $A_{0}$ appear in the equations of motion without the multiplication of the powers of $A_{0}^{\prime}$. In such cases, $A_{0}$ is generally forced to vanish.

The coupling $G_{2}(X)$ is a rather specific case in which a non-vanishing effective mass term $G_{2,X}A_{0}$ is present in Eq. (12) even at spatial infinity. This fact does not allow the existence of hairy BH solutions. For other derivative interactions $G_{3,4,5,6}$ and $g_{5}$ the terms containing $A_{0}$ are typically multiplied by the powers of $A_{0}^{\prime}$ or by the inverse powers of $r$, so there exists the solution whose asymptotic behavior for $r\to\infty$ is $A_{0}\to P\neq 0$. In such cases, it is possible to realize hairy BH solutions with $A_{0}\neq 0$ outside the horizon. We also note that, in the presence of general derivative interactions, there are branches of solutions where the longitudinal mode $A_{1}$ does not vanish.

III Exact BH solutions for the quartic coupling $G_{4}$

In this section, we first revisit the exact hairy BH solution with $A_{1}\neq 0$ which are known to exist for the theory given by the coupling $G_{4}(X)=M_{\rm pl}^{2}/2+X/4$ Chagoya and then derive another solution corresponding to the branch $A_{1}=0$. The exact BH solutions of Ref. Chagoya obeys the two conditions

[TABLE]

where $X_{c}$ is a constant. These two conditions are imposed to search for exact solutions in this section and Sec. IV. The condition (34) translates to

[TABLE]

where $\epsilon=\pm 1$, and we used Eq. (33).

The longitudinal mode (35) exhibits the divergence at the horizon where $f=0$ for $A_{0}\neq 0$. However, this divergence simply comes from the choice of coordinates. To see this, we introduce the tortoise coordinate $dr_{*}\equiv dr/f(r)$ and consider the scalar product Minami

[TABLE]

Since $A_{1}\simeq\epsilon A_{0}/f$ around the horizon, the product (36) reduces to

[TABLE]

where $u_{+}\equiv t+r_{*}$ and $u_{-}\equiv t-r_{*}$ are the advanced and retarded null coordinates, respectively. The coordinates $u_{+}$ and $u_{-}$ are regular at the future and past event horizons, respectively. Thus, the regularity of the vector field is ensured at the corresponding (future or past) horizon.

For the general quartic coupling $G_{4}(X)$ the vector-field Eqs. (12) and (13) reduce, respectively, to

[TABLE]

From Eq. (39) there are two branches characterized by $A_{1}\neq 0$ and $A_{1}=0$. In the following, we will consider the two cases separately.

III.1 $A_{1}\neq 0$

If the second derivative of $G_{4}$ with respect to $X$ obeys the condition

[TABLE]

then Eq. (39) can be satisfied for

[TABLE]

with $G_{4,X}\neq 0$. Under the condition (33), the solution to Eq. (41) is given by the Schwarzschild metric

[TABLE]

where $M$ is an integration constant. Then, Eq. (38) is satisfied for

[TABLE]

whose integrated solution is

[TABLE]

where $P$ and $Q$ are constants. Now, we search for solutions obeying the condition (34). On using Eqs. (35), (42), (44), and their derivatives with respect to $r$, we find that Eqs. (14)-(16) can be satisfied for

[TABLE]

with the longitudinal vector component

[TABLE]

Since the constant $P$ in Eq. (44) does not depend on $M$ and $Q$, the Proca hair is of the primary type Herdeiro . The function $G_{4}(X)$ obeying the two conditions (40) and (45) is given by

[TABLE]

where $X_{c}=P^{2}/2$, and $b_{n}$’s are constants. The model $G_{4}(X)=M_{\rm pl}^{2}/2+X/4$ of Ref. Chagoya is the special case of Eq. (48), i.e., $G_{4}(X_{c})=M_{\rm pl}^{2}/2+X_{c}/4$ and $b_{n}=0$ for $n\geq 3$. The above solution is a stealth Schwarzschild solution with a non-vanishing longitudinal vector component.

III.2 $A_{1}=0$

Let us proceed to another branch characterized by $A_{1}=0$. Imposing the condition (34), we have $A_{0}^{2}(r)=2f(r)X_{c}$ from Eq. (35). Under the condition (33), Eq. (12) reduces to

[TABLE]

If we consider the case in which the relation $rf^{\prime}+f-1=0$ holds, then the resulting solution $f=1-2M/r$ does not obey Eq. (49). Then, we search for solutions satisfying

[TABLE]

with

[TABLE]

Integration of Eq. (50) leads to

[TABLE]

where $C$ and $M$ are constants. The solution (52) is consistent with Eqs. (14)-(16) for $C=1$ and

[TABLE]

An explicit model satisfying the conditions (51) and (53) is given by

[TABLE]

The resulting exact solution reads

[TABLE]

where $P=\epsilon\sqrt{2X_{c}}$. This corresponds to the extremal RN BH solution.

IV Exact BH solutions for general couplings

We proceed to the derivation of exact BH solutions in the presence of the couplings $G_{3}(X),G_{5}(X),G_{6}(X),g_{5}(X)$ and $G_{2}(X,F)=-2g_{4}(X)F$. Throughout the analysis we take into account the Einstein-Hilbert term $M_{\rm pl}^{2}/2$ in $G_{4}(X)$. Analogous to the derivation of exact solutions given in Sec. III, we will impose the two conditions (33) and (34) in the following discussion.

IV.1 Cubic coupling $G_{3}(X)$

For the cubic interaction $G_{3}(X)$, Eq. (13) reduces to

[TABLE]

Since there are two branches satisfying (i) $G_{3,X}(X_{c})=0$ and (ii) $G_{3,X}(X_{c})\neq 0$, we will discuss each case separately.

IV.1.1 $G_{3,X}(X_{c})=0$

For the branch (i), Eqs. (12) and (14) reduce to Eqs. (21) and (22), respectively, so we obtain the RN solutions (23) and (24). From Eq. (35) the longitudinal mode reduces to

[TABLE]

Since the constant $P$ is independent of $M$ and $Q$, it can be regarded as the primary hair. A concrete example realizing this exact solution is given by

[TABLE]

IV.1.2 $G_{3,X}(X_{c})\neq 0$

The branch (ii) includes the case of vector Galileons ($G_{3}=\beta_{3}X$). On using the conditions (33) and (34) in Eq. (56), it follows that

[TABLE]

Substituting this relation into Eq. (12), we obtain $rA_{0}^{\prime\prime}+2A_{0}^{\prime}=0$. Hence the integrated solution is $A_{0}=P+Q/r$ with two constants $P$ and $Q$. Then, Eq. (59) is integrated to give

[TABLE]

where $C$ is a constant. To satisfy the asymptotically flat boundary condition $f\to 1$ as $r\to\infty$, we require that $P^{2}=2X_{c}$. The above solutions are consistent with Eqs. (14)-(16) for $C=0$ and $X_{c}=M_{\rm pl}^{2}$. On defining $M=\pm Q/(\sqrt{2}M_{\rm pl})$ for $P=\mp\sqrt{2}M_{\rm pl}$, we obtain the extremal BH solution (55) with $P=\epsilon\sqrt{2}M_{\rm pl}$. The longitudinal mode $A_{1}$ vanishes for this exact solution.

IV.2 Quintic coupling $G_{5}(X)$

For the quintic interaction $G_{5}(X)$, combining Eq. (12) with Eq. (13) leads to

[TABLE]

The solution to Eq. (61) is given by $A_{0}=P+Q/r$. If

[TABLE]

then Eq. (62) is satisfied either for (i) $A_{0}A_{0}^{\prime}=X_{c}f^{\prime}$ or (ii) $A_{0}^{2}=2fX_{c}$.

For the branch (i) we have

[TABLE]

which is integrated to give $f=C-2M/r+Q^{2}/(2X_{c}r^{2})$ with $M=-PQ/(2M_{\rm pl}^{2})$. For the consistency with Eqs. (14)-(15) we require that $C=1$ and $X_{c}=M_{\rm pl}^{2}$, so we obtain the RN solution

[TABLE]

with the vector components

[TABLE]

The existence of this solution requires the condition $2M^{2}M_{\rm pl}^{2}>Q^{2}$. Since the Proca hair $P=-2MM_{\rm pl}^{2}/Q$ is fixed by $M$ and $Q$, it is of the secondary type.

From Eq. (35) the branch (ii) corresponds to $A_{1}=0$. In this case, the integration of Eqs. (14)-(15) gives rise to the RN solutions (23) and (24). On using the property $A_{0}^{2}=2fX_{c}$, the metric $f$ reduces to the extremal RN solution $f=(1-M/r)^{2}$ with the particular relation $Q^{2}=2M^{2}M_{\rm pl}^{2}$. Indeed, this case can be regarded as the special case of the solutions (66) with $A_{1}=0$.

A concrete mode realizing the above solutions is given by

[TABLE]

where $X_{c}=M_{\rm pl}^{2}$.

IV.3 Sixth-order coupling $G_{6}(X)$

In the presence of the sixth-order coupling $G_{6}(X)$, Eq. (13) reduces to

[TABLE]

Let us search for exact solutions satisfying either $A_{0}^{\prime}=0$ or $A_{1}=0$.

IV.3.1 $A_{0}^{\prime}=0$

In this case we have

[TABLE]

under which Eq. (12) is trivially satisfied. From Eqs. (14) and (15) we obtain $rf^{\prime}+f-1=0$, so the integrated solution is given by the Schwarzschild metric

[TABLE]

In fact, this solution exists for general couplings $G_{6}(X)$ with any value of $A_{1}$. In the present case the longitudinal mode is subject to the constraint (35), so it is given by

[TABLE]

Since $A_{1}$ approaches the constant $\epsilon\sqrt{P^{2}-2X_{c}}$ as $r\to\infty$, we require the condition $P^{2}>2X_{c}$ for the existence of this solution.

IV.3.2 $A_{1}=0$

We proceed to the case in which the longitudinal mode obeys

[TABLE]

In this case we have $A_{0}^{2}(r)=2f(r)X_{c}$, so we take the $r$-derivative of this relation and substitute them into Eqs. (12)-(14). Then, Eqs. (16) and (12) reduce, respectively, to

[TABLE]

Let us search for exact solutions satisfying the two conditions

[TABLE]

From Eqs. (73)-(74) we obtain the integrated solution $f=(C-M/r)^{2}$ with $X_{c}=M_{\rm pl}^{2}$. Since the integration constant is fixed to be $C=1$ from Eq. (14), we obtain the extremal RN solution (55) with $P=\epsilon\sqrt{2}M_{\rm pl}$. This is equivalent to the solution derived for the quintic coupling $G_{5}(X)$ with the branch $A_{0}^{2}=2fX_{c}$.

A concrete model realizing this solution is

[TABLE]

where $X_{c}=M_{\rm pl}^{2}$.

IV.4 Quartic intrinsic vector-mode coupling $g_{4}(X)$

Let us consider the coupling given by

[TABLE]

where $g_{4}(X)$ is a function of $X$. This corresponds to the intrinsic vector mode originally introduced in ${\cal L}_{4}$ as a form $g_{4}(X)(\nabla_{\rho}A_{\sigma}\nabla^{\rho}A^{\sigma}-\nabla_{\rho}A_{\sigma}\nabla^{\sigma}A^{\rho})$ with $g_{4}(X)=c_{2}G_{4,X}$ Heisenberg . Then, Eq. (13) reduces to

[TABLE]

Let us consider the case in which the relation

[TABLE]

is satisfied. From Eqs. (12) and (14) it follows that

[TABLE]

For $g_{4}(X_{c})\neq 1/2$ these equations are integrated to give

[TABLE]

with the longitudinal mode (35). The metric is of the RN type with the effective charge $Q_{\rm eff}=\sqrt{1-2g_{4}(X_{c})}Q$, which is different from $Q$ unless $g_{4}(X_{c})=0$.

A concrete model realizing this solution is given by

[TABLE]

If $g_{4}(X_{c})=1/2$, then the Schwarzschild solution $f=h=1-2M/r$ follows from Eq. (81) with $A_{0}$ undetermined. This comes from the fact that, for $g_{4}(X_{c})=1/2$, the Lagrangian $F$ is compensated by the term $-2g_{4}F$.

It is also possible to satisfy Eq. (78) either for (i) $A_{0}^{\prime}=0$ or (ii) $A_{1}=0$. For the branch (i), we obtain the Schwarzschild solution $f=h=1-2M/r$ with $A_{0}={\rm constant}$ and $A_{1}$ given by Eq. (35) for general couplings $g_{4}(X)$. For the branch (ii), there exists an exact solution under the conditions $g_{4,X}(X_{c})=0$ and $g_{4}(X_{c})=1/2$. In this case, the resulting solution reads

[TABLE]

This solution exists for the function (84) with $g_{4}(X_{c})=1/2$.

IV.5 Quintic intrinsic vector-mode coupling $g_{5}(X)$

Let us finally proceed to the exact solution for the quintic coupling $g_{5}(X)$. Then, Eq. (13) reduces to

[TABLE]

For the branch $A_{0}^{\prime}=0$, the Schwarzschild solution $f=h=1-2M/r$ follows with $A_{1}$ given by Eq. (35) for general couplings $g_{5}(X)$.

For the other branch $fg_{5}=(A_{0}^{2}-2fX_{c})g_{5,X}$, there exists an exact solution under the condition

[TABLE]

From Eqs. (12) and (14) we obtain the equations same as Eqs. (21) and (22), respectively, so the integrated solutions to $f,h,A_{0}$ yield the RN solutions (23) and (24) with $A_{1}$ given by Eq. (35). Since $g_{5}(X_{c})=0$ in this case, the quintic interaction in the form

[TABLE]

gives rise to the RN solution with the non-vanishing longitudinal mode.

V Power-law cubic couplings $G_{3}(X)$

In Secs. III and IV, we have imposed the two conditions (33) and (34) for the purpose of deriving exact BH solutions. We now focus on the solutions where $f\neq h$ and $X$ is not constant. Numerical works are generally required to find such non-exact solutions. In this section, we first study the model in which the function $G_{3}$ is given by the power-law function $X^{n}$, where $n$ is assumed to be a positive integer. In the subsequent sections, we will study the models in which the functions $G_{4},G_{5},G_{6},g_{4},g_{5}$ contain the function $X^{n}$. In the whole analysis by the end of Sec. IX, we include the Einstein-Hilbert term $M_{\rm pl}^{2}/2$ in $G_{4}$. We will focus on the asymptotically flat solutions and not take into account the vector-field mass and the cosmological constant.

We begin with the power-law cubic coupling model given by

[TABLE]

with $G_{4}=M_{\rm pl}^{2}/2$, where $\beta_{3}$ is a dimensionless constant. In the following, we will discuss the cases of $n=1$ (vector Galileons) and $n\geq 2$, separately.

V.1 $n=1$

From (13) the longitudinal component for $n=1$ is related to $A_{0},f,h$, as

[TABLE]

We substitute Eq. (90) and the $r$-derivative of it into Eqs. (12), (14), and (15) to eliminate the $A_{1}$ dependence.

Around the BH horizon characterized by the distance $r_{h}$, we expand $f,h,A_{0}$ in the following forms

[TABLE]

where $f_{i},h_{i},a_{0},a_{i}$ are constants. The effect of the coupling $\beta_{3}$ works as corrections to the RN metrics given by

[TABLE]

The constant $\mu$ is in the range $0<\mu<1$, so that $r_{h}$ corresponds to the outer horizon. Compared to Eq. (23), there is the correspondence $Q^{2}=2r_{h}(2M-r_{h})M_{\rm pl}^{2}$ with the inner horizon $\tilde{r}_{h}=2M-r_{h}$. Hence the constant $\mu$ is given by $\mu=2M/r_{h}-1$ with $M<r_{h}<2M$. To derive the coefficients $f_{i},h_{i},a_{0},a_{i}$ in Eq. (91) iteratively, we assume that $f_{1},h_{1},a_{0}$ are positive and choose the positive branch of Eq. (90) for $r>r_{h}$. We also take the contributions up to linear order in $\beta_{3}$ under the assumption that the coupling $\beta_{3}$ works as a correction to the RN solutions. Up to the order of $(r-r_{h})^{2}$, the coefficients are given by

[TABLE]

and

[TABLE]

where

[TABLE]

From Eq. (93) and the condition $0<\mu<1$, the quantity $h_{1}r_{h}$ is in the range $0<h_{1}r_{h}<1$. The corrections to the RN solutions from the coupling $\beta_{3}$ arise at the order of $(r-r_{h})^{2}$ for $f,h,A_{0}$.

Taking the positive branch of Eq. (90), the behavior of the longitudinal mode around the horizon is given by

[TABLE]

which exhibits the divergence at $r=r_{h}$. Analogous to the discussion given after Eq. (35), the scalar product (36) reduces to $A_{\mu}dx^{\mu}\simeq a_{0}du_{+}$ around $r=r_{h}$ for the solution (96). Hence the regularity of the vector field is ensured at the future horizon.

We also derive asymptotic flat solutions satisfying $f,h\to 1$, and $A_{0}\to P$ as $r\to\infty$, where $P$ is a constant. To obtain the solutions at spatial infinity, we expand $f,h,A_{0}$ as the power series of $1/r$, as

[TABLE]

For the cubic coupling model (89) with $n=1$, there exists an asymptotic solution where the longitudinal mode is given by $A_{1}=\sum_{i=1}^{\infty}\tilde{b}_{i}/r^{i}$. Substituting this expression of $A_{1}$ and Eqs. (97) into Eqs. (12)-(16), we obtain the following iterative solutions

[TABLE]

where we have set $\tilde{f}_{1}=\tilde{h}_{1}=-2M$. The metric $f$ does not contain the term proportional to $P^{2}M^{2}/(M_{\rm pl}^{2}r^{2})$ unlike the metric $h$, so there is the difference between $f$ and $h$ for $P\neq 0$ at the order of $1/r^{2}$. The leading-order solutions to temporal and longitudinal vector components are given, respectively, by $A_{0}\simeq P(1-M/r)$ and $A_{1}\simeq\tilde{b}_{2}/r^{2}$, which are also consistent with the solutions derived under the weak gravity approximation far outside a spherically symmetric body DeFeliceVain . The effects of the constants $\tilde{b}_{2}$ and $\beta_{3}$ start to appear in the metrics $f$ and $h$ at the orders of $1/r^{4}$ and $1/r^{5}$, respectively.

To study whether the near-horizon solutions (91) connect to the large-distance solutions (97), we numerically solve Eqs. (12)-(16) outside the horizon by using the boundary conditions (93)-(94) around $r=r_{h}$. In Fig. 1 we plot the numerically integrated solutions of $f,h,A_{0},A_{1},h-f$ for the coupling $\beta_{3}=7\times 10^{-3}/(r_{h}M_{\rm pl})$. The two asymptotic solutions smoothly join each other without any singular behavior. While the metric $h$ approaches 1 as $r\to\infty$, the existence of the coupling $\beta_{3}$ in Eq. (94) around the horizon leads to the value of $f$ slightly different from 1 in the large-distance limit. We exploit the freedom of performing a time rescaling to shift $f$ to 1 by solving the equations of motion up to $r=10^{7}r_{h}$. After this rescaling, the difference between $f$ and $h$ is of the order of $10^{-2}$ around the horizon for the coupling $\beta_{3}$ chosen in Fig. 1. For increasing $|\beta_{3}|$, this difference tends to be larger. The deviation from GR is most significant in the regime of strong gravity and it starts to decrease for $r\gtrsim 10r_{h}$. This signature can be potentially probed by the future measurements of gravitational waves in the nonlinear regime of gravity.

From Fig. 1 we find that the longitudinal vector component behaves as $A_{1}\propto(r-r_{h})^{-1}$ around the horizon. As explained already, the apparent divergence of $A_{1}$ at $r=r_{h}$ does not spoil the regularity of the vector field. In the limit $r\to\infty$, we also numerically confirm that $A_{1}$ decreases in proportion to $r^{-2}$. The temporal vector component approaches the constant $P$ as $r\to\infty$.

Due to the continuity of two asymptotic solutions outside the horizon, the model parameters $P,M,\tilde{b}_{2}$ in Eqs. (98)-(101) should be related to the parameters $\mu,r_{h},a_{0}$ in Eqs. (93)-(94), as $P=P(\mu,r_{h},a_{0})$, $M=M(\mu,r_{h},a_{0})$, and $\tilde{b}_{2}=\tilde{b}_{2}(\mu,r_{h},a_{0})$. Since the constant $P$ cannot be fixed by other two parameters $M$ and $\tilde{b}_{2}$, it can be regarded as the primary Proca hair.

V.2 $n\geq 2$

Let us proceed to the cubic coupling (89) with the powers $n\geq 2$. In this case, Eq. (13) reduces to

[TABLE]

Hence, there are two branches characterized by (i) $fh(rf^{\prime}+4f)A_{1}^{2}+rA_{0}(2fA_{0}^{\prime}-f^{\prime}A_{0})=0$ or (ii) $A_{0}^{2}-fhA_{1}^{2}=0$, which we will discuss separately below.

V.2.1 Branch (i)

For this branch, the longitudinal mode $A_{1}$ satisfies the relation same as Eq. (90). Taking the positive branch of $A_{1}$ for $r>r_{h}$ and expanding $f,h,A_{0}$ around the horizon as Eq. (91), we obtain the coefficients $f_{1},h_{1},a_{1},f_{2},h_{2},a_{2}$ whose forms are the same as Eqs. (93)-(94) but with different values of ${\cal F}_{2},{\cal H}_{2},\alpha_{2}$. For $n=2$, they are given by

[TABLE]

where ${\cal F}_{2}^{(n=1)}$, ${\cal H}_{2}^{(n=1)}$, and $\alpha_{2}^{(n=1)}$ correspond to the values of ${\cal F}_{2},{\cal H}_{2},\alpha_{2}$ for $n=1$, and

[TABLE]

The effect of the coupling $\beta_{3}$ works in a similar way to that discussed for $n=1$. For larger $\beta_{3}$, the difference between the metrics $f$ and $h$ tends to be more significant around the horizon. The longitudinal mode around the horizon is given by Eq. (96), so it behaves as $A_{1}\propto(r-r_{h})^{-1}$.

Far outside the horizon ($r\gg r_{h}$) we can also derive iterative solutions by using the expansion (97). Up to the order $1/r^{4}$, the solutions to $f,h,A_{0}$ for $n\geq 2$ are exactly the same as Eqs. (98), (99), and (100), respectively. At leading order, the longitudinal component decreases as $A_{1}=\tilde{b}_{2}/r^{2}$. Thus, the behavior of large-distance solutions is similar to that for $n=1$. The asymptotic solutions of $f,h,A_{0},A_{1}$ for $r\gg r_{h}$ contain the three parameters $P,M,\tilde{b}_{2}$, which are related to other three parameters $\mu,r_{h},a_{0}$ appearing in the solutions expanded around the horizon. Since $P$ is not solely fixed by $M$ and $\tilde{b}_{2}$, the vector hair is of the primary type.

V.2.2 Branch (ii)

The branch (ii) satisfies the relation

[TABLE]

Substituting Eq. (105) into Eqs. (12), (14), and (15), we obtain the differential equations same as Eqs. (18), (19), and (20), respectively. Imposing the boundary conditions $f=h=1$ at $r\to\infty$, we obtain the RN metrics (23) with the temporal vector component (24). Indeed, the relation (105) corresponds to the special case of Eq. (35) with $X_{c}=0$ and $f=h$. As we discussed in Sec. IV.1.1, the RN solutions (23) and (24) with the longitudinal mode (57) exist for the theory given by the function (58) under the condition $G_{3,X}(X_{c})=0$. Since we are now considering the coupling $G_{3}(X)\propto X^{n}$ with $n\geq 2$, the condition $G_{3,X}(X_{c})=0$ is satisfied for $X_{c}=0$. Setting $X_{c}=0$ in Eq. (57), it follows that $A_{1}$ is equivalent to Eq. (105) with $f,h,A_{0}$ given by Eqs. (23) and (24).

VI Power-law quartic couplings $G_{4}(X)$

Let us proceed to the model of quartic power-law interactions given by

[TABLE]

where $\beta_{4}$ is a dimensionless constant. From Eq. (13) the longitudinal mode obeys

[TABLE]

For $n=1$, Eq. (107) reduces to

[TABLE]

so we have two branches satisfying (i) $hrf^{\prime}+(h-1)f=0$ or (ii) $A_{1}=0$. As we showed in Sec. III, there exists an exact BH solution for the branch (i) with $\beta_{4}=1/4$. In Ref. Chagoya2 the solutions for general $\beta_{4}$ were discussed for the two branches (i) and (ii), so we will not repeat the analysis here.

In the following we will study the $n=2$ (vector Galileons) and $n\geq 3$ cases, separately.

VI.1 $n=2$

In this case, Eq. (107) yields

[TABLE]

Since there are two branches characterized by $A_{1}\neq 0$ or $A_{1}=0$, we will discuss such two cases in turn.

VI.1.1 Branch with $A_{1}\neq 0$

This branch corresponds to the longitudinal mode satisfying

[TABLE]

We first derive the solutions around the horizon by expanding the functions $f,h,A_{0}$ in the forms (91) with the constraint (110). We take the positive branch of Eq. (110) for $r>r_{h}$, assume that $f_{1}>0$ and $h_{1}>0$, and pick up terms linear in $\beta_{4}$. The resulting coefficients up to the order of $(r-r_{h})^{2}$, which recover the RN metrics (92) in the limit $\beta_{4}\to 0$, read

[TABLE]

and

[TABLE]

where

[TABLE]

Apart from the appearance of the $\beta_{4}$-dependent term in $a_{1}$, the structure of solutions around the horizon is similar to that of the power-law cubic models studied in Sec. V. The coupling $\beta_{4}$ works as corrections to the leading-order RN solutions characterized by the first terms on the r.h.s. of Eq. (112). Taking the positive branch, the behavior of the longitudinal mode (110) is given by

[TABLE]

The leading-order contribution to $A_{1}$ is similar to that in Eq. (96) of the cubic-coupling case, so the regularity of solutions is ensured around the horizon.

For the distance $r$ much lager than $r_{h}$, we perform the expansions of $f,h,A_{0}$ given by Eq. (97). In doing so, we take the $r$ derivative of Eq. (110) and eliminate the terms $A_{1}$ and $A_{1}^{\prime}$ from Eqs. (12), (14), and (15). Picking up the leading-order terms of $\beta_{4}$, we obtain the iterative solutions

[TABLE]

where we have set $f_{1}=h_{1}=-2M$ and $\tilde{a}_{1}=Q$. We recover the RN solutions (23)-(24) by taking the limit $\beta_{4}\to 0$ in Eqs. (118)-(120). The existence of the coupling $\beta_{4}$ leads to the difference between the two metric components $f$ and $h$ at the order of $1/r^{2}$. The leading-order longitudinal mode decreases as $A_{1}\propto 1/\sqrt{r}$, whose property is different from that in the cubic power-law models (in which case $A_{1}\propto r^{-2}$).

In Fig. 2 we show the numerically integrated solutions to $f,h,A_{0},A_{1},f-h$ for $\beta_{4}=10^{-2}$ derived by choosing the boundary conditions (91) with Eqs. (111)-(112) near the horizon. The solutions smoothly connect to those in another asymptotic regime $r\gg r_{h}$. As in the numerical simulation of Fig. 1, the metric $f$ is rescaled to 1 at a sufficiently large distance ($r=10^{7}r_{h}$). For $\beta_{4}=10^{-2}$ the maximum difference between $f$ and $h$ is of the order of $10^{-3}$ around the horizon, but it decreases for larger $r$ according to Eqs. (118)-(119). Numerically, we also confirm that the longitudinal mode behaves as $A_{1}\propto 1/(r-r_{h})$ around $r\simeq r_{h}$ and $A_{1}\propto 1/\sqrt{r}$ for $r\gg r_{h}$. The parameters $\mu,r_{h},a_{0}$ around the horizon are related to $P,Q,M$ in Eqs. (118)-(121). Since the constant $P$ can not be fixed by $Q$ and $M$, it is regarded as a primary hair.

VI.1.2 Branch with $A_{1}=0$

From Eq. (109) there exists the other branch satisfying $A_{1}=0$. Substituting $A_{1}=0$ and $A_{1}^{\prime}=0$ into Eqs. (12)-(15) and expanding the functions $f,h,A_{0}$ as Eq. (91) around the horizon, the solutions recovering the RN metric (92) in the limit $\beta_{4}\to 0$ are given by

[TABLE]

where

[TABLE]

and we have chosen the branch $a_{1}>0$. The temporal component $A_{0}$ exactly vanishes at $r=r_{h}$. The effect of $\beta_{4}$ on $f,h,A_{0}$ arises at the order of $(\Delta x)^{3}$, which is higher than the order of the branch $A_{1}\neq 0$.

For the distance $r\gg r_{h}$, the expansions of $f,h,A_{0}$ in the forms (97) give rise to the following iterative solutions (up to linear order in $\beta_{4}$):

[TABLE]

where we have set $\tilde{a}_{1}=Q$ and $\tilde{h}_{1}=-2M$. The coupling $\beta_{4}$ works as corrections to the RN solutions (23)-(24), which induce the difference between $f$ and $h$ even at the order of $1/r$. If we consider the case in which the conditions $|Q|\ll|MP|$ and $|\beta_{4}|\ll 1$ are satisfied, the horizon radius $r_{h}$ is close to $2M$. For $|P|$ of the order of $M_{\rm pl}$, Eqs. (VI.1.2) and (127) show that $|f-h|$ is of the order of $(r_{h}/r)\beta_{4}$ for $r\gg r_{h}$. This means that, even if the corrections induced by $\beta_{4}$ are small in the very vicinity of the horizon, they are not necessarily negligible for $r$ larger than the order of $r_{h}$ (say, $r=10r_{h}$). At spatial infinity the effects of $\beta_{4}$ on $f$ and $h$ are suppressed.

Besides the coupling $\beta_{4}$, we have two parameters $\mu$ and $r_{h}$ in Eqs. (122)-(124), whereas there are three parameters $P,M,Q$ in Eqs. (VI.1.2)-(128). Numerically we confirmed that the two asymptotic solutions in the regimes $r\simeq r_{h}$ and $r\gg r_{h}$ smoothly join each other. Hence $P,Q,M$ are related to $\mu,r_{h}$ according to $P=P(\mu,r_{h})$, $Q=Q(\mu,r_{h})$, and $M=M(\mu,r_{h})$. In this case $P$ depends on $Q$ and $M$, so the vector hair is of the secondary type.

VI.2 $n\geq 3$

For the theories with $n\geq 3$, there are the three branches of solutions:

[TABLE]

The first two branches are similar to those discussed for $n=2$, but the branch (iii) arises only for $n\geq 3$.

For the branch (i), the solution around the horizon is given by Eq. (91) with the coefficients in the forms (111) and (112) but with different values of $\alpha_{1},{\cal F}_{2},{\cal H}_{2},\alpha_{2}$. As in the case of $n=2$, the effect of the coupling $\beta_{4}$ induces the difference between $f$ and $h$ at the order of $(r-r_{h})^{2}$. At spatial infinity the modifications to the RN solutions arise at the order of $1/r^{2}$, whose property is similar to Eqs. (118)-(120) of the $n=2$ case. If $n=3$, for example, the leading-order corrections to $f,h,A_{0}$ of the RN solutions (23)-(24) are given, respectively, by

[TABLE]

The longitudinal mode decreases as $A_{1}\propto 1/\sqrt{r}$ for $r\gg r_{h}$. Our numerical simulations show that the solutions in two asymptotic regimes $r\approx r_{h}$ and $r\gg r_{h}$ smoothly join each other.

For the branch (ii), the effect of the coupling $\beta_{4}$ on the solution (91) expanded around the horizon appears at higher orders with increasing $n$, e.g., at the order of $(r-r_{h})^{4}$ for $n=3$. At spatial infinity, the coupling $\beta_{4}$ leads to the difference between $f$ and $h$ at the order of $1/r$. If $n=3$, for example, the two metric components are given by

[TABLE]

so that $|f-h|$ can be of the order of $(r_{h}/r)\beta_{4}$ for $|Q|\ll|MP|$ and $|P|={\cal O}(M_{\rm pl})$.

For the branch (iii), it follows that Eqs. (12), (14), and (15) reduce to Eqs. (18), (19), and (20), respectively. Hence this branch corresponds to the RN solutions with $A_{1}=\epsilon A_{0}/f$. Since $X=0$ in this case, we have $G_{4,X}=0$ and $G_{4,XX}=0$ for $n\geq 3$. This is the reason why all the $X$-dependent terms arising from $G_{4}(X)$ vanish from the equations of motion. Note that the exact solution discussed in Sec. III is different from the above non-exact solution, because the former satisfies the condition $G_{4,X}=1/4$.

VII Power-law quintic couplings $G_{5}(X)$

In this section, we consider power-law quintic couplings given by

[TABLE]

with $G_{4}=M_{\rm pl}^{2}/2$, where $\beta_{5}$ is a dimensionless constant. From Eq. (13) it follows that

[TABLE]

If $n=1$, then this reduces to

[TABLE]

so there are two physically equivalent branches corresponding to $\epsilon=\pm 1$.

For $n\geq 2$, it follows that

[TABLE]

where

[TABLE]

We require the two conditions $\xi_{1}^{2}\geq\xi_{2}$ and $\xi_{1}[1\pm\sqrt{1-\xi_{2}/\xi_{1}^{2}}]\geq 0$ for the existence of the solutions (137).

For $n\geq 3$, Eq. (135) admits the following solutions

[TABLE]

besides the solutions satisfying Eq. (137). As in the cases of cubic and quartic power-law couplings, the branches (139) correspond to the RN solutions (23) and (24).

In what follows, we will focus on the branches (136) and (137) for $n=1$ and $n\geq 2$, respectively. At spatial infinity, we expand $f,h,A_{0}$ in the forms (97) and also assume the asymptotic behavior $A_{1}=\sum_{i=1}^{\infty}\tilde{b}_{i}/r^{i}$. Then, we obtain the iterative solutions

[TABLE]

where $\mu_{P}=P^{2}-2M_{\rm pl}^{2}$, and we set $\tilde{a}_{0}=P>0$ and $\tilde{f}_{1}=\tilde{h}_{1}=-2M$. The existence of the above solutions requires the condition $\mu_{P}\geq 0$, i.e., $P^{2}\geq 2M_{\rm pl}^{2}$. On using the large-distance solutions (140)-(142), we find that the quantities $\xi_{1}$ and $\xi_{2}$ defined by Eq. (138) behave as

[TABLE]

For the branch with the positive sign inside the square root of Eq. (137), we have that $A_{1}^{2}\simeq 2\xi_{1}\simeq P^{2}/n={\rm constant}$ at spatial infinity. The branch with the negative sign inside the square root of Eq. (137) gives rise to the solution $A_{1}^{2}\simeq\xi_{2}/(2\xi_{1})\simeq M^{2}P^{2}\mu_{P}/(2M_{\rm pl}^{2}r^{2})$, so this corresponds to the leading-order solution of Eq. (143).

Outside the horizon the metric $h$ is in the range $0<h<1$. For $n=1$, Eq. (136) shows that $A_{1}^{2}$ exhibits the divergence at

[TABLE]

For the $\epsilon=-1$ branch of Eq. (136) we can derive the solutions in the form (91) expanded around the horizon. On using such analytic boundary conditions and solving the equations of motion numerically, we find that $A_{1}$ indeed diverges as $h$ approaches $1/3$. Hence the solutions in the strong-gravity regime ($h<1/3$) are disconnected to the large-distance solutions (140)-(143).

From Eqs. (137)-(138) we find that the similar divergence of $A_{1}$ occurs at

[TABLE]

for $n\geq 2$. Since the metric (146) is in the range $0<h\leq 1/5$, the solutions outside the horizon cannot avoid passing through this divergent point. Our numerical simulations show that there are no regular exterior BH solutions that smoothly connect to Eqs. (140)-(143).

VIII Power-law sixth order couplings $G_{6}(X)$

We proceed to the case of power-law sixth-order interactions given by

[TABLE]

with $G_{4}=M_{\rm pl}^{2}/2$, where $\beta_{6}$ is a dimensionless constant. The $U(1)$-invariant gravitational coupling advocated by Horndeski Horndeski76 corresponds to $n=0$, so we will also include such a case in the analysis. From Eq. (13) the longitudinal mode obeys

[TABLE]

where

[TABLE]

Equation (148) admits the solution $A_{0}^{\prime}=0$, but this corresponds to the stealth Schwarzschild solution (70). This is analogous to the discussion given in Sec. IV.3.1, but the difference is that $A_{1}$ is arbitrary in the present case (since we are not imposing the condition that $X$ is constant). There exist other two branches satisfying ${\cal G}=0$ or $A_{1}=0$. Let us first discuss the possibility for the realization of the branch ${\cal G}=0$.

For $n=0$ the branch ${\cal G}=0$ is realized for $A_{1}^{2}/A_{0}^{2}=(3h-1)/[fh(h-1)]$, but the real solutions to $A_{1}$ do not exist for $1/3<h<1$. When $n=1$ the quantity (149) simply reduces to ${\cal G}=1-3h$, so there is no consistent branch satisfying ${\cal G}=0$ in the whole region outside the horizon. For $n=2$ we have ${\cal G}=A_{1}^{2}fh(5h-1)-A_{0}^{2}(3h-1)$, so the real solutions to ${\cal G}=0$ are not present for $1/5<h<1/3$. For $n\geq 3$ there are two solutions to ${\cal G}=0$, i.e.,

[TABLE]

The branch (i) does not exist in the region $1/(2n+1)<h<1/3$ outside the horizon. For the branch (ii) the solutions are described by the RN solutions (23)-(24). Therefore, apart from the trivial branch (ii) present for $n\geq 3$, there are no consistent solutions satisfying ${\cal G}=0$.

Since the remaining possibility is the branch

[TABLE]

we will focus on this case in the following discussion. Substituting $A_{1}=0$ and $A_{1}^{\prime}=0$ into Eqs. (12), (14), (15) and expanding $f,h,A_{0}$ in the forms (97) at large distances ($r\gg r_{h}$), the iterative solutions are given by

[TABLE]

For $n=0$ these results match with those derived by Horndeski in Ref. HorndeskiBH . The coupling $\beta_{6}$ works as corrections to the leading-order RN solutions. The difference between $f$ and $h$ arises at the order of $1/r^{4}$.

We expand the solutions around the horizon as Eq. (91) and pick up the terms linear in $\beta_{6}$. The resulting solutions, which recover the RN metrics in the limit $\beta_{6}\to 0$, are given by

[TABLE]

for $n=0$, and

[TABLE]

for $n=1$. If $n\geq 2$, the effect of the coupling $\beta_{6}$ arises at the order of $(\Delta x)^{n+1}$ in $f,h,A_{0}$.

For $n=0$ there exists the $U(1)$ gauge symmetry, so the constant $P$ in Eq. (154) has no physical meaning with the value of $a_{0}$ unconstrained. In this case, we have two physical hairs $M$ and $Q$ related to the parameters $\mu$ and $r_{h}$ around the horizon.

For $n\geq 1$ we have that $a_{0}=0$, so the parameters $M,Q,P$ are related to the two parameters $\mu,r_{h}$ appearing for the solutions around the horizon. Then, the Proca hair is of the secondary type. This situation is analogous to what happens for the quartic power-law interactions with the branch $A_{1}=0$.

In Fig. 3 we plot the numerically integrated solutions to $f,h,A_{0},h-f$ for $n=1$ with the branch $A_{1}=0$. On choosing the boundary conditions (158)-(160) around $r=r_{h}$, the solutions smoothly connect to those in the regime $r\gg r_{h}$, i.e., Eqs. (152)-(154) with $n=1$. Since $a_{0}=0$ in this case, the temporal component $A_{0}$ vanishes on the horizon. By normalizing the metric $f$ to be 1 at $r\to\infty$, the coupling $\beta_{6}$ induces the difference between $f$ and $h$ around the horizon. Compared to the cases of cubic and quartic couplings plotted in Figs. 1 and 2, $|f-h|$ decreases faster for increasing $r$. Thus, the future precise measurements for the deviation from GR in the strong-gravity regime may allow us to distinguish between hairy solutions with different couplings. For $n=0$ we have also confirmed that the numerical solutions are regular outside the horizon with the difference between $f$ and $h$. For $n\geq 2$, the effect of the coupling $\beta_{6}$ arises at higher order in metrics around the horizon.

IX Power-law intrinsic vector mode couplings

$g_{4}(X)$ and $g_{5}(X)$

Let us finally study the models of power-law intrinsic vector-mode couplings given by

[TABLE]

where $g_{4}(X)$ is given in Eq. (77), with $G_{4}=M_{\rm pl}^{2}/2$, where $\gamma_{4}$ and $\gamma_{5}$ are dimensionless constants.

IX.1 $\gamma_{4}\neq 0$ and $\gamma_{5}=0$

In this case, there is the relation (78) with $g_{4,X}=n\gamma_{4}(A_{0}^{2}-fhA_{1}^{2})^{n-1}/[(2f)^{n-1}M_{\rm pl}^{2n}]$. The branch satisfying $A_{0}^{\prime 2}=0$ corresponds to the stealth Schwarzschild BH solution (70). For $n\geq 2$ there exists the branch $A_{1}=\epsilon\sqrt{A_{0}^{2}/(fh)}$, in which case the solutions are described by the RN solutions (23)-(24).

In what follows we will focus on the last branch of Eq. (78) with the vanishing longitudinal mode ($A_{1}=0$). Expanding the functions $f,h,A_{0}$ in the forms (97), the resulting large-distance solutions are given by

[TABLE]

For $n\geq 1$, the coupling $\gamma_{4}$ induces the difference between the two metrics $f$ and $h$ at the order of $1/r^{2}$.

For $n=1$, the solutions expanded around $r=r_{h}$, which recover the RN metrics in the limit $\gamma_{4}\to 0$, are given by

[TABLE]

For $n\geq 2$, the effect of the coupling $\gamma_{4}$ appears at the order of $(\Delta x)^{n+1}$ in the expansions of $f,h,A_{0}$. Since $a_{0}=0$ for $n\geq 1$, there are two parameters $\mu$ and $r_{h}$ around the horizon. Since these two parameters are related to $P,Q,M$ in Eqs. (162)-(164), the Proca hair $P$ is of the secondary type.

Numerically, we have confirmed that the two asymptotic solutions given above smoothly connect to each other. The coupling $\gamma_{4}$ gives rise to the difference between $f$ and $h$ in the strong-gravity regime, whose effect tends to be smaller for larger $n$.

IX.2 $\gamma_{5}\neq 0$ and $\gamma_{4}=0$

For the quintic intrinsic vector-mode interaction, Eq. (13) reduces to

[TABLE]

Again, the branch $A_{0}^{\prime 2}=0$ corresponds to the stealth Schwarzschild solution (70). For $n\geq 2$ there exists the branch $A_{1}=\epsilon\sqrt{A_{0}^{2}/(fh)}$, which leads to the RN solutions (23)-(24).

For $n\geq 1$ we have the last branch satisfying

[TABLE]

where we will take the $\epsilon=+1$ branch in the following. Differentiating Eq. (169) with respect to $r$ and substituting $A_{1}$ and $A_{1}^{\prime}$ into Eqs. (12), (14), and (15), the iterative solutions at spatial infinity are given by

[TABLE]

The longitudinal mode behaves as

[TABLE]

which approaches the constant $P/\sqrt{2n+1}$ as $r\to\infty$. The coupling $\gamma_{5}$ induces the difference between $f$ and $h$ at the order of $1/r^{3}$.

Expanding the solutions around $r=r_{h}$ as Eq. (91) for $n\geq 1$, it follows that $a_{0}=0$. The effect of the coupling $\gamma_{5}$ arises at the order of $(r-r_{h})^{n+2}$ in the expansions of $f,h,A_{0}$ as corrections to the leading-order RN solutions. When $n=1$, for example, the resulting solutions are given by

[TABLE]

There are two parameters $\mu$ and $r_{h}$ in the expansions (174)-(176), which are related to the three parameters $P,Q,M$ in Eqs. (170)-(172). Hence the Proca hair $P$ is of the secondary type.

In Fig. 4 we plot the numerically integrated solutions for $n=1$ and $\gamma_{5}=0.1r_{h}M_{\rm pl}$ derived by using the boundary conditions (174)-(176) around $r=r_{h}$. The solutions in two asymptotic regimes smoothly join each other. Substituting Eqs. (174)-(176) into Eq. (169), it follows that $A_{1}$ approaches the constant $M_{\rm pl}\sqrt{2\mu/[(1+2n)(1-\mu)^{2}]}$ as $r\to r_{h}$. Indeed, the numerical simulation of Fig. 4 shows that $A_{1}$ starts to decrease from this finite value with increasing $r$ and it approaches another constant $P/\sqrt{2n+1}$ in the limit $r\to\infty$. From Eqs. (174)-(176) the effect of the coupling $\gamma_{5}$ on $f,h,A_{0}$ does not arise up to the order of $(\Delta x)^{3}$, but the difference between $f$ and $h$ still remains in the strong-gravity regime.

X Conclusions

In this paper, we have studied the static and spherically symmetric BH solutions in second-order generalized Proca theories with nonlinear derivative vector-field interactions. In Sec. II we derived the full background equations of motion for the action (1) and revisited the non-existence of hairy BH solutions for a massive Proca field given by the Lagrangian $G_{2}=m^{2}X$. More generally, we found that the bare $X$-dependent coupling $g_{2}(X)$ in $G_{2}$ without being multiplied by derivative terms like $F$ is generally the obstacle for the existence of hairy BH solutions. On the other hand, other derivative couplings like those appearing in the Lagrangians (3)-(6) can give rise to a variety of hairy BH solutions.

In Sec. III we reviewed the exact stealth BH solution which is known to exist for the specific coupling $X/4$ in $G_{4}$, and also obtained the extremal RN solution with the vanishing longitudinal mode (55) present for the model (54). Imposing the two conditions of the two identical metric components (33) and the constant norm of the vector field (34), we also constructed a family of exact BH solutions for other interactions $G_{3},G_{5}$ and intrinsic vector-mode couplings $G_{6},G_{2}=-2g_{4}(X)F,g_{5}$ in Sec. IV. The models allowing for their existence are given, respectively, by Eqs. (58), (67), (76), (84), and (88). The corresponding metrics are described by either RN, extremal RN, or Schwarzschild types.

In Sec. V we explored the existence of non-exact BH solutions for the power-law cubic interaction (89) with the non-vanishing longitudinal mode (90). Expanding $f,h,A_{0}$ around the horizon $r\simeq r_{h}$ as Eq. (91) for $n\geq 1$, it follows that the coupling $\beta_{3}$ arises as corrections to the RN solutions. On using Eqs. (93)-(95) as the boundary conditions around $r=r_{h}$ and numerically solving the equations of motion outside the horizon, we showed that the solutions smoothly connect to those at spatial infinity given by Eqs. (98)-(101). The constant $P$ appearing in Eq. (100) cannot be fixed by other two parameters $M$ and $\tilde{b}_{2}$, so it corresponds to a primary hair. As seen in Fig. 1, the difference between two metrics $f$ and $h$ is most significant in the regime of strong gravity ($r\lesssim 10r_{h}$).

In Sec. VI we showed that the power-law quartic coupling (106) with $n\geq 2$ gives rise to two branches characterized by (i) $A_{1}\neq 0$ with a primary Proca hair, and (ii) $A_{1}=0$ with a secondary Proca hair. In both cases, the solutions are regular throughout the horizon exterior with the difference between $f$ and $h$ induced by the coupling $\beta_{4}$, see Fig. 2 for the branch (i). For $n\geq 3$ there is the branch satisfying $A_{1}=\epsilon\sqrt{A_{0}^{2}/(fh)}$, but this merely corresponds to the RN solutions. For the power-law quintic coupling (134), it was shown in Sec. VII that the two asymptotic solutions around the horizon and at spatial infinity are discontinuous due to the divergence of $A_{1}$ at $h=1/(2n+1)$.

In Sec. VIII we studied BH solutions for the power-law sixth-order coupling (147) and found that the consistent branch for $n\geq 0$ corresponds to $A_{1}=0$. For $n=0$, there exists the $U(1)$ gauge symmetry, in which case the integration constant $P$ in Eq. (154) does not have physical meaning. For $n\geq 1$, the Proca hair $P$ is of the secondary type by reflecting the fact that the near-horizon expansions (158)-(160) contain only two parameters $\mu$ and $r_{h}$. In Fig. 3 the solutions in two asymptotic regimes $r\simeq r_{h}$ and $r\gg r_{h}$ smoothly join each other with the largest difference between $f$ and $h$ around the horizon. Compared to the cases of cubic and quartic couplings plotted in Figs. 1 and 2, $|f-h|$ decreases faster for increasing $r$.

In Sec. IX we also discussed the role of intrinsic vector modes given by the two power-law interactions (161). For the coupling $g_{4}(X)=\gamma_{4}(X/M_{\rm pl}^{2})^{n}$, the branch leading to the difference between the two metric components corresponds to $A_{1}=0$. The coupling $\gamma_{4}$ induces corrections to the RN solutions with the secondary Proca hair $P$ arising for $r\to\infty$. The coupling $g_{5}(X)=(\gamma_{5}/M_{\rm pl}^{2})(X/M_{\rm pl}^{2})^{n}$ gives rise to the non-vanishing $A_{1}$ branch (169). This is rather a specific case in which $A_{1}$ approaches finite constants for both the limits $r\to r_{h}$ and $r\to\infty$, with the secondary hair $P$ at spatial infinity. The numerical simulation of Fig. 4 for $n=1$ shows that, unlike the power-law quintic interaction (134), the solutions are regular throughout the horizon exterior.

In summary, for the power-law models with cubic and quartic couplings $G_{3}(X)$ and $G_{4}(X)$, we showed the existence of regular BH solutions with a primary hair related to the longitudinal vector propagation. The power-law couplings $G_{6}(X),g_{4}(X),g_{5}(X)$ associated with intrinsic vector modes generally give rise to regular BH solutions with a secondary hair. In both cases the deviation from GR is most significant in the strong-gravity regime, with the recovery of GR at spatial infinity. The deviation from GR can be potentially probed in future measurements of gravitational waves in the nonlinear regime of gravity.

There are several issues we did not address in this paper. While we focused on the static and spherically symmetric configurations as a first step, we can extend the analysis to asymptotically non-flat BH solutions and rotating BH solutions, see Refs. Minami ; Cisterna ; Babichev17 for such solutions in the specific models. In particular the existence of hairy Kerr BH solutions and stars with gravitational solitons was recently found for a complex Proca field Herdeiro2 ; Brito , so it is of interest to study what happens in the presence of vector-field derivative couplings. In addition, the stability analysis of BH solutions against odd- and even-parity perturbations along the line of Refs. Suyama may constrain the strength of derivative couplings studied in this paper. It will also be interesting to investigate solutions of NSs and other (exotic) compact objects in generalized Proca theories, where the deviations from GR may be more evident than the BH case, see Ref. Chagoya2 for NS solutions in the model with nonminimal coupling to the Einstein tensor. These interesting issues will be left for future works.

Appendix: Coefficients in the gravitational

equations

In Eqs. (14)-(16) the coefficients $c_{1,2,\cdots,19}$ are given by

[TABLE]

Acknowledgements

LH thanks financial support from Dr. Max Rössler, the Walter Haefner Foundation and the ETH Zurich Foundation. RK is supported by the Grant-in-Aid for Young Scientists B of the JSPS No. 17K14297. ST is supported by the Grant-in-Aid for Scientific Research Fund of the JSPS No. 16K05359 and MEXT KAKENHI Grant-in-Aid for Scientific Research on Innovative Areas “Cosmic Acceleration” (No. 15H05890). MM was supported by FCT-Portugal through Grant No. SFRH/BPD/88299/2012. We also wish to acknowledge the hospitality of the Yukawa Institute for Theoretical Physics where part of this work was conducted during the workshop on gravity and cosmology for young researchers (YITP-X-16-10).