{'fraction_non_alphanumeric': 0.03198003198003198, 'fraction_numerical': 0.011349011349011348, 'mean_word_length': 4.407423028258119, 'pattern_counts': {'":': 0, '<': 3, '<?xml version=': 0, '>': 1, 'https://': 0, 'lorem ipsum': 0, 'www.': 0, 'xml': 0}, 'pii_count': 0, 'substrings_counts': 0, 'word_list_counts': {'cursed_substrings.json': 4, 'profanity_word_list.json': 0, 'sexual_word_list.json': 0, 'zh_pornsignals.json': 0}}
Self-interacting scalar field distributions around Schwarzschild black holes

Long-lived configurations of massive scalar fields around black holes may form if the coupling between the mass of the scalar field and the mass of the black hole is very small. In this work we analyze the effect of self-interaction in the distribution of the long-lived cloud surrounding a static black hole. We consider both attractive and repulsive self-interactions. By solving numerically the Klein Gordon equation on a fixed background in the frequency domain, we find that the spatial distribution of quasi stationary states may be larger as compared to the non interacting case. We performed a time evolution to determine the effect of the self-interaction on the life time of the configurations our findings indicate that the contribution of the self-interaction is subdominant.


I. INTRODUCTION
Scalar fields are promising candidates to explain the nature of Dark Matter. In a Scalar Field Dark Matter model it is assumed that a classical field makes up the main fraction of the dark matter of the Universe. The description relies on the fact that a coherently, weakly interacting, oscillating light scalar field formed a Bose Einstein Condensate (BEC) in a non-relativistic and low-momentum state [1][2][3][4][5][6][7][8]. The BEC has a phase space density that enables it to describe the density profile of the galactic dark matter halo for a convenient approximations. Amongst the many proposed candidates, Axions, particles introduced by Peccei and Quinn [9], are scalar fields with a non-zero vacuum expectation value and keep the CP invariance of the strong interactions in the Lagrangian involving all Yukawa couplings. The effects of such scalar field and other axion-like particles (with a lower mass) in astrophysical and cosmological scenarios have been investigated extensively; see, e.g., Refs. [10][11][12][13][14][15][16][17][18][19].
Since the detection of gravitational waves by the LIGO-Virgo-Kagra collaboration [20][21][22][23] and the recent developments in electromagnetic observations, in particular the image of the central black hole and its shadow of M87 and Sagittarius A* captured by the Event Horizon Telescope [24][25][26], there is a renewed interest in the physical processes that may occur in the vicinity of compact objects.
Massive bosonic fields may form quasi-bound states around a black holes [27][28][29][30][31][32][33]. In a Schwarzschild background, all of these states are unstable to decay, leaking part of the field towards the black hole [34,35]. These scalar configurations are characterized by instability timescales which are much longer than the timescale set by the mass of the central black hole. Because the low rate of decay, the scalar field configurations may remain surrounding a black hole for large time scales depending on the values of the parameters involved [36]. Therefore, scalar fields around black holes represent a very convenient set up to model the super massive black hole surrounded by a dark matter halo in galaxies. Quasi stationary solutions to the Klein-Gordon equation on a Schwarzschild background in that context, were described in detail by Barranco et al. in Refs. [36][37][38]. The results obtained in those references, indicate that it is possible that scalar field halos may last for cosmological time scales around super massive black holes. These long lived configurations composed by a complex, massive and non-self-interacting scalar field, were found to be characterized essentially by two parameters, namely the integer , associated with the angular distribution of the field and the dimensionless quantity formed by the mass coupling between the mass of the black hole M , and the mass of the scalar particle m, GM m/ c, which can also be interpreted as one half the ratio between the black hole horizon radius r h = 2GM/c 2 and the characteristic wavelength of the scalar field λ φ = c/m.
For a scalar particle, the simplest non gravitational interaction is a quartic self-interaction. This generalization can be achieved by expanding a potential about a symmetric minimum, and realizing that the quartic term is the most important interaction term for small amplitudes. From the point of view of a field theory the quartic self-interaction is the largest value in the exponent that allows a renormalizable theory. In an astrophysical scenario, Colpi et al. [39] show that self-interactions in boson stars can produce significant phenomenological changes. In particular, they show that the upper limit on the mass of a boson stars increases notably compared to the noninteracting case. Such findings show that boson stars can have masses even larger than a solar mass. Furthermore, in some regimes self-interacting boson stars are compact enough to be considered black hole mimickers [40][41][42][43][44][45]. From a cosmological perspective, the role of self-interaction is of great importance in dark matter models, particularly within complex scalar field models as in Refs. [46,47] where it has been shown that with a repulsive self-interaction, the scalar field goes through a radiation-like stage in the Friedman evolution which in turn increases the effective number of relativistic degrees of freedom that matches the estimates of big bang nucleosynthesis (see also [48] for further discussion). Quasi stationary states of attracting self-interacting scalar fields where found in the context of gravitational collapse of unstable boson stars in Ref. [49]. The evolution of self-interacting boson stars in spherical symmetry was performed by solving the Einstein Klein Gordon system numerically. When describing an unstable configuration, boson stars collapsed forming a black hole surrounded by a remnant of scalar field leaving long lasting states around the newly formed black hole.
In this work, we solve the Klein-Gordon equation describing a self-interacting scalar field in the background of a Schwarzschild black hole. We assume the field oscillates coherently and analyse the system in the frequency domain to find resonant states. We analyze in detail the consequences of a quartic self-interaction in the distribution of the field and present a thorough analysis of the quasi stationary solutions. We observe that the self-coupling have important consequences in the phenomenology associated with scalar field distributions. We further analyse the behaviour in time of the scalar distribution.
The paper is organized as follows: In Sec.II we present the formulation, the background spacetime and method of construction of solutions. In Sec. III we present some examples of solutions of the Klein Gordon equation and describe the properties of resonant modes. In Sec. IV we discuss the time development of resonant modes. Finally, we summarize our results and present concluding remarks in Sec.V. Throughout the paper we use geometric units such that c = 1 = G.

II. SET UP
In this work, we will focus on the regime where self-interaction effects on the field become significant before gravitational backreaction and therefore, we shall restrict to a fixed black hole spacetime. We consider a complex massive scalar field Φ with mass m described by the action where µ = m/ and for convenience we consider λ > 0. The attractive or repulsive nature of the self-interaction is set by η = 1 or η = −1 respectively. The stress-energy tensor of the scalar field reads The conservation of the stress energy tensor ∇ ν T µν = 0 provides the corresponding equation of motion for the field, the Klein-Gordon equation: and its complex conjugate. We consider the background space-time as the non-rotating Schwarzschild black hole in Boyer-Lindquist coordinates (t, r, θ, ϕ) with element of line given by with N (r) = 1 − 2M/r, and dΩ 2 = dθ 2 + sin 2 θdϕ 2 . We assume that the scalar field is spherically symmetric and thus can be written as where ω is a real number. The resulting Klein-Gordon equation takes the form where, We will use the fact that the static spacetime at hand has a timelike Killing vector to define a diagnostic quantity with respect to the stress energy of the system. The energy density of the scalar field measured by a static observer, ρ = −T 0 0 , provides a useful measure of the distribution of the field in the spacetime, explicitly is given as This quantity will be used to characterize the scalar field distribution in the following sections. Another quantity of interest is the total mass-energy storaged by the scalar field obtained by the integration of the energy density This quantity gives an estimation of the amount of scalar field in the spacetime that can be compared with the mass of the black hole.

A. Scaling properties
In order to explore the space of solutions of the Klein Gordon equation (6), we define a new function v(r) as such that the function U(r) becomes The function v encodes the self-interaction parameter and allow us to explore the solution space of Eq. (6) providing an infinite set of solutions for each value of λ. Additionally, in order to better resolve the behaviour of the field in the region close to the horizon, we use the tortoise coordinate, r * = r + 2M ln(r/2M − 1), to rewrite Eq. (6) as where r is given implicitly in the definition of r * . One can thus eliminate λ from the numerical task because always enters as a scale. The equation (12), posses a further rescaling property, which is specified in Table I. This rescaling allows us to use dimensionless quantities such as r/M , M µ, M ω, etc. to specify any quantity in terms of the mass of the black hole. The dependence on the parameters in the Klein Gordon equation is in M µ, and the value of η (that is +1 or −1).   where the amplitude A and the phase δ are real numbers. Notice however, that the function (13) together with the ansatz (5), contains no physical solutions since it represents both out-coming and in-going modes. The physical situation we shall consider, a field distribution around a black hole, excludes outgoing modes. In the following section we will address this issue in more detail. Let us consider values of ω such that ω < µ. Then, in the limit r * → ∞ Eq. (12) admits exponentially decaying or growing solutions of the form v(r * ) ≈ C exp(−kr * ) + D exp(+kr * ) with k = µ 2 − ω 2 . Since we are interested in describing localized solutions, consistent with the asymptotically flat metric, we will focus on the solutions with exponential decay and set D = 0.
The previous descriptions concerns only the asymptotic behaviour of the solutions of Eq. (12). With only this information, the possible values that ω can take are infinite and the spectra is continuous. In the next section we shall show that there is a set of discrete values of ω that allows solutions in which the in-going modes have a larger amplitude than the out-going modes. These solutions are thus closer to the physical description of having an event horizon as the left boundary.
In order to find the solutions of Eq. (12) we integrate it numerically as follows. We fix η, µ and choose some value ω. Then we set a value of v at the outermost point of the numerical grid r * max and integrate the equation inwards up to the left boundary r * min . Fig. 1 shows the plot of v for M µ = 0.14 and M ω = 0.1370 obtained by this procedure. There is a region in which the behaviour (13) is reached and the function v is completely characterized by its maximum amplitude in that region A, and the asymptotic behaviour. Furthermore, we have found that for each solution v(r), the quantity v max defined as the local maximum attained by the function v in the far horizon region, can be used to characterize that solution. In the next section we will show that comparing the amplitude in the far region, v max , with the amplitude near the horizon, A, it is possible to find a discrete set of values of ω that allows quasi bound or resonant configurations.
Decaying solution are characterized by the amplitudes A and v max for a given value of ω. Figs. 2 display sets of solutions with the desired asymptotic behaviour for some values of ω. We will refer solutions with v max → 0 as solutions in the weak self-interacting regime. As the value of v max increases, the term v 2 in Eq. (12) becomes important and defines solutions in the strong self-interacting regime. These functions however, are not physical in the sense that we have not yet imposed the boundary condition at the left. That is, these functions contain outgoing modes coming from the left which is not physically possible.

III. RESONANT STATES
It is possible to quantify the amount of the scalar field that accumulates in the outer region with respect to the amplitude near the horizon by comparing the values of v max and A. In Fig. 3 we show a projection map of the ratio A/v max using several values of ωM for µM = 0.14. The left panel corresponds to η = −1 and the right panel corresponds to η = 1. In order to better visualize the frequencies of the resonant states and get some understanding on their change in v max , Fig, 4 shows constant v max cuts of the projected surface Fig. 3. In these plots, resonant states  Other results related to the scalar field distribution as well as the issue about the lifetime of these states will be addressed below.

A. Properties of resonant states
As stated before, the effect of the magnitude of λ on the radial profile of the scalar field can be described in terms of v max . The first scenario we examine is the one where η = −1 and ω is fixed. As a particular example, by taking the frequency M ω = 0.13837 (that corresponds to the first resonant state with λ = 0), configurations in the weak self-interacting regime are almost indistinguishable from the case with λ = 0. In the strong regime the profiles differ slightly because these configurations are no longer resonant states. Fig. 6 displays the profiles v(r * ) and the radial energy density ρ e = 4πr 2 ρ, with η = −1 and M µ = 0.14 for some representative values of v max in both weak and strong self-interacting regimes, the case with λ = 0 is also shown for comparison. The configurations in Fig. 6 correspond to the asterisk symbols in Fig. 3. For large values of v max , the ratio A/v max becomes larger and it may happen that for large enough values of v max along the fixed frequency, a profile with one or more nodes is found. Like in the λ = 0 case, overtones correspond to solutions with increasing number of nodes as ω approaches µ. In the negative self-interaction case, η = −1, the effect of increasing λ in the mass density is that the scalar field distribution spreads over a larger region and the maximum moves away from the horizon. The normalized mass density is shown in the right panel of Fig. 6. Given the behaviour (13), solutions of Eq. (12) produce an infinite energy density at the horizon. A positive self-interaction η = +1 produces analogous changes in the radial profile of the scalar field, in this case however, the field concentrates closer the horizon and the changes are smaller in magnitude. Another important difference is that for large values of v max resonant solutions cease to exist. Fig. 7 displays the radial scalar field profile and the radial energy density for η = +1 In the following we shall focus on the fundamental resonant states as v max increases. Fig. 8 shows the scalar field profile of the fundamental resonant state for both η = −1 and η = +1. Each resonant state is characterized by its frequency and v max . As v max increases the amplitude of the radial profile v(r) increases for both η = −1 and η = +1. The radial density presents the same behaviuor as is displayed in Fig. 9. For η = +1 the effective size of the configuration increases in the strong interacting regime. Another important effect of the self-interaction on the distribution of the scalar field is in the mass of the configuration. Self-interacting configurations are heavier than their non interacting partners. In the weak self-interacting regime, the mass increases almost linearly with v max for both η = −1 and η = +1. In the strong regime however, the growth slows down for η = −1. Fig. 10 shows the mass M Φ as a function of v max for some values of the gravitational coupling M µ considering η − 1 and η = +1.
In the non self-interacting case, an effective potential of the time-independent Schrödinger-like equation (12) (with λ = 0) can be defined. In the self-interacting case, such effective potential has the form: The effective potential interpretation may help in the characterization of the solutions (see [36,37] and references therein) however, for λ = 0, V eff can only be obtained a posteriori since the knowledge of v(r) (and consequently u(r)) is required. In any case, the effective potential may be useful to determine the resonance band formed by the  parameters M µ and M ω for which the solutions have the possibility of being concentrated. The existence of the potential well in the λ = 0 case is granted by the condition M µ < 1/4 [36]. A similar bound however, can not be found for λ = 0.

IV. TIME DOMAIN DESCRIPTION
In this section we are interested in describing the behaviour in time of the radial profiles of resonant states once the harmonic time dependence (5) is relaxed. In particular we are interested in the life-halftime of the configuration. To proceed, we solve the Klein-Gordon equation in the time domain using the in-going Kerr-Schild coordinate system. These coordinates are more convenient for numerical analysis because a constant time hypersurface is nonsingular and horizon penetrating. This characteristic is important in order to impose convenient boundary conditions.
The relation between the Kerr-Schild and Boyer Lindquist coordinates is given through the time transformatioñ where r * is the tortoise coordinate defined above. In these coordinates, the Schwarzschild metric takes the form In terms of Arnowit Desser Misner (ADM) variables, see [50] for an introduction, the lapse function α, the shift vector β r , and the induced 3-metric can be read from Eq. (16) as The complex scalar field can be split into two real scalar fields according to Φ = φ + iφ c . However, due to the spherical symmetry of the problem at hand, both fields evolve following the same dynamical equation. Consequently, in the following we shall refer to the dynamical properties of φ. In order to solve numerically the Klein Gordon Equation, in the background metric (16), we introduce the auxiliary first order functions where we have drop the tilde in the time coordinate. Equation (18) is thus equivalent to the following system of evolution equations for φ, π and ψ ∂ t φ = r r + 2M π + 2M r ψ , As a diagnostic quantity we use the total energy of the scalar field that in Kerr Schild coordinates is written as where ρ(r) = 1 2 α 2 π 2 + 2πψβ 2 + γ rr ψ 2 + µ 2 φ 2 + 1 2 λφ 4 .
The evolution equations for the radial components were solved with the 1+1 dimensional PDE solver described in [51] and used in several scenarios [52,53]. The time evolution is done via the method of lines with a third order TVD Runge-Kutta algorithm. The spatial derivatives are approximated with a second order symmetric finite difference stencil. A standard fourth order dissipation term was added in order to guarantee the stability of the scheme. The evolution scheme is complemented imposing an outgoing wave boundary condition.  Table II.

A. Initial data: Resonant modes
One may want to evolve directly a resonant state constructed in the frequency domain as described above in order to determine the time rate of decay, however, those states diverge at the horizon. In order to use the radial function of resonant states as initial data we employ a regularization procedure as described in [37] to mimic a quasi-resonant state and smooth its behaviour in the region close the horizon. Such procedure allows us to set an initial data close enough to the resonant states found in previous sections. We have used as initial data at t = 0 the radial profile of the field φ, obtained in III with M µ = 0.14. For our analysis we consider states in both weak and strong self-interacting regimes. Fig. 12 shows the initial radial distribution φ(r) with η = +1 and η = −1.
As a result of the evolution, the scalar field present an oscillating behaviour and a slow rate of decay. By measuring the field amplitude at a fixed point r = r 1 one obtains a time series for the amplitude at that point. We may thus perform a Fast Fourier Transform to obtain the frequency of oscillation. The power spectrum obtained from the time evolution, shows that the field oscillates with a frequency that corresponds to the frequency obtained in III. In Fig. 13 we also plot the total energy of the scalar configuration as a function of time E(t), defined in Eq. (21). Fig. 13 shows that the decay rate of energy is exponential after an initial transient state. We assume that the energy has the form the form E(t) ∼ exp(−st) and perform a linear fit of ln(E(t)/E 0 ) with t in order to find the value of the rate of s. The half-lifetime, of the configuration is thus given by t 1/2 /M = ln(2)/s. We have further considered resonant modes with other values of the coupling M µ < 1/4 and compute the rate of decay of the energy. The results are summarized in Table II. Despite the fact the initial morphology of the scalar field state is different for configurations with η = +1 and η = −1 the rate of decay is of the same order of magnitude for a given value of M µ. We conclude that the effect of λ on the rate of decay is subdominant as compared with the effect of M µ. Notice however, that the effective size of the modes with η = +1 and large enough values of λ is slightly larger than the non interacting case,   thus if the leaking rate is of the same order of magnitude in both cases, self-interacting clouds will last longer since the initial distribution in the interacting case is larger. Values of the boson mass µ motivated by dark matter scalar field models correspond to µ ∼ 10 −24 eV, and for a black hole with mass M ∼ 10 8 M , the gravitational coupling is of the order of M µ ∼ 10 −6 . This value is slightly smaller than the ones considered in this work, such small values cannot be reached due to the limitations of the numerical integration. However if we extrapolate our findings as done in Ref. [38] resonant modes may last about 10 9 years.

V. CONCLUDING REMARKS
In this paper we studied the role played by the self-interaction in quasi resonant states of test scalar fields around a Schwarzschild black hole. We assumed that self-interaction is mediated by a term ∼ λ Φ 4 and investigated the properties of the field distribution that forms in the vicinity of black holes in the case of both attractive and repulsive self-interactions. We characterized these cases by means of a parameter η, leaving the parameter λ as a non negative quantity. We rewrite the Klein Gordon equation in terms of a function v(r), allowing us to eliminate the explicit dependence in λ and solve it numerically. With this function one can characterize each solution with its amplitude in the near and far horizon regions. We focus our analysis on the resonant modes, which are exponentially decaying solutions at infinity concentrated well outside the event horizon. Furthermore, like in the non interacting case, resonant modes posses a definite frequency of oscillation.
The first conclusion that can be drawn regarding the role of self-interaction in resonant states is the fact that the values of the discrete frequencies change with respect to the non-interacting case. For η = −1 the frequency decreases and for η = +1 increases compared with the corresponding values of the frequencies without self-interaction.
We have found that in the presence of interactions the size and distribution of the scalar field changes depending on the value of λ. In the regime of strong self-interaction with η = +1 the size of the resonant scalar cloud distribution increases considerably and the scalar field tend to concentrate more in a region far from the horizon as compared to the case with no self-interaction. Even more, we have found that for large enough values of λ, resonant solutions do not exist. Regarding the self-interaction with η = −1, the field concentrates closer the horizon with almost no change in size. In this way, we can conclude that, for η = +1 and large self-interaction the size of the distribution is significantly larger than the η = −1 and the non-interacting cases.
From a classical perspective this happens because the particles of the configuration tend to concentrate outside the horizon while gravity and attractive self-interactions tend to shrink the configuration towards the black hole. We further investigate the life time of resonant states by means of a numerical time evolution. We found that despite the important role played by the self-interaction in the spatial distribution of the scalar field around the black hole, its role in the life time is negligible as compared to the effect of the mass term. In the scenario described in this work we focus on the regime where self-interaction dominates over self gravitation effects and the parameter λ entered in the equations as a global scale. However, when the spacetime reacts to the presence of the scalar field, the term with λ enters explicitly on the energy density and thus a stronger effect on the properties of the whole configuration is found. A further study, considering the back reaction of the spacetime is under way.