1. Introduction
Black holes are known by their mass, angular momentum and their charge. Choptuik [
1] illustrates that there may be another parameter explaining the gravitational collapse solutions. Christodolou [
2] elaborated on the spherically symmetric collapse of the scalar field. Choptuik [
1] empirically showed the critical behaviour of the solutions that illustrated the discrete self-similarity for the real scalar field gravitational collapse. The gravitational collapse solution shows space-time self-similarity so that dilations occur. The critical solution provides the scaling law. We can illustrate the initial condition of the real scalar field by parameter
p, which is related to the field amplitude. Let
define the critical solution. The black hole is then formed when
p gets bigger than
. If
, then one finds out the mass of the black hole or the Schwarzschild radius by a scaling law as follows
The critical exponent for a real scalar field in four dimensions is found to be
[
1,
3,
4], while for the other dimensions (
), the mass of black hole is given by [
5,
6] as
The numerical optimization for different matter content can be read from [
7,
8,
9,
10,
11,
12]. The collapse solutions of the perfect fluid had been studied in [
5,
13,
14,
15] and its critical exponent
was also found in [
14]. It has also been discussed in [
16] that
may have a universal value for all the fields that are coupled to gravity in four dimensions. One method for obtaining the critical exponent is based on perturbations of self-similar solutions. This method was established in [
5,
15,
17]. Let
h be any field, its perturbation is given by the following
where
has the scaling
which is related to different modes. The most relevant mode
is related to the biggest value of
. The minus sign demonstrates a growing mode near the black hole formation time
. It was proved in [
5,
15,
17] that
can be related to Choptuik exponent by the following relation
The axial symmetry was studied in [
18]. On the other hand, Alvarez-Gaume et al. [
19] explore other solutions as well as shock waves. The critical exponent
for the axion-dilaton system in four dimensions was first discovered in [
20]. The analysis of [
3,
21] had been done in four dimensions for the elliptic case; however, their techniques can be used in hyperbolic and parabolic cases as well as in other dimensions, for further results see [
22,
23].
As the first motivation for the critical collapse of the axion-dilaton system, one can consider the AdS/CFT correspondence [
24,
25,
26,
27], relating the Choptuik exponent, the imaginary part of quasinormal modes and also the dual conformal field theory [
28]. In the context of AdS/CFT, it is desirable to consider the spaces which tend asymptotically to
. From the bosonic fields of the theory, one can employ the system of axion-dilaton and the self-dual 5-form to analyze the critical collapse solutions in different dimensions. The other motivation of the axion-dilaton system is actually related to the holographic description of black hole formation [
13,
29]. The next motivation could potentially be the application of the system to black hole physics [
30,
31] as well as the important role of S-duality in self similar solutions [
32].
The families of continuous self-similar solutions were explored by [
33] in four and five dimensions for all classes of SL(2,R). These solutions are the extended results of [
34,
35]. More specifically, [
36] studied the perturbations and reproduced the existing value [
20] of
in four dimensions. Some other critical exponents for four and five dimensions were proved in [
23]. Hatefi and Hatefi [
37] proposed the Fourier-based regression models for the nonlinear critical collapse functions. To address the challenges of [
37], Hatefi and Hatefi [
38] employed truncated power basis, natural spline and penalized B-spline regression models to estimate the non-linear critical collapse functions. Recently, Hatefi et al. [
39] applied artificial neural networks to analyze the black hole solutions in parabolic class in higher dimensions in detail.
In this paper, we propose a new methodology, which allows us to approach the problem of the complexity of elliptical black holes in 4 dimensions by means of Bayesian statistical modeling. Technically, we address the problem of modeling the self-similar solutions for the described spherical gravitational collapse [
3,
21] in four dimensions, using Hamiltonian Monte Carlo with Stacked Neural Networks. The modeling of black holes in the 4-dimensional elliptic case is reduced to a set of differential equations, which have always been approached by numerical methods that try to find a unique solution that satisfies the initial conditions. On the contrary, we argue that it is possible to attack the problem with Bayesian statistical modeling of the differential equations. In order to simplify the equations and the parameters of the equations of motion, researchers often have to use a variety of numerical approaches due to the complicated and highly nonlinear form of the equations of motion in the physics of black holes (e.g. [
22,
38,
39]). For this reason, while calculating the parameters of the equations of motion, they must disregard the measurement errors imposed by the numerical calculations.
To the best of our knowledge, for the first time in the literature on the axion-dilation black hole, we investigate the equations of motion of the critical collapse functions within a Bayesian framework. This Bayesian framework enables us to address the measurement errors imposed by the numerical calculations in the estimation process. The parameters of the equations of motion are handled as random variables due to this statistical modelling. Because of this, we can derive the posterior distributions for these parameters, which we can then utilize as legitimate priors for a collection of stacked neural network-based solvers that will calculate the critical collapse function. In a nutshell, this paper describes a novel pipeline that starts solving the elliptic black hole equations in 4 dimensions using Hamiltonian Monte Carlo statistical modelling. Then, the posterior distributions discovered, for the parameters of the equations, are used with stacked neural networks to estimate the critical collapse functions. Unlike methods in the literature, this probabilistic perspective enables us to recover the deterministic solution in the literature and find all possible solutions that may occur due to numerical measurement errors in the parameter estimation phase. This paper paves the path to view the problem from a stochastic rather than a deterministic perspective.
This paper is organized as follows.
Section 2 discusses the effective action of the axion-dilaton system in four dimensions, its equations of motion, and the initial conditions for the axion-dilaton system that come from the requirement of continuous self-similarity.
Section 3 describes our novel approach for modelling the complexity of elliptic black hole solution in 4d from Hamiltonian Monte Carlo with stacked neural networks. Through extensive numerical studies, we develop the Bayesian posterior mean and median neural networks and the 95% credible intervals for all the critical functions in the elliptic class 4d in
Section 4. Finally, in
Section 5, we present the summary and concluding remarks.
2. The Relevant System
The axion-dilaton fields are combined into a single complex scalar field
. The four-dimensional axion-dilaton (
) action is given by
where
R is the scalar curvature. By taking variations from the metric and
, one finds all the equations of motion as follows
Note that the theory is classically invariant under SL(2,R) transformations, that is, if
gets exchanged by
where
,
and
as well as the action remain to be invariant. This group can be broken to SL(2,Z) as a consequence of duality in String Theory (for review of duality see for instance [
40,
41,
42,
43]).
The spherically symmetric metric is given by [
21] as
One can implement time re-definitions for (
9), hence as in [
21] we set
and regularity condition implies
. The continuous self-similarity means the existence of a killing vector
that produces global scale transformation as follows
In spherical coordinates, one chooses . Now by defining the scale invariant variable , self-similarity of the metric reflects the fact that all the functions can be expressed in terms of z.
The action in (
5) is SL(2,R)-invariant, so one compensates a scale transformation of
by an SL(2,R) transformation. Hence, by considering the change of variables to
, one derives a differential condition (as originally pointed out in [
34]) for
as
where all
are real numbers. The above equation has two roots which can be either two complex conjugate numbers, two real numbers or a double real root where they are related to compensating the scaling transformation in SL(2,R). In the elliptic class, the general form of the ansatz is given by
where under a scaling transformation
,
changes by a SL(2,R) rotation, which means that all equations are invariant under the following transformation
where
is a complex function that satisfies |
f(
z)| < 1 to be determined by solutions to non-linear ordinary differential equations, and
is a real constant.
In the following, we describe briefly the derivation of the equations of motion for the elliptic class in four dimensions. Applying continuous self-similarity ansätze (
12) to all the equations of motion (
6) and (
7), one finds the ordinary differential equations for
,
,
. By taking the spherical symmetry,
can be expressed in terms of
and
as
Note that the first derivative of
can be eliminated from all equations of motion by the following constraint
Hence, the equations of motion in the elliptic class in four dimensions are then given by
The first equation of motion includes
where
b is a first-order linear in-homogeneous equation where its initial condition is given by
. The initial conditions for
are known by applying the smoothness of the critical solution. Applying polar coordinate representation
, one can easily show that all equations are invariant under a global redefinition of the phase of
. This implies that
. From the regularity condition at
and the residual symmetries in the equations of motions (
17), one can derive the initial boundary conditions of the equations of motion by
4. Numerical Studies
In this section, we estimate the critical collapse functions of the equations of motion in the elliptic class of 4 dimensions in a Bayesian framework. In this framework, the parameters of the equations are treated as random variables to account for the measurement error involved in numerically estimating the parameters of the equations. We also explore the effect of this complexity on the estimation of the critical collapse functions.
Various research articles in the literature [
36,
39] investigated the non-linear equations of motion of the axion-dilaton system in different dimensions and for different ansatz. Hatefi et al. [
37,
38] used the properties of the Fourier-based and spline nonlinear smothers in predicting the critical collapse functions, respectively. Recently, Hatefi et al. [
39] developed an NN-based solver to show that there is no solution in the parabolic class of black holes and [
36] used a grid search root-finding method to estimate the parameters of the equations of motion numerically. According to the complex and highly nonlinear pattern of the equations, all these methods in the literature use almost a series of numerical techniques to make the equations tractable and approximate the parameters of the models. These techniques include, for instance, imposing some constraints on the equations of motions such as finiteness of
as
. Another constraint that one employs is the vanishing of the divergent part of
that produces a complex-valued constraint at
. As the equations are not tractable analytically in elliptic equations of motion, for example, [
33] used a numerical profile-root finding method and applied numerically a grid search discrete optimization on the coordinates of the extended parameter space of the equations of motion. To do so, the equations are simplified by approximating the target functions by the first two orders of the Taylor expansions. Applying a grid search optimization on the coordinates of the extended parameter space in the equation of motion, [
33] could estimate the required parameters of the equations of motion. Removing the spurious roots of the equations, they used parameter estimates to find the solution to the equations of motion and estimate the critical collapse functions.
Unlike [
36,
37,
38,
39], to our best knowledge, for the first time in the literature, through a probabilistic viewpoint, we develop a Bayesian method to model the measurement errors into the estimation of the critical collapse functions. Hence the framework provides a robust approach against measurement errors in parameter estimation. These Bayesian estimates incorporate information from all possible parameters’ outcomes, that may occur in the numerical experiments, from the posterior distribution and accommodate this complexity in estimating critical collapse functions.
Here all the critical collapse functions are treated as DE variables of DE equations of motion
16 and
17. Since all the DE variables are numerically solved simultaneously, it is reasonable to assume that the observed DE variables have the same standard deviation
parameter representing the variability in the numerical experiments. Therefore
will be the set of all unknown parameters of the likelihood function under equations of motion. We first use the HMC sampler, described in Sub
Section 3.2, to find the posterior distribution of the likelihood parameter
. We used the Python package Pymc3 [
52] to implement the HMC for the equations of motion. To evaluate the effect of the prior information of likelihood parameters in estimating the critical functions, we considered two different prior distributions for the parameter
. We set the first prior distribution as Gaussian distribution with a mean
and standard deviation
; that is
. In contrast to the unimodal Gaussian, we assumed the second prior followed non-informative Uniform distribution between (1,1.5); that
. We also used the half-Cauchy distribution with scale parameter
as the prior distribution for
to incorporate the uncertainty of the observed DE variables in the likelihood function.
Figure 1 and
Figure 6 show the posterior distributions
and
under the Gaussian and Uniform prior distributions, respectively. Comparing
Figure 1 and
Figure 6, under both prior distributions, we observe that the posterior distribution of
converges to a roughly uni-modal, bell-shaped distribution where the 94% highest density interval (HDI) roughly is (1.165,1.205). In other words, considering the measurement errors from numerical analysis of the parameter estimation, the true value of
in the equations of motion, with probability 94%, will vary between (1.165,1.205). Interestingly this is compatible with the solution for the elliptic case for four dimensions, obtained using various numerical experiments [
21,
34], which as reported as
. We see that the 94% HDI contains (under both prior distributions)
as one possible outcome of the numerical analysis. It can be concluded that given the numerical measurement errors involved in the DE variables, the probability that the true value of
appears
is roughly 13%; that is
.
In the second part of the analysis, we incorporated the posterior information of the parameters into the NNs, described in Sub
Section 3.3, to obtain the NN-based posterior mean and posterior median estimators of the critical collapse functions
and
. To do so, we randomly took
realizations
from the posterior distribution
. We then applied these 200 realizations iteratively and obtained the NN-based estimates of
and
given that
at
equally spaced
in the domain of the equations of motion. These NN-based estimators are denoted by
,
and
for
and
. In each iteration, we implement the NNs solver using the Python Package NeuroDiffEq [
53] configured with four hidden layers, each of 16 nodes, following [
39]. We ran the solvers for 1000 epochs using the initial boundary conditions in Eq.
18.
Before establishing the Bayesian stacked NN estimators, we would like to explain the probabilistic perspective in sampling from the posterior distribution obtained. The posterior distribution
lists all possible values of
and how often they are observed as the posterior estimate of the parameter of the equations of motion. The posterior candidate
, as a sample from the posterior distribution, changes from one sample to another. Accordingly, given the posterior candidate
, the posterior NN candidates, as estimators for the critical collapse functions, will change from one sample to another. The probability that each of these posterior NN candidates is observed in estimating the true critical collapse functions corresponds to the posterior distribution
. To show this probabilistic perspective,
Figure 2 and
Figure 7 represent 10 posterior NN candidates of the critical collapse functions based on 10 randomly selected
from
realizations of
using the Gaussian and Uniform prior distributions, respectively.
Because we have
L Bayesian NNs candidates corresponding to
L realizations from the posterior distribution, we applied the Bayesian model averaging. We obtained the two mean and median stacked NNs for all the critical functions at each space-time
. In addition to posterior mean and median estimates, we also computed the 95% Bayesian credible intervals for the critical functions on the entire domain. To compute the Bayesian 95% credible intervals, we sorted the NN estimates from the smallest value to the largest value at each space-time. We then computed the 2.5 and 97.5 percentiles of the estimates as lower and upper bands of the interval, respectively. To compare our method with the non-Bayesian approach of [
39], we computed the NNs-based estimates of the critical collapse functions similar to [
39], treating
with the same configuration as described above. This estimate is, henceforth, called Fixed NNs.
Figure 3 and
Figure 8 show the Bayesian posterior mean and median estimates and their 95% credible intervals using Gaussian and Uniform prior distributions, respectively. Moreover,
Figure 3 and
Figure 8 compare the Bayesian estimates with their corresponding fixed NN-based estimates of the critical functions. One can conclude from the figures that, given all possible values for
in the equations of motion, the true critical collapse functions, with probability 0.95, should be in the reported Bayesian credible intervals. Also, we observe that all 95% credible intervals contain the Fixed NNs estimates of [
39] as one possible candidate for estimating the critical functions. We also observe that the posterior mean and median NN estimator represent the form and curvature of the critical collapse functions as the average and median of all the Bayesian NN-based candidates.
Note that the fixed NN-based method of [
39] is not robust against the variability and measurement error in calculating the parameter. In other words, the Fixed NN estimates will change slightly from
in estimating the parameter’s value. Unlike the fixed NNs, the Bayesian credible intervals, posterior mean and median NN estimators will remain robust against changes in measurement errors in estimating
w because these Bayesian estimates have already considered all possibilities of
w in the domain from the posterior distribution.
Figure 4 and
Figure 9 show the posterior mean and median stacked NNs and 95% credible intervals of all the critical collapse functions together under Gaussian and Uniform prior distributions, respectively. It is apparent from
Figure 4 and
Figure 9 that the NN-based Bayesian credible intervals confirm that there must be a single self-similar solution, on average, in the area
to the elliptic equations of motion in 4 dimensions where the imaginary and real parts of the
G function intersects, where the explicit form of the
G function for the elliptic case in 4 dimension was given in [
33]. Interestingly, this finding is compatible with the finding in the literature where [
36,
39] showed that there was a single self-similar solution at
when the parameter of the equation is estimated to be
.
In order to evaluate the convergence of the proposed NN-based estimates using the posterior distribution
, we computed the train and validation losses in the last epoch of the NN solvers for each
.
Figure 5 and
Figure 10 show the trace and box plots of the loss differences for
realizations from the posterior distribution
under Gaussian and Uniform prior distributions, respectively. We clearly observe that the loss differences, on average, are less than
which validates the convergence of NN-based estimates in estimating the critical collapse functions in the elliptic class of 4 dimensions.
5. Summary and Concluding Remarks
This paper proposes a Bayesian methodology to model the self-similar solutions for the spherical gravitational collapse for an elliptic case in four dimensions. In the literature on the axion-dilaton system, due to the complex and highly nonlinear equations of motion, one needs to employ various numerical techniques and constraints to make the equations of motion tractable. These constraints include, for instance, the finiteness of
as
and the vanishing of the divergent part of
producing a complex-valued constraint at
. Because the equations of motion are not tractable analytically, for instance, [
33] used a series of numerical methods to find the solutions to the equations. Their numerical profile root-finding includes, for example, approximating the main effects of the equations by the first two orders of the Taylor expansion, applying a numerical grid search discrete optimization on the coordinates of the extended parameter space of the equations and removing the spurious roots, to name a few. This enables them to estimate the parameters of the equations and then find self-similar black hole solutions. According to the complex form of the equations and the vital role of the parameters, researchers have to ignore the measurement errors imposed in finding the parameters through numerical techniques.
To fix this problem, we have proposed a novel pipeline to take the measurement errors into the statistical models in finding the solution to the elliptic equations of motion. To do so, unlike the methods in the literature, we developed a Bayesian estimation method where the parameter is treated as a random variable. We proposed Hamiltonian Monte Carlo to derive the posterior distribution of the parameter efficiently. The posterior distribution determines all the possible outcomes in estimating the parameter of the equations of motion. The posterior distribution interestingly shows that the deterministic solution available in the literature to the elliptic equations of motion, namely , is the true value for the parameter with a probability of almost 0.13. In addition, the posterior distribution lists all possible values for the parameter of elliptic equations of motion along with their probabilities.
Next, we utilized neural networks (NNs) to incorporate the posterior information in solving the equations of motion and estimating the critical collapse functions in the entire domain of the equations. From a probabilistic viewpoint, one can find the NN estimate of the critical function as an estimator, which would be the true trajectory of the critical collapse function whose probability distribution corresponds to the posterior distribution of the parameter. Accordingly, we applied a Bayesian model averaging strategy and found the posterior mean and median stacked NNs and the 95% Bayesian credible intervals in estimating the critical collapse functions.
Given all variability in estimating the parameter of the equations of motion, the developed Bayesian credible intervals will contain the true critical collapse functions with a probability of 0.95. Comparing our Bayesian method with the deterministic NN approach of [
39], the developed Bayesian credible intervals contain the deterministic estimate as one possible candidate in estimating the critical functions. Unlike the estimate of [
39], the Bayesian proposals will remain robust against measurement errors in estimating
because these Bayesian estimates have already considered all possibilities of the parameter in the domain of the posterior distribution. It is worth mentioning that our approach will help us to reveal the value and range of the critical exponent of axion-dilaton system in a wider range in different dimensions that we would like to pursue in the near future.