Quasibound states of scalar fields in the consistent 4D Einstein-Gauss-Bonnet-(Anti-)de Sitter gravity

We examine the interaction between massless scalar fields and the gravitational field generated by a black hole solution that was recently obtained in the consistent well-defined 4-dimensional Einstein-Gauss-Bonnet gravity with a cosmological constant. In order to do this, we calculate quasibound state frequencies of scalar fields for the spherically symmetric black hole in the consistent 4-dimensional Einstein-Gauss-Bonnet-de Sitter and Anti-de Sitter theories. The expression for the quasibound states is obtained by using the polynomial condition associated to the Heun functions, and their values are overdamped. We also demonstrate the stability of the systems.


I. INTRODUCTION
In the present paper, we extend these results by studying the quasibound states for the conformally coupled massless scalar field in both 4D Einstein-Gauss-Bonnet-de Sitter (4DEGBdS) and 4D Einstein-Gauss-Bonnet-Anti-de Sitter (4DEGBAdS) backgrounds. These calculations are carried out in the stability region of the GB coupling constant [15]. It is worth emphasizing that the cosmological constant should be, in principle, a negligible factor in the studies of black holes, since its current value is quite small and hence it does not play any role in the black hole dynamics. However, an effective cosmological constant can appear when regions with high curvature in alternative/modified theories of gravity are considered. In this case, such term could, in principle, describe the interaction of the black hole with a dark energy component of the universe. Therefore, it is interesting to investigate the role played by the cosmological constant within the 4DEGB theory.
In this work, we will calculate the spectrum of quasibound states, as well as the corresponding angular and radial wave eigenfunctions, for massless scalar particles in the 4DEGB black hole solutions with a cosmological constant by using the polynomial condition of the Heun functions, the so-called Vieira-Bezerra-Kokkotas (VBK) approach (for details, see Refs. [16,17]). We show that the quasibound states depend on the GB coupling constant, a, as well as on the cosmological constant, Λ. We also investigate the stability of the systems.
The quasibound states [18,19] are a kind of wave phenomena occurring near the black hole exterior event horizon. They are localized in the black hole potential well, which means that there exist a flux of particles crossing into the black hole. Thus, the spectrum of quasibound states is constituted by complex frequencies, which can be expressed as ω = ω R + iω I , where ω R and ω I are the real and imaginary parts, respectively. The real part of the resonant frequency is the oscillation frequency, while the imaginary part determines the stability of the system. Thus, the wave solution is said to be stable when the imaginary part of the resonant frequency is negative (ω I < 0), which signalizes the existence of a decay rate with the time. Otherwise, the wave solution is unstable when the imaginary part of the resonant frequency is positive (ω I > 0)and therefore, contrarily to the previous case, there is a growth rate with the time.
The paper is organized as follows. In Sec. II, we introduce the general metric corresponding to the 4DEGB black hole spacetimes with a cosmological constant. In Sec. III, we separate the conformally coupled massless Klein-Gordon equation in the background under consideration. In Sec. IV, we discuss the aforementioned metric for the case of a positive cosmological constant (4DEGBdS black hole) and find an exact solution for the radial equation in this background. Then, we obtain an expression for the quasibound states. In Sec. V, we repeat this analysis for the case of a negative cosmological constant (4DEGBAdS black hole). Finally, in Sec. VI, we summarize the obtained results. Here we adopt the natural units where G ≡ c ≡ ≡ 1.

II. THE CONSISTENT WELL-DEFINED 4-DIMENSIONAL EINSTEIN-GAUSS-BONNET GRAVITY WITH A COSMOLOGICAL CONSTANT
In the studies of black hole radiation, which means the investigation of emission/transmission/reflection of any kind of quasispectrum, it is crucial that the black hole solutions should be solutions of both truly 4-dimensional AGM (4DAGM) [7][8][9] and extra scalar degrees of freedom [20][21][22][23][24][25][26][27] theories. Thus, we can call it as a consistent well-defined 4-dimensional Einstein-Gauss-Bonnet theory of gravity, which includes the black hole solutions with a cosmological constant.
Then, by taking into account all of the aforementioned approaches, we will briefly review the basic ideas behind it in order to present the 4DEGB black hole solutions with a cosmological constant.
Let us start by stating that the consistent well-defined 4DEGB gravity is based on the ADM formalism, where the 4D metric is defined by with N being the so-called lapse function, N i is a shift vector, and γ ij is the spatial metric. In a 4D spacetime, the following gauge condition is valid: where D k is the covariant derivative and π ij is the canonical momentum conjugate to γ ij . Then, the gravitational action, S, is defined as with where κ is the gravitational coupling (strength) constant, a is the (rescaled) GB coupling constant, R ij is the Ricci tensor, R is the Ricci scalar, λ GF is a gauge-fixing (GF) parameter, and M = M i i . Here, the dot denotes the derivative with respect to the time t. The final step is to choose the GF parameter such that λ GF = 0 and appropriately rescale both the gravitational and GB coupling constants, as well as to take M as the total mass centered at the origin of the system of coordinates. Therefore, an exact solution describing the 4DEGB black hole spacetime with a cosmological constant has the following form where the metric function, f ± (r), is given by The metric function f + (r) corresponds to an asymptotically de Sitter spacetime in the limit when Λ → 0. On the other hand, the metric function f − (r) corresponds to an asymptotically flat spacetime when Λ → 0. Here, we focus on the "minus" case and hence we set f (r) ≡ f − (r), such that the surface equation can be written as In general, Eq. (10) has four solutions, which can be complex numbers, as well as positive and negative real numbers.
In the present work, we will consider only the real solutions, including the negative ones, which may describe the (Cauchy) "interior" cosmological horizon. It is worth emphasizing that we will obtain these solutions directly from Eq. (9) by using the built-in symbol Solve on Wolfram Mathematica 12.3. Next, we will consider massless scalar fields in the 4DEGB black hole spacetime with a cosmological constant given by Eq. (8).

III. THE KLEIN-GORDON EQUATION
We are interested in some basic characteristics of these 4DEGB black hole spacetime with a cosmological constant, in particular the ones related to their interaction with quantum fields, including the Hawking radiation, and classical scalar wave scattering such as the quasibound states (QBSs). To this end, we want to consider the conformally coupled massless scalar field as a probe, whose covariant equation of motion is given by In order to obtain solutions of the conformally coupled massless Klein-Gordon equation (11), and due to stationarity and axisymmetry, we use the following separation ansatz where U (r) is the radial function, P (θ) is the polar angle function, m (∈ Z) is the magnetic or azimuthal quantum number, and ω is the frequency (energy, in the natural units). By substituting Eq. (12) into Eq. (11), we obtain two ordinary differential equations, namely, where λ is the separation constant and the prime denotes differentiation of the polar and radial functions with respect to ϑ and r, respectively. The Ricci scalar R can be written in terms of the metric function f (r) as The final (parametrized) expression for the Ricci scalar R will depend on the parametrization adopted for the metric function f (r). The general solution of the polar equation (13) is given in terms of the associated Legendre functions P (θ) = P m ν (cos θ) with general degree ν (which can be a complex number C) and m ≥ 0 ∈ Z, such that λ = ν(ν + 1). Now, let us turn our attention to the solution of the radial equation. Firstly, we will analyze the solutions of the surface equation f (r) = 0, that is, the number of event horizons, and then write the metric function f (r) in a more convenient form to write Eq. (14) as a Heun-type equation that is useful to obtain exact solutions, which is useful to study the Hawking radiation and QBS spectra, which will be presented in the next sections.

IV. THE 4DEGB-DE SITTER BLACK HOLE SPACETIME
For a 4DEGBdS black hole spacetime, we fix M = 1/2 and then calculate the black hole event horizons according to the values of the cosmological constant (Λ > 0), as well as for the GB coupling constant in the stability region −0.20 ≤ a ≤ 0.20 [15].
The behavior of the metric function f (r) is shown in Fig. 1. In addition, the horizons are shown in the inlay plots as functions of the cosmological constant Λ. We can see that the surface equation f (r) = 0 has three real solutions when a < 0, otherwise, there are four real solutions when a > 0.
A. (dS) Case 1: a < 0 In this case, we adopt the following parametrization for the metric function, where r c is the (positive) "exterior" de Sitter cosmological horizon, r + is the (positive) "exterior" black hole event horizon, and r −− is the (negative) "interior" de Sitter cosmological horizon. Then, with the metric function given by Eq. (16), we will show that Eq. (14) is totally appropriate to study QBSs with purely ingoing boundary conditions at the exterior black hole event horizon and vanishing boundary conditions at infinity, since it is a (general) Heun-type equation [28] with three finite regular singularities (associated to the three event horizons) and one regular singularity at (spatial) infinity. Now, we follow the steps described in the VBK approach [16,17] to obtain the analytical solution of Eq. (14) in the 4DEGBdS black hole spacetime with a < 0 (without the assumption of specific boundary conditions). First of all, we need to define a new radial coordinate, x, as These definitions move the three singularities (r −− , r + , r c ) to the points (0, 1, ∞); the exterior de Sitter cosmological horizon is now a regular singularity at (spatial) infinity. The next step is to perform an F-homotopic transformation U (x) → y(x) given by where Thus, by substituting Eqs. (16)- (22) into Eq. (14), we get where Equation (23) has the (canonical) form of a general Heun equation, where y(x) ≡ HeunG(η, q; α, β, γ, δ; x) denotes the general Heun function, which is the solution corresponding to the exponent 0 at x = 0 and assumes the value 1 there. If γ is not a negative integer, then from the Fuchs-Frobenius theory it follows that the HeunG(η, q; α, β, γ, δ; x) exists, is analytic in the disk |x| < 1, and has a Maclaurin expansion given by with −qc 0 + ηγc 1 = 0, P n c n−1 − (Q n + q)c n + X n c n+1 = 0 (n ≥ 1), (28) where P n = (n − 1 + α)(n − 1 + β), Q n = n[(n − 1 + γ)(1 + η) + ηδ + ǫ], X n = (n + 1)(n + γ)η.
Here, the normalization c 0 = 1 was adopted by Karl Heun [29]. On the other hand, if γ is not a positive integer, the solution of Eq. (23) corresponding to the exponent 1 − γ at x = 0 is x 1−γ HeunG(η, (ηδ + ǫ)(1 − γ) + q; α + 1 − γ, β + 1 − γ, 2 − γ, δ; x). Therefore, the analytical solution for the radial part of the conformally coupled massless Klein-Gordon equation, in the 4DEGBdS black hole spacetime with a < 0, can be written as where C 1,j and C 2,j are constants to be determined, and j = {η, 0, 1, ∞} labels the solution at each singular point. Thus, the pair of linearly independent solutions at x = 0 (r = r −− ) is given by Similarly, the pair of linearly independent solutions corresponding to the exponents 0 and 1 − δ at x = 1 (r = r + ) is given by The pair of linearly independent solutions corresponding to the exponents 0 and 1 − ǫ at x = η (r = ∞) is given by Finally, the pair of linearly independent solutions corresponding to the exponents α and β at x = ∞ (r = r c ) is given by The assumption of a specific asymptotic behavior on the aforementioned analytical solutions near the exterior black hole event horizon (r → r + ) and spatial infinity (r → r c ) can lead to various physical solutions. Next, we will use the VBK approach [16,17] to derive the characteristic resonance equation and then find the spectrum of resonant frequencies related to quasibound states. In order to compute this quasispectrum, we need to impose two boundary conditions on the radial solution: it should describe an ingoing wave at the exterior black hole event horizon and tend to zero far from the black hole at spatial infinity.

Hawking radiation
In the limit when r → r + , which implies that x → 1, the radial solution given by Eq. (30) behaves as where C 1,1 and C 2,1 are constants to be determined, in which all remaining constants were included. Then, we can include the time dependence and hence this solution is written as Here, the ingoing, Ψ in,1 , and outgoing, Ψ out,1 , scalar wave solutions at the exterior black hole event horizon are given by In addition, the gravitational acceleration, κ + , on the exterior black hole event horizon is Thus, by using the analytic continuation described in the VBK approach [16,17], we can obtain the relative scattering probability, Γ + , and the exact spectrum of Hawking-Unruh radiation,N ω , given by where k B is the Boltzmann constant, and T + (= |κ + |/2πk B ) is the Hawking temperature at the exterior black hole event horizon. This spectrum has a thermal character and hence it is analogous to the spectrum of black body radiation.

Quasibound states
In order to satisfy purely ingoing boundary conditions, from the asymptotic behavior of the radial solution at the exterior black hole event horizon described by Eq. (40), we must impose that C 2,1 = 0. Now, let us analyze under which circumstances the radial solution vanishes at spatial infinity. To do this, we have to use the two linearly independent solutions of the general Heun equation at spatial infinity, given by Eqs. (37) and (38), and then obtains the following asymptotic behavior for the radial solution where with Thus, the radial solution fully satisfies the second boundary condition for QBSs if the sign of the real part of σ is such that Re[σ] > 0; whereas the radial solution diverges at spatial infinity if Re[σ] < 0. Then, the final asymptotic behavior of the radial solution at spatial infinity will be determined when we know the values of the coefficient σ, which depends on the frequencies ω, and on the parameter α. In order to determine this, we use a polynomial condition for the general Heun functions to match the two asymptotic solutions of the scalar radial equation in their common overlap region. It is known that the general Heun function becomes a (class I) polynomial of degree n if and only if the following two conditions are satisfied [28]: where n = 0, 1, 2, . . .. Equation (50) is the so-called α condition, which provides the expression for the frequencies ω n . On the other hand, the accessory parameter q must be appropriately chosen so that Eqs. (28) and (29) are consistent, which means that it is also necessary for the accessory parameter q to be an eigenvalue of the general Heun equation, calculated via Eq. (51). Furthermore, the accessory parameter q contains the separation constant λ, which indicates that we could obtain the eigenvalues of the separation constant λ n (ω n ), corresponding to the appropriate eigenvalues of the frequencies ω n , from the polynomial solution for the radial equation and then use it to show the (regular) angular behavior of massless scalar QBSs in the background under consideration (for details, see Ref. [30]).
In the present case, by imposing the condition given by Eq. (50), with the parameter α given by Eq. (24), we obtain a first order equation for ω, whose solution is the exact spectrum of QBSs: where n can be called the principal quantum number. As a result, the quasispectrum of frequency eigenvalues of the conformally coupled massless Klein-Gordon equation (11) is discrete, that is, the frequency quasispectrum consists of an infinite number of discrete frequency levels corresponding to quasibound states. Furthermore, note that it is a function of the horizons r j , and therefore the QBSs are also function of continuous variables, such as the cosmological constant, GB coupling constant, and cosmological radius.
In Table I we present the QBSs ω n as functions of the GB coupling constant a, for some values of the cosmological constant Λ. In addition, we also present the behavior of the QBSs ω n in Fig. 2. We can conclude that the spectrum given by Eq. (52) is physically admissible, since the real part of σ is positive, and therefore it represents QBS frequencies for conformally coupled massless scalars in the 4DEGBdS black hole spacetime with a < 0. Furthermore, we can see that the decay is overdamped, since this quasispectrum is purely imaginary.

Radial wave eigenfunctions
The radial wave eigenfunctions, which are related to QBSs of conformally coupled massless scalars propagating in the 4DEGBdS black hole spacetime with a < 0, can be obtained by using the condition given by Eq. (51), and imposing the condition given by Eq. (50), as it is described in the VBK approach [16,17].
As we explained before, these polynomial radial eigenfunctions are related to the appropriate determination of the eigenvalue q. Since it is calculated via Eq. (51), we index its solutions by a parameter s(= 1, 2, 3, . . .), which can be conveniently denoted by q n;s . Thus, the corresponding general Heun polynomials are now denoted as HeunGp n;s (x). Therefore, the QBS radial wave eigenfunctions for conformally coupled massless scalars propagating in a 4DEGBdS black hole spacetime with a < 0 are given by HeunGp n;s (η, q n,s ; −n, β, γ, δ; x), where C n;s is a constant to be determined. Next, we calculate the general Heun polynomials related to the fundamental and first excited modes. The general Heun polynomial for the fundamental mode n = 0 is given by where the eigenvalue q 0;1 must obey whose unique solution (s = 1) is On the other hand, the general Heun polynomials for the first excited mode n = 1 are given by  where the eigenvalues q 1;s must obey whose two solutions (s = 1, 2) are where ∆ = [γ(1 + η) + ηδ + ǫ] 2 + 4ηαβγ. It is worth noticing that the polynomial radial eigenfunctions for the first excited mode n = 1 are degenerate, since there exist two solutions for the eigenvalue q. However, we can choose one of these solutions by analyzing their behaviors. Indeed, the first solution (s = 1), and its corresponding eigenvalue (q 1;1 ) as well, is not suitable for our description and therefore it does not describe quasibound states, that is, it does not fully satisfy the two boundary conditions, since it has a false (nonphysical) singular point in the region between the exterior black hole event horizon and spatial infinity. In Fig. 3, we plot the first two squared QBS radial wave eigenfunctions. We observe that the radial solution tends to zero at spatial infinity (x → ∞ or r c = 4.86958, for a = −0.20, and Λ = 0.1) and diverges at the exterior black hole event horizon, which represents QBSs. Note that the radial wave eigenfunctions reach a maximum value at the exterior black hole event horizon (x = 1 or r + = 1.13694, for a = −0.20, and Λ = 0.1), and then cross this surface, as shown in the log-scale plots.

Angular wave eigenfunctions
The radial equation (14) and the angular equation (13) are coupled by a discrete set of angular eigenvalues λ, which are solutions due to the regularity requirements imposed to the angular functions at the two boundaries θ = 0 and θ = π.
For the QBSs, we can obtain an expression for the angular eigenvalues λ n;s (ω n ) from the polynomial equation (51) that determined the accessory parameter q n;s . Let us focus on the fundamental mode n = 0. Thus, from Eqs. (26) and (56), we obtain the following expression for the angular eigenvalues with a < 0: In this case, the angular eigenvalues λ 0;1 are real and positive, and hence the general degree ν 0;1 will be positive real number, which can be numerically evaluated from Eq. (61).
In Table I, we present the angular eigenvalues λ 0;1 , as well as the corresponding general degree ν 0;1 , as functions of the GB coupling constant a, for some values of the cosmological constant Λ. In addition, we also present the behavior of the QBS angular wave eigenfunctions P (θ) in Fig. 4, as functions of the new polar coordinate z = cos θ, for some values of the magnetic quantum number m. It is worth emphasizing that the two numerically satisfactory solutions of the associated Legendre equation of general degree are given in terms of the Ferrers function of the first kind P −m ν (−z) and P −m ν (z) in the interval −1 < x < 1. Note that these angular solutions are regular at the two boundaries θ = π (z = −1) and θ = 0 (z = 1).
In this case, we adopt the following parametrization for the metric function, where r c is the (positive) "exterior" de Sitter cosmological horizon, r + is the (positive) "exterior" black hole event horizon, r − is the (positive) "interior" black hole event horizon, and r −− is the (negative) "interior" de Sitter cosmological horizon. Then, with the metric function given by Eq. (62), it is easy to see that Eq. (14) is also a (general) Heun-type equation with four finite regular singularities (associated to the four event horizons) and one regular singularity at (spatial) infinity. Therefore, by following straightforwardly the steps described in the previous section, the analytical solution for the radial part of the conformally coupled massless Klein-Gordon equation, in the 4DEGBdS black hole spacetime with a > 0, can be written as where C 1,j and C 2,j are constants to be determined, and j = {x 4 , 0, 1, ∞} labels the solution at each singular point. Thus, the pair of linearly independent solutions corresponding to the exponents 0 and 1 − ǫ at x = x 4 (r = r −− ) is given by Similarly, the pair of linearly independent solutions at x = 0 (r = r − ) is given by Furthermore, the pair of linearly independent solutions corresponding to the exponents 0 and 1 − δ at x = 1 (r = r + ) is given by Finally, the pair of linearly independent solutions corresponding to the exponents α and β at x = ∞ (r = r c ) is given by Here, the new radial coordinate, x, is defined as In addition, we set a new parameter, x 4 , which is associated with the four finite regular singularities through the relation The parameters α, β, γ, δ, ǫ, and q are given by

Hawking radiation
In the limit when r → r + , which implies that x → 1, the radial solution given by Eq. (63) behaves as where and Therefore, these wave solutions give a thermal spectrum analogous to the spectrum of black body radiation.
In order to satisfy purely ingoing boundary conditions, from the asymptotic behavior of the radial solution at the exterior black hole event horizon described by Eq. (81), we must impose that C 2,1 = 0.
On the other hand, in the limit when r → r c , which implies that x → ∞, the radial solution given by Eq. (63) behaves as where with In the present case, by imposing the condition given by Eq. (50), with the parameter α given by Eq. (75), we obtain a first order equation for ω, whose solution is the exact spectrum of QBSs: In Table I we present the QBSs ω n as functions of the GB coupling constant, a, for some values of the cosmological constant Λ. In addition, we also present the behavior of the QBSs ω n in Fig. 5. We can conclude that the spectrum given by Eq. (89) is physically admissible, since the real part of σ is positive, and therefore it represents QBS frequencies for conformally coupled massless scalars in the 4DEGBdS black hole spacetime with a > 0. Furthermore, we can see that the decay is overdamped, since this quasispectrum is purely imaginary.

Radial wave eigenfunctions
In the present case, the QBS radial wave eigenfunctions for conformally coupled massless scalars propagating in a 4DEGBdS black hole spacetime with a > 0 are given by U n;s (x) = C n;s x γ−1 2 (x − η) HeunGp n;s (x 4 , q n,s ; −n, β, γ, δ; x).
In Fig. 6, we plot the first two squared QBS radial wave eigenfunctions. We observe that the radial solution tends to zero at spatial infinity (x → ∞ or r c = 4.899130, for a = 0.20, and Λ = 0.1) and diverges at the exterior black hole event horizon, which represents QBSs. Note that the radial wave eigenfunctions reach a maximum value at the exterior black hole event horizon (x = 1 or r + = 0.916563, for a = 0.20, and Λ = 0.1), and then cross this surface, as shown in the log-scale plots.

Angular wave eigenfunctions
In the present case, from Eqs. (80) and (56), we obtain the following expression for the fundamental mode angular eigenvalues with a > 0: In this case, the fundamental mode angular eigenvalues λ 0;1 are real and negative, and hence the general degree ν 0;1 will be complex, which can be numerically evaluated from Eq. (91). In fact, one of the first application of complex angular momentum techniques to atomic and molecular scattering was reported in the end of the 1970s [31], where it was evaluated the Legendre functions of the first kind of complex degree ν.
In Table I, we present the fundamental mode angular eigenvalues λ 0;1 , as well as the corresponding general degree ν 0;1 , as functions of the GB coupling constant a, for some values of the cosmological constant Λ. In addition, we also present the behavior of the QBS angular wave eigenfunctions P (θ) in Fig. 7, as functions of the new polar coordinate z = cos θ, for some values of the magnetic quantum number m. Note that these angular solutions are regular at the two boundaries θ = π (z = −1) and θ = 0 (z = 1).

V. THE 4DEGB-ANTI-DE SITTER BLACK HOLE SPACETIME
For a 4DEGBAdS black hole spacetime, we will use the relation Λ/3 = −1/R 2 , where R is the AdS cosmological radius, and then fix R = 1, measure the black hole (interior and exterior) radius in units of R, and calculate the corresponding value of the total mass M for each fixed value of the GB coupling constant in the stability region −0.20 ≤ a ≤ 0.20.
We show the behavior of the metric function f (r) in Fig. 8. In addition, the horizons are shown in the inlay plots as functions of the AdS cosmological radius R. We can see that the surface equation f (r) = 0 has one real solution when a < 0, otherwise, there are two real solutions when a > 0.
A. (AdS) Case 1: a < 0 In this case, we adopt the following parametrization for the metric function, where r + is the (positive) "exterior" black hole event horizon. Then, with the metric function given by Eq. (92), we will show that Eq. (14) is totally convenient to study QBSs with purely ingoing boundary conditions at the exterior black hole event horizon and vanishing boundary conditions at infinity, since it is a (biconfluent) Heun-type equation [28] with one finite regular singularity (associated to the exterior black hole event horizon) and one irregular singularity of rank 2 at (spatial) infinity. Now, we also follow the steps described in the VBK approach [16,17], as well as in the VB description of the biconfluent Heun functions [32], to obtain the analytical solution of Eq. (14) in the 4DEGBAdS black hole spacetime with a < 0 (without the assumption of specific boundary conditions). First of all, we need to define a new radial coordinate, x, as where κ(= −1 1/4 /ω 1/2 ) is a (dimensionless) constant. This definition moves the singularity r + to the point 0. The next step is to perform a special case of the s-homotopic transformation U (x) → y(x) given by where Thus, by substituting Eqs. (92)-(96) into Eq. (14), we get where γ = 2ir + ω, Equation (97) has the (canonical) form of a biconfluent Heun equation, where y(x) ≡ HeunB(α, β, γ, δ; x) denotes the biconfluent Heun function, which is the solution corresponding to the exponent 0 at x = 0 and assumes the value 1 there. If α is not a relative integer number, then from the Fuchs-Frobenius theory it follows that the HeunB(α, β, γ, δ; x) exists, and can be expanded as with ρ(ρ + α) = 0, where the indicial equation has two solutions, namely, ρ = 0 and ρ = −α. Let us focus in the ρ = 0 solution. Thus, the HeunB(α, β, γ, δ; x) is a entire function and given by with A 0 = c 0 = 1 , where (1 + α) n = Γ(1 + α + n)/Γ(1 + α). On the other hand, for ρ = −α, the solution of Eq. (97) corresponding to the exponent −α at x = 0 is x −α HeunB(−α, β, γ, δ; x). In addition, when α is a negative integer number, such that α = −m with m ≥ 1, it is possible to define the function HeunB(−m, β, γ, δ; x) = x m HeunB(m, β, γ, δ; x). Therefore, the analytical solution for the radial part of the conformally coupled massless Klein-Gordon equation, in the 4DEGBAdS black hole spacetime with a < 0, can be written as where C 1,j and C 2,j are constants to be determined, and j = {0, ∞} labels the solution at each singular point. Thus, the pair of linearly independent solutions at x = 0 (r = r + ) is given by y 1,0 = HeunB(α, β, γ, δ; x), (105) y 2,0 = x −α HeunB(−α, β, γ, δ; x).
Furthermore, the solutions of Eq. (97) at infinity (x → ∞) have the following asymptotic developments where | arg x| ≤ π/2 − ǫ. Next, we also use the VBK approach [16,17] to derive the characteristic resonance equation and then find the spectrum of resonant frequencies related to quasibound states by imposing two boundary conditions on the radial solution: it should describe an ingoing wave at the exterior black hole event horizon (r → r + ) and tend to zero far from the black hole at spatial infinity (r → ∞).

Hawking radiation
In the limit when r → r + , which implies that x → 0, the radial solution given by Eq. (104) behaves as where with and Therefore, these wave solutions give a thermal spectrum analogous to the spectrum of black body radiation.

Quasibound states
In order to satisfy purely ingoing boundary conditions, from the asymptotic behavior of the radial solution at the exterior black hole event horizon described by Eq. (109), we must impose that C 2,0 = 0.
On the other hand, in the limit when r → ∞, which implies that x → ∞, the radial solution given by Eq. (104) behaves as In the present case, it is easy to see that the radial solution tends to zero far from the black hole at spatial infinity, since the exponential (negative) term dominates its asymptotic behavior. However, the final asymptotic behavior of the radial solution at spatial infinity will be determined when we know the values of the parameters β and γ, which depends on the frequencies ω. In order to determine this, we use a polynomial condition for the biconfluent Heun functions to match the two asymptotic solutions of the scalar radial equation in their common overlap region. It is known that the biconfluent Heun function becomes a polynomial of degree N if and only if the following two conditions are satisfied [28]: where N = 0, 1, 2, . . . is the principal quantum number. Equation (115) will be called as the γ condition, which provides the expression for the frequencies ω N . On the other hand, the parameter δ must be appropriately chosen so that Eqs. (101) and (103) where D ν (z) are the parabolic cylinder (or Weber-Hermite) functions [33][34][35]. Thus, by imposing the γ condition given by Eq. (115), with the parameters α and γ given by Eqs. (95) and (98), respectively, we obtain a first order equation for ω, whose solution is the exact spectrum of QBSs: In Table II we present the QBSs ω N as functions of the exterior black hole event horizon r + . In fact, the exterior black hole event horizon r + is measured in terms of the AdS cosmological radius R, such that it is a fixed parameter In addition, we also present the behavior of the QBSs ω N in Fig. 9. We can conclude that the spectrum given by Eq. (118) is physically admissible, and therefore it represents QBS frequencies for conformally coupled massless scalars in the 4DEGBdS black hole spacetime with a < 0. Furthermore, we can see that the decay is overdamped, since this quasispectrum is purely imaginary.

Radial wave eigenfunctions
The radial wave eigenfunctions, which are related to QBSs of conformally coupled massless scalars propagating in the 4DEGBAdS black hole spacetime with a < 0, can be obtained by using the δ condition given by Eq. (116), as well as by imposing the γ condition given by Eq. (115) to such a parameter, as it is described in the VBK approach [16,17].
As we explained before, these polynomial radial eigenfunctions are related to the appropriate determination of the eigenvalue δ. Since it is calculated via Eq. (116), we index its solutions by a parameter s (= 1, 2, 3, . . .), which can be denoted by δ N ;s . Thus, the corresponding biconfluent Heun polynomials are now denoted as HeunBp N ;s (x). Therefore, the QBS radial wave eigenfunctions for conformally coupled massless scalars propagating in a 4DEGBAdS black hole spacetime with a < 0 are given by where C N ;s is a constant to be determined. Next, we calculate the biconfluent Heun polynomials related to the fundamental and first excited modes. In the present case, we have that α = −N − 1, and hence we will write the biconfluent Heun polynomials by using Eq. (117).
The biconfluent Heun polynomial for the fundamental mode N = 0 is given by where the eigenvalue δ 0;1 must obey whose unique solution (s = 1) is On the other hand, the biconfluent Heun polynomials for the first excited mode N = 1 are given by where the eigenvalues δ 1;s must obey whose two solutions (s = 1, 2) are It is worth noticing that the polynomial radial eigenfunctions for the first excited mode N = 1 are degenerate, since there exist two solutions for the eigenvalue δ. Here, we choose the first solution s = 1, and hence the eigenvalue δ 1;1 . In Fig. 10, we plot the first two squared QBS radial wave eigenfunctions. We observe that the radial solution tends to zero at spatial infinity (x → ∞) and diverges at the exterior black hole event horizon, which represents QBSs. Note that the radial wave eigenfunctions reach a maximum value at the exterior black hole event horizon (x = 0 or r + = 0.5), and then cross this surface, as shown in the log-scale plots.

Angular wave eigenfunctions
In the present case, from Eqs. (99) and (122), we obtain the following expression for the fundamental mode angular eigenvalues with a < 0: In this case, the fundamental mode angular eigenvalues λ 0;1 are real and negative, and hence the general degree ν 0;1 will be complex, which can be numerically evaluated from Eq. (127).
In Table II, we present the fundamental mode angular eigenvalues λ 0;1 , as well as the corresponding general degree ν 0;1 , as functions of the exterior black hole event horizon r + . In addition, we also present the behavior of the QBS angular wave eigenfunctions P (θ) in Fig. 11, as functions of the new polar coordinate z = cos θ, for some values of the magnetic quantum number m. Note that these angular solutions are regular at the two boundaries θ = π (z = −1) and θ = 0 (z = 1).

B. (AdS) Case 2: a > 0
In this case, we adopt the following parametrization for the metric function f (r), where r + is the (positive) "exterior" black hole event horizon, and r − is the (positive) "interior" black hole event horizon. Then, with the metric function given by Eq. (128), we will show that Eq. (14) is totally suitable to study QBSs with purely ingoing boundary conditions at the exterior black hole event horizon and vanishing boundary conditions at infinity, since it is a (confluent) Heun-type equation [28] with two finite regular singularities (associated to the black hole event horizons) and one irregular singularity of rank 1 at (spatial) infinity. Now, we also follow the steps described in the VBK approach [16,17], to obtain the analytical solution of Eq. (14) in the 4DEGBAdS black hole spacetime with a > 0 (without the assumption of specific boundary conditions). First of all, we need to define a new radial coordinate, x, as This definition move the singularities (r + ,r − ) to the points (0,1). The next step is to perform a special case of the s-homotopic transformation U (x) → y(x) given by where Thus, by substituting Eqs. (128)-(132) into Eq. (14), we get where with δ = 2(r 2 Equation (134) has the (canonical) form of a confluent Heun equation, where y(x) ≡ HeunC(α, β, γ, δ, η; x) denotes the confluent Heun function, which is the solution corresponding to the exponent 0 at x = 0 and assumes the value 1 there. If β is not a relative integer number, then this standard confluent Heun function is defined via the convergent Taylor series expansion in the disc |x| < 1, where it is assumed the normalization HeunC(α, β, γ, δ, η; 0) = 1. The coefficients c n are given by with the initial conditions c −1 = 0 and c 0 = 1. The expressions for P n , T n , and X n are given by Therefore, the analytical solution for the radial part of the conformally coupled massless Klein-Gordon equation, in the 4DEGBAdS black hole spacetime with a > 0, can be written as where C 1,j and C 2,j are constants to be determined, and j = {0, 1, ∞} labels the solution at each singular point. Thus, the pair of linearly independent solutions at x = 0 (r = r − ) is given by where the parameters α, β, γ, δ, and η are given by Eqs. (131), (132), (133), (137), and (138), respectively. Furthermore, by following the same steps, the pair of linearly independent solutions at x = 1 (r = r + ) is given by where the new parameters α, β, γ, δ, and η are given by Finally, the solutions of Eq. (134) at infinity (x → ∞) have the following asymptotic developments Next, we also use the VBK approach [16,17] to derive the characteristic resonance equation and then find the spectrum of resonant frequencies related to quasibound states by imposing two boundary conditions on the radial solution: it should describe an ingoing wave at the exterior black hole event horizon (r → r + ) and tend to zero far from the black hole at spatial infinity (r → ∞).

Hawking radiation
In the limit when r → r + , which implies that x → 1, the radial solution given by Eq. (144) behaves as where with and Therefore, these wave solutions give a thermal spectrum analogous to the spectrum of black body radiation.
In order to satisfy purely ingoing boundary conditions, from the asymptotic behavior of the radial solution at the exterior black hole event horizon described by Eq. (156), we must impose that C 2,1 = 0.
On the other hand, in the limit when r → ∞, which implies that x → ∞, the radial solution given by Eq. (144) behaves as where In the present case, when ω = −iω I , it is easy to see that the radial solution tends to zero far from the black hole at spatial infinity, since p > 0 and q > 0. Thus, in order to determine the frequencies ω, we use a polynomial condition for the confluent Heun functions to match the two asymptotic solutions of the scalar radial equation in their common overlap region. It is known that the confluent Heun function becomes a polynomial of degree N if and only if the following two conditions are satisfied [28]: where N = 0, 1, 2, . . . is the principal quantum number. Equation (164) will be called as the δ condition, which provides the expression for the frequencies ω N . On the other hand, the parameter ξ must be appropriately chosen so that Eqs. (140)-(143) are consistent, which means that it is also necessary for the parameter ξ to be an eigenvalue of the confluent Heun equation, calculated via Eq. (165). Furthermore, the parameter ξ, that is, the parameter η, contains the separation constant λ, which indicates that we could obtain the eigenvalues of the separation constant λ N (ω N ), corresponding to the appropriate eigenvalues of the frequencies ω N , from the polynomial solution for the radial equation and then use it to show the (regular) angular behavior of massless scalar QBSs in the background under consideration. Note that the δ condition is equivalent to X N +2 = 0 in Eq. (143), and the ∆ condition is equivalent to the requirement c N +1 = 0 in Eq. (140). Thus, by imposing the δ condition given by Eq. (164), we obtain a first order equation for ω, whose solution is the exact spectrum of QBSs: In Table II we present the QBSs ω N as functions of the exterior black hole event horizon r + for some values of the GB coupling constant a. In fact, the exterior black hole event horizon r + is measured in terms of the AdS cosmological radius R, such that it is a fixed parameter which does not depend on the GB coupling constant a. On the other hand, the interior black hole event horizon r − depends on the GB coupling constant a.
In addition, we also present the behavior of the QBSs ω N in Fig. 12. We can conclude that the spectrum given by Eq. (166) is physically admissible, and therefore it represents QBS frequencies for conformally coupled massless scalars in the 4DEGBdS black hole spacetime with a < 0. Furthermore, we can see that the decay is overdamped, since this quasispectrum is purely imaginary.

Radial wave eigenfunctions
The radial wave eigenfunctions, which are related to QBSs of conformally coupled massless scalars propagating in the 4DEGBAdS black hole spacetime with a > 0, can be obtained by using the δ condition given by Eq. (164), as well as by imposing the ∆ condition given by Eq. (165) to parameter ξ, as it is described in the VBK approach [16,17].
As we explained before, these polynomial radial eigenfunctions are related to the appropriate determination of the eigenvalue ξ. Since it is calculated via Eq. (165), we index its solutions by a parameter s(= 1, 2, 3, . . .), which can be conveniently denoted by ξ N ;s . Thus, the corresponding confluent Heun polynomials are now denoted as HeunCp N ;s (x). Therefore, the QBS radial wave eigenfunctions for conformally coupled massless scalars propagating in a 4DEGBAdS black hole spacetime with a > 0 are given by where C N ;s is a constant to be determined. Next, we calculate the confluent Heun polynomials related to the fundamental and first excited modes. The confluent Heun polynomial for the fundamental mode N = 0 is given by where the eigenvalue ξ 0;1 must obey whose unique solution (s = 1) is On the other hand, the confluent Heun polynomials for the first excited mode N = 1 are given by where the eigenvalues ξ 1;s must obey whose two solutions (s = 1, 2) are with ∆ = (α − 2 − β − γ) 2 + 4α(1 + β). It is worth noticing that the polynomial radial eigenfunctions for the first excited mode N = 1 are degenerate, since there exist two solutions for the eigenvalue ξ. Here, we choose the first solution s = 1, and hence the eigenvalue ξ 1;1 . In Fig. 13, we plot the first two squared QBS radial wave eigenfunctions. We observe that the radial solution tends to zero at spatial infinity (x → ∞) and diverges at the exterior black hole event horizon, which represents QBSs. Note that the radial wave eigenfunctions reach a maximum value at the exterior black hole event horizon (x = 0 or r + = 0.5), and then cross this surface, as shown in the log-scale plots.

Angular wave eigenfunctions
In the present case, from Eqs. (135) and (170), we obtain the following expression for the fundamental mode angular eigenvalues with a > 0: In this case, the fundamental mode angular eigenvalues λ 0;1 are real and positive, and hence the general degree ν 0;1 will be real, which can be numerically evaluated from Eq. (175).
In Table II, we present the fundamental mode angular eigenvalues λ 0;1 , as well as the corresponding general degree ν 0;1 , as functions of the exterior black hole event horizon r + . In addition, we also present the behavior of the QBS angular wave eigenfunctions P (θ) in Fig. 14, as functions of the new polar coordinate z = cos θ, for some values of the magnetic quantum number m. Note that these angular solutions are regular at the two boundaries θ = π (z = −1) and θ = 0 (z = 1).

VI. FINAL REMARKS
In this work, we obtained the exact analytical general solutions of the conformally coupled massless Klein-Gordon equation in the 4DEGB de Sitter and Anti-de Sitter black hole spacetimes. The general solution of the angular part is given in terms of the associated Legendre functions with arbitrary degree, while the general solutions for the radial part are given in terms of the Heun functions.
On the radial solutions, we imposed the two boundary conditions related to the quasibound state phenomena. Near the exterior event horizon, the radial solutions describe the Hawking radiation spectrum. In particular, the ingoing waves reach a maximum value and then crosses into the black hole. On the other hand, far from the black hole at asymptotic (spatial) infinity, the radial solutions tend to zero, that is, the probability of finding any scalar particles in the spatial infinity is null.
The spectrum of quasibound states for massless scalar fields was obtained by using the polynomial condition of the Heun functions. In fact, that is a new (analytical) approach developed by Vieira, Bezerra, and Kokkotas [16,17].
Finally, we can discuss the stability of the systems. All the systems are stables, and present an over-damped motion, since their resonant frequencies is purely imaginary and negative. In addition, it is worth pointing out that these quasibound state frequencies were obtained directly from the Heun functions by using a polynomial condition, and, to our knowledge, there is no similar result in the literature for the backgrounds under consideration.
We hope that our results, which describe an unquestionably phenomenon associated with purely quantum effects in gravity, may be used to fit some astrophysical data in the near future and hence shed some light on the physics of black holes, and alternative/modified theories of gravity as well.

DATA AVAILABILITY
The data that support the findings of this study are available from the corresponding author upon reasonable request.