Holomorphic Distributions and Connectivity by Integral Curves of Distributions

Holomorphic Distributions and Connectivity
by Integral Curves of Distributions††This paper is a contribution to the
Special Issue on Algebra, Topology, and Dynamics in Interaction in honor of
Dmitry Fuchs. The full collection is available at
https://www.emis.de/journals/SIGMA/Fuchs.html

Vladimir A. ZORICH

V.A. Zorich

Division of Mathematical Analysis, Faculty of Mechanics and Mathematics,
Lomonosov Moscow State University, GSP-1, Leninskie Gory, Moscow, 119991,
Russia<EMAIL_ADDRESS>
Received July 19, 2019, in final form September 14, 2019; Published online
September 19, 2019

It is known that the classical Frobenius theorem on conditions of
integrability for distributions of planes can be extended to the case of
complex holomorphic distributions. We show that an alternative criterion for
integrability, namely, non-connectivity, discovered (or at least, marked and
explicitly formulated) by Carathéodory in relation to classical
thermodynamics, also admits a holomorphic formulation.

holomorphic distribution; integral curve; connectivity; thermodynamic states;
adiabatic transitions; Carathéodory’s theorem

53D10; 32B15; 80A05

##  1 Introduction

We consider a holomorphic version of a theorem of Carathéodory initiated by
his studies related to mathematical aspects of classical thermodynamics. The
generalization which we are going to discuss, like the Carathéodory theorem
itself, is easily provable by modern tools of analysis and geometry. However,
an interesting, useful and important feature is rather a remarkable duality of
connectivity of the space of states of a thermodynamic system, or in any
abstract space, by admissible (tangent, contact) paths for a distribution of
hyperplanes in the space, and integrability of this distribution itself.

In more detail, in classical equilibrium thermodynamics one often has to
consider adiabatic transitions of a thermodynamic system from one
thermodynamic state to another. The adiabatic condition in mathematical
language means the requirement that the corresponding transition path goes
along the kernels (zero spaces) of a special 1-form, i.e., it is constrained
to be tangent to a certain distribution of hyperplanes. Justifying
(mathematically formalizing) the so-called second law of thermodynamics,
generating the entropy function of a thermodynamic state, Carathéodory [5, 6]
investigated the possibility of adiabatic transitions between thermodynamic
states. In mathematical terms he showed that the absence of local connectivity
of points in space by paths tangent to the distribution of hyperplanes is
equivalent to integrability of the distribution. Thus, the question of
connectivity (non-connectivity) turned out to be related with the issue of
integrability of the distribution given by a 1-form (and if we omit some
details, it turned out to be equivalent to it). But the problem of
integrability of a distribution was already solved by the well-known Frobenius
theorem. Thus, the problem of connectivity reduces to the already solved
problem.

Note that the first significant observations in this direction were already
made by Poincaré [12], also in connection with thermodynamics.

The reference [4] is a nice extended exposition by Born of Carathéodory’s
work, possibly more accessible to physicists. A mathematical generalization of
Carathéodory’s theorem for distributions of arbitrary codimension (not only to
distributions of hyperplanes) is now known as the Rashevskiy–Chow theorem [8,
13].

##  2 Integrability of distributions

Consider a smooth field (distribution) of hyperplanes in the space
${\mathbb{R}}^{n}$. Such a distribution can have integral surfaces, which are
surfaces whose tangent planes at each point coincide with the corresponding
planes of the distribution. If the space is locally foliated by such surfaces,
then such a distribution of hyperplanes is called an integrable distribution.

Generic vector fields are always integrable (this is one of fundamental
initial facts of the theory of ordinary differential equations). For
distributions of higher dimension this is not always the case. For hyperplane
distribution defined by kernels (zeros) of a non-degenerate 1-form $\omega$ a
necessary and sufficient condition of local integrability of the distribution
is the condition $\omega\wedge{\rm d}\omega\equiv 0$. This is the Frobenius
theorem [7], or rather one of its versions.

Below we will be interested in holomorphic distributions of complex
hyperplanes in the complex space. It can be verified that the local
integrability condition expressed by the Frobenius theorem is also applicable
when the distribution is given by kernels (zeros) of a holomorphic 1-form
$\omega$. (See, for example, [11] or proofs of the classical Frobenius theorem
in [7] or [3].)

##  3 Connectivity by paths tangent to distributions

Given a plane distribution in space, we are allowed to move from one to
another point of the space subject to the constraint that at each moment the
velocity vector of motion has to lie in the corresponding plane of the
distribution. In this case the path is said to be tangent (or subordinated) to
the distribution, or that it is an admissible path or an integral curve for
the distribution.

In mechanics this is a typical example of motion in the presence of
nonholonomic constraints. In thermodynamics it may be a problem of the
adiabatic transfer of a thermodynamic system from one equilibrium state to
another. In control theory one looks for (im)possibility of reaching one point
of space from another by means of a controlled admissible path for a given
distribution, determined by the problem itself.

If the distribution is integrable, then obviously no admissible path can leave
the integral surface of the distribution to which the starting point of the
path belongs. No point lying outside of the integral surface can be connected
by an admissible path with any point of this integral surface.

So, if there is integrability of a distribution of hyperplanes, then,
certainly, there is no local connectivity of points of space by admissible
paths. (For instance, in a thermodynamic system there cannot be an adiabatic
transition through equilibrium states between two equilibrium states having
different entropy values.)

A nontrivial observation of Carathéodory, initiated by his work on the
mathematical formalization of classical thermodynamics, is the above-mentioned
statement that in the case of a hyperplane distribution the inverse theorem is
valid: if there is no local connectivity, then there is local integrability.

In short, although not quite accurately, there is an alternative: either there
is local connectivity and there is no local integrability, or there is no
local connectivity, and then there is local integrability.

Since the question of the local integrability of the distribution of
hyperplanes is solved by the Frobenius theorem, Carathéodory’s theorem
provides an effective tool to verify the possibility of connecting points
(states of the system) by means of paths admissible for a given distribution.

##  4 Main statement

Now we can turn to the main subject of this note. We are going to show that
not only the holomorphic analogue of the Frobenius theorem is valid, which was
already mentioned above, but the holomorphic analogue of the Carathéodory
theorem takes place as well.

We will consider the simplest nontrivial case, the distribution of complex
planes in the complex space ${\mathbb{C}}^{3}$.

Let our distribution be given as zeros of the standard 1-form $x{\rm d}y+{\rm
d}z$. Instead of the form $x{\rm d}y+{\rm d}z$ at the moment it is technically
more convenient for us to consider the linearly equivalent form $x{\rm
d}y-y{\rm d}x-{\rm d}z$.

Let $x=t$, $y=t^{2}$, and $z=\frac{1}{3}t^{3}$. Since ${\rm d}z=x{\rm
d}y-y{\rm d}x$, we have an admissible holomorphic curve parametrized by the
complex parameter $t$. The function $z=\frac{1}{3}t^{3}$ for some value of
$t=t_{1}$ can take any prescribed value $z_{1}$. At the complex moment $t_{1}$
the curve $\big{(}x=t,y=t^{2}$, $z=\frac{1}{3}t^{3}\big{)}$ passes through the
point
$(x_{1},y_{1},z_{1})=\big{(}t_{1},t^{2}_{1},\frac{1}{3}t^{3}_{1}\big{)}$. By
an additional transition along the admissible complex line
$(x=x_{1}t,y=y_{1}t,z=z_{1})$ one can move to the point $(0,0,z_{1})$. But
once we managed to get from the point $(0,0,0)$ to any point of the form
$(0,0,z_{1})$ by an admissible holomorphic curve, then, using the same
construction, we can reach any point of the space ${\mathbb{C}}^{3}$.

Thus, in principle, it is possible to pass from any point of the space
${\mathbb{C}}^{3}$ to any other point along holomorphic curves subjected to
the distribution defined in ${\mathbb{C}}^{3}$ as $\\{\ker(x{\rm d}y-y{\rm
d}x-{\rm d}z)\\}$ or as $\\{\ker(x{\rm d}y+{\rm d}z)\\}$.

In the above discussion a holomorphic broken line has no more than three
links. One can carry out this reasoning more carefully, and instead of a
piece-wise holomorphic curve use a single holomorphic curve connecting a
selected pair of points. For instance, in order to connect in
${\mathbb{C}}^{3}$ the origin $(0,0,0)$ with an arbitrary point
$(x_{1},y_{1},z_{1})$ by a holomorphic curve admissible for the distribution
$\\{\ker(x{\rm d}y-y{\rm d}x-{\rm d}z)\\}$, one can consider the curve
$\big{(}x=x_{1}t,y=y_{1}t^{2}-ct(t-1),z=\frac{1}{3}x_{1}y_{1}t^{3}-\frac{1}{3}cx_{1}t^{3}\big{)}$
by choosing the constant $c$ satisfying the condition
$z_{1}=\frac{1}{3}x_{1}y_{1}-\frac{1}{3}cx_{1}$ (for $t=1$). This linear
equation is solvable if $x_{1}\neq 0$. If $x_{1}=0$ and $y_{1}\neq 0$, then we
can do this procedure by interchanging the variables. If $x_{1}=y_{1}=0$ the
task becomes extremely concrete: starting from the origin we have to reach the
point $(0,0,z_{1})$ by holomorphic curve admissible for the given
distribution.

Note that the quasi-homogeneous dilation $(x,y,z)\mapsto\big{(}\alpha x,\alpha
y,\alpha^{2}z\big{)}$ preserves our distribution, therefore, it suffices to
check local connectivity in the neighbourhood of the origin.

For the general holomorphic distribution of complex codimension 1 the Darboux
theorem holds: by the suitable holomorphic change of local coordinates the
form generating the distribution can be reduced to the normal form $x{\rm
d}y+{\rm d}z$ (see [3], for the proof of the classical Darboux theorem, and
[2], for the proof in the case of a holomorphic distribution). Thus, “generic”
holomorphic distributions are completely non-integrable in the sense that
locally they allow transitions between points of space along holomorphic
curves admissible for the distribution.

In the case of the space ${\mathbb{C}}^{3}$, in contrast to the general
multidimensional case, holomorphic curves integral (admissible) to the form
$x{\rm d}y+{\rm d}z$ are holomorphic Legendrian integral manifolds (of maximal
possible dimension). But all such varieties (in ${{\mathbb{C}}^{2n+1}}$ or in
${{\mathbb{R}}^{2n+1}}$) are described explicitly as follows

$\displaystyle y_{I}=\frac{\partial S}{\partial x_{I}},\qquad
x_{J}=-\frac{\partial S}{\partial y_{J}},\qquad z=S-x_{I}\frac{\partial
S}{\partial x_{I}},$

by means of the generating function $S$ (see [3, Appendix 4] on contact
structures). Here $S=S(x_{I},y_{J})$, and $I\cup J$ is any partition of the
indices $(1,\dots,n)$ when we consider the space ${\mathbb{R}}^{2n+1}$ or the
space ${\mathbb{C}}^{2n+1}$. In the case of ${\mathbb{C}}^{3}$ the function
$S$ depends only on one of the variables $x,y$. This allows one (by selecting
the function $S$ as a polynomial $ax^{2}+bx^{3}$ or as a polynomial
$ay^{2}+by^{3}$) directly verify that any point in the neighborhood of the
origin can be connected to the origin by a holomorphic curve integral for the
distribution $\\{\ker(x{\rm d}y+{\rm d}z)\\}$.

Now we will try to obtain an inverse theorem, that of the Carathéodory theorem
type: if a holomorphic distribution of complex codimension 1 does not allow
connections between any points in any neighborhood of the space by means of
holomorphic curves, admissible for the distribution, then the distribution is
integrable.

So, we will assume that in any neighborhood of every point $p$ there are
points inaccessible from $p$ by admissible holomorphic curves. We fix some
starting point $p$ together with its certain neighborhood. The idea of the
further heuristic reasoning could be as follows.111This idea can be seen
already in Poincaré [12, Section 193e]. Starting with the adiabatic
inaccessibility of some equilibrium states of the thermodynamic system,
Poincaré (and then Carathéodory) derives integrability of the corresponding
differential 1-form of heat influx and obtains one of the most important
characteristics of an equilibrium state of a thermodynamic system, the
entropy. Levels of the entropy function (adiabats, isoentropes) integrate the
distribution given by the heat influx 1-form.

The set of points of this neighborhood inaccessible (unreachable) from the
point $p$ is an open set (since if accessible points converge to a point, then
the limiting point is also accessible; of course, it requires further
mathematical justification, although from a physical point of view it is
almost a tautology).

Consider the boundary $\Gamma$ of the set of inaccessible points. Our starting
point $p$ (the center of the neighborhood where events occur) must lie on
$\Gamma$. The boundary $\Gamma$ is expected to be regular, due to the
smoothness of the distribution (for example, it cannot contain conic points
and corners).

The specific plane of the distribution attached at the point $p$ cannot be
transversal to $\Gamma$, as otherwise it would be possible to enter the region
of inaccessible points along this plane. Therefore, the considered plane of
the distribution turns out to be tangent to $\Gamma$.

We would like to show that in the case of the holomorphic distribution we
started with, $\Gamma$ has complex codimension 1, and therefore $\Gamma$ turns
out to be an integral surface of the initial holomorphic distribution.

Note that the holomorphic property of the distribution must be used
substantially. Indeed, consider, for example, the foliation of the space
${{\mathbb{C}}^{n}}$ by concentric spheres and look at the distribution of
complex tangents to the spheres. This distribution is not holomorphic (it is
parameterized by spheres). If it were holomorphic, then in the case of
${{\mathbb{C}}^{2}}$ it would even be holomorphically straightened and, like a
vector field, it would be integrable. Within each of the spheres, the
distribution of their complex tangents is completely non-integrable: there is
possible to connect any two points of the sphere along real one-dimensional
paths admissible for the distribution. However, it is impossible to leave the
sphere along such a path. Therefore in this setting in the neighborhood of
each point of the whole space there are points that cannot be reached even
along real one-dimensional paths admissible for the distribution.

Now we return to a holomorphic distribution in ${\mathbb{C}}^{3}$. A somewhat
more detailed study of the boundary set $\Gamma$ leads to the desired
conclusion that $\Gamma$ is a complex hypersurface and it is integral for the
initial holomorphic distribution of hyperplanes. Instead, we will replace such
a description of the set $\Gamma$ by the following shorter argument222It was
proposed by F. Forstnerič and is published with his kind permission. using
that distribution is holomorphic, while the space is three-dimensional. In the
space ${\mathbb{C}}^{3}$ the 3-form $\omega\wedge{\rm d}\omega$ is determined
by only one functional coefficient. If the coefficient is identically equal to
zero, then $\omega\wedge{\rm d}\omega\equiv 0$ and the Frobenius condition for
integrability of the distribution is fulfilled.

On the other hand, the non-connectivity assumption implies that there are no
points where this coefficient differs from zero. Indeed, if $\omega\wedge{\rm
d}\omega\not\equiv 0$, then the equation $\omega\wedge{\rm d}\omega=0$
determines a complex hypersurface in ${\mathbb{C}}^{3}$ outside of which the
holomorphic 2-form ${\rm d}\omega$ must be non-degenerate on $\ker\omega$.

But in a neighborhood of points where the form ${\rm d}\omega$ is
nondegenerate, the form $\omega$ can be written in the form $\omega=x{\rm
d}y+{\rm d}z$, thanks to the Darboux theorem. Hence, as we saw above, there
appears accessibility of all points of that neighborhood by admissible paths,
while this contradicts the assumption. Thus, the following theorem holds.
Consider a holomorphic distribution of complex hyperplanes in the space
${\mathbb{C}}^{3}$. If in any neighborhood of each point of the space there
are points inaccessible from the center of the neighborhood by holomorphic
curves admissible for the distribution, then the distribution is integrable.
Considering holomorphic distributions of complex hyperplanes we confine
ourselves to the first nontrivial case.

To complete the general case consider the form $\omega\wedge({\rm
d}\omega)^{k}$, where $k$ is minimal of the natural numbers such that
$\omega\wedge({\rm d}\omega)^{k}\equiv 0$. If $k=1$, we have the integrability
condition. On the other hand, $k$ must be equal to 1, otherwise the non-
connectivity condition will be violated. It follows from the Darboux theorem
applied to the form $\omega\wedge({\rm d}\omega)^{k-1}$ in the neighbourhood
of the point where the form is not equal to zero. (The following example can
explain the above reasoning. If in the space
${\mathbb{C}}^{5}={\mathbb{C}}^{3}\times{\mathbb{C}}^{2}$ instead of the form
$x_{1}{\rm d}y_{1}+x_{2}{\rm d}y_{2}+{\rm d}z$ one considers the form
$\omega=x_{1}{\rm d}y_{1}+{\rm d}z$, then $\omega\wedge({\rm
d}\omega)^{2}\equiv 0$ but $\omega\wedge{\rm d}\omega\neq 0$.)

Note that since the theorem is actually local, it suggests the similar
considerations for an arbitrary complex manifold.

##  5 Final comments

The idea of Poincaré mentioned above (very natural in the context of
thermodynamics) is particularly attractive, since it directly and explicitly
relates the problem of connectivity of the points of the space by admissible
paths and (non)integrability of the corresponding distribution, presenting the
integral surface $\Gamma$.

However, the detailed proof of the fact that the manifold $\Gamma$ is a
holomorphic hypersurface turns out to be rather intricate. Thus, it is worth
saying a few words about another approach to the relation between
integrability of the distribution and connectivity of points of the space by
integral paths of the distribution. We mean the following proof (essentially
being an illustration to the Frobenius theorem), which works both in real and
complex cases and for different dimensions of the space under consideration.

Above we have already mentioned Frobenius criterion $\omega\wedge{\rm
d}\omega\equiv 0$ for the distribution $\\{\ker\omega\\}$ to be an integrable
one. In turn, this condition is equivalent to the fact that the differential
${\rm d}\omega$ of the 1-form $\omega$ vanishes on any pair of vector fields
tangent to the planes of that distribution.

But according to Cartan’s formula, known in calculus of differential forms
(see, e.g., [7]), if $X$ and $Y$ are any smooth vector fields and $\omega$ is
a 1-form on a manifold, then ${\rm
d}\omega(X,Y)=X\omega(Y)-Y\omega(X)-\omega([X,Y])$, where $X\omega(Y)$ and
$Y\omega(X)$ are Lie derivatives of the functions $\omega(Y)$ and $\omega(X)$
along the fields $X$ and $Y$, respectively, and $[X,Y]$ is the Lie bracket
(commutator) of these vector fields.

If the vector fields belong to the planes of the distribution, then
$\omega(Y)\equiv 0$ and $\omega(X)\equiv 0$, and the equality ${\rm
d}\omega(X,Y)=0$ depends on whether the commutator $[X,Y]$ belongs to the
planes of the distribution.

If the commutator is not tangent to the distribution, then, as is well known,
by going along such fields one can move in directions transversal to the
planes of distribution (see [3, 9, 10]).

Thus, Cartan’s formula provides a direct connection between integrability of a
hyperplane distribution in space and connectivity of points of the space by
integral paths of the distribution. Namely, in the case of a typical
distribution the brackets of vector fields spanning the distribution of
hyperplanes can be transversal to the hyperplanes of distribution. This
provides local connectivity of points of the space by paths admissible for the
distribution.

For distributions of higher codimension, and even in degenerate cases of
hyperplane distributions, instead of integral surfaces of distribution, there
could appear such foliations of the space that the distribution is completely
non-integrable and the connectivity of points by admissible paths is possible
within each leaf of the foliation, but the transfer between the leaves is
impossible. We have already seen such a phenomenon in examples above. Singular
cases, of course, require special consideration and description. They were not
presented in the initial papers by Poincaré and Carathéodory mentioned above
and cited in the references.

In the introduction we briefly explained a relation of the issue discussed
above, as well as of Carathéodory’s theorem, with classical thermodynamics.
One can read more about that relation, e.g., in [14] and [15].

### Acknowledgements

I am very grateful to S.Yu. Nemirovski for fruitful discussions and for the
reference [2], where the authors not only report their impressive results, but
provide full proofs of several facts of mathematical folklore (e.g.,
holomorphic versions of the Darboux theorem), which they use, and which we
also needed above. I am also indebted to F. Forstnerič for his argument used
above, and for pointing out the new reference [1]. Special thanks to B. Khesin
who carefully read my source text and corrected not only my English but my
math as well. Additional thanks go to the referees for valuable remarks and
suggestions.

## References

  * [1] Alarcón A., Forstnerič F., Darboux charts around holomorphic Legendrian curves and applications, Int. Math. Res. Not. 2019 (2019), 893–922, arXiv:1702.00704.
  * [2] Alarcón A., Forstnerič F., López F.J., Holomorphic Legendrian curves, Compos. Math. 153 (2017), 1945–1986, arXiv:1607.00634.
  * [3] Arnold V.I., Mathematical methods of classical mechanics, 2nd ed., Graduate Texts in Mathematics, Vol. 60, Springer-Verlag, New York, 1989.
  * [4] Born M., Kritische Betrachtungen zur traditionellen Darstellung der Thermodynamik, Phys. Z. 22 (1921), 282–286.
  * [5] Carathéodory C., Untersuchungen über die Grundlagen der Thermodynamik, Math. Ann. 67 (1909), 355–386.
  * [6] Carathéodory C., Untersuchungen über die Grundlagen der Thermodynamik, in Gesammelte Mathematische Schriften, Band 2, C.H. Beck, München, 1955, 131–166.
  * [7] Cartan H., Calcul différentiel, Hermann, Paris, 1967.
  * [8] Chow W.-L., Über Systeme von linearen partiellen Differentialgleichungen erster Ordnung, Math. Ann. 117 (1939), 98–105.
  * [9] Gromov M., Carnot-Carathéodory spaces seen from within, in Sub-Riemannian Geometry, Progr. Math., Vol. 144, Birkhäuser, Basel, 1996, 79–323.
  * [10] Montgomery R., A tour of subriemannian geometries, their geodesics and applications, Mathematical Surveys and Monographs, Vol. 91, Amer. Math. Soc., Providence, RI, 2002.
  * [11] Nirenberg L., A complex Frobenius theorem, Sem. Analytic Functions 1 (1958), 172–189.
  * [12] Poincaré H., Thermodynamique, Gautier-Villar, Paris, 1908.
  * [13] Rashevskiy P.K., On the connectivity of any two points of a completely nonholonomic space with an admissible line, Uch. Zap. Moskov. Gos. Ped. Univ. im. K. Libknekhta Ser. Fiz.-Mat. 3 (1958), no. 2, 83–94.
  * [14] Zorich V., Mathematical analysis of problems in the natural sciences, Springer, Heidelberg, 2011.
  * [15] Zorich V., Mathematical aspects of classical thermodynamics, MCCME, Moscow, 2019\.


# Normalizing Flows for Bayesian Posteriors: Reproducibility and Deployment

Yukari Yamauchi<EMAIL_ADDRESS>Institute for Nuclear Theory, University of
Washington, Seattle, WA 98195, USA  Landon Buskirk  Facility for Rare Isotope
Beams, Michigan State University, East Lansing, Michigan 48824, USA  Pablo
Giuliani  Facility for Rare Isotope Beams, Michigan State University, East
Lansing, Michigan 48824, USA  Department of Statistics and Probability,
Michigan State University, East Lansing, Michigan 48824, USA  Kyle Godbey
Facility for Rare Isotope Beams, Michigan State University, East Lansing,
Michigan 48824, USA

###### Abstract

We present a computational framework for efficient learning, sampling, and
distribution of general Bayesian posterior distributions. The framework
leverages a machine learning approach for the construction of normalizing
flows for the general probability distributions typically encountered in
Bayesian uncertainty quantification studies. This normalizing flow can map a
trivial distribution to a more complicated one and can be stored more
efficiently than the empirical distribution samples themselves. Once the
normalized flow is trained, it further enables parallelized and uncorrelated
sampling of the learned distribution. We demonstrate our framework with three
test distributions with strong non-linear correlations, multi-modality, and
heavy tails, as well as with a realistic posterior distribution obtained from
a Bayesian calibration of a nuclear relativistic mean-field model. The
performance of the framework, as well as its relatively simple implementation,
positions it as one fundamental cornerstone in the development and deployment
of continuous calibration pipelines of physical models and as a key component
of future reproducible science workflows.

††preprint: INT-PUB-23-041

##  I Introduction

Uncertainty quantification (UQ) is a crucial aspect of scientific inquiry
across various domains, such as physics, biology, engineering, and social
sciences. A driving motivation of UQ studies is to assess the credibility of
experimental methods and measurements, numerical simulations, and physical
models, and to rigorously model and propagate the uncertainties coming from
all such sources. UQ efforts often directly inform decision making based on
the estimated uncertainties, such as optimizing for the best design of an
apparatus, selecting the most suitable parameter region of a theoretical
model, or evaluating the potential impact of an experimental campaign.

Given such an important and widespread application spectrum, a careful
consideration of reproducibility and accessibility should reside at the heart
of any UQ study. Reproducibility means that the same results can be obtained
by different researchers using the same data and methods, while accessibility
means that the data and methods are available and understandable to others. In
practice, achieving reproducibility and accessibility in UQ studies can be
challenging, especially when dealing with high-dimensional data sets and
Bayesian posterior distributions, commonly featured in many scientific
applications. One possible way to compress and share the relevant probability
distributions is to approximate them by multivariate Gaussian laws (see for
example covariance analysis for nuclear Density Functional Theory models [1,
2, 3]). This reduction, or simplification, comes at the cost of possible
misrepresentation of the underlying distributions, ignoring features such as
multi-modality, heavy tails, and non-linear relations between variables.
Variational Bayes methods are another approach for representing posterior
distributions (see [4, 5, 6, 7, 8] for reviews and examples of applications),
expanding the Gaussian law approximation by including a family of
distributions with parameters to be determined by the data. Normalizing flows
-realized through a machine learning approach- is another possibility, and the
method we explore in this paper.

Normalizing flows can model arbitrary probability distributions under some
reasonable assumptions [9] by transforming a simple base distribution, such as
a standard multivariate normal distribution, into a more complex target
distribution, using a sequence of invertible and differentiable functions. By
learning this mapping from the base distribution to the target distribution,
normalizing flows can capture the complex dependencies and variations in the
data, and provide a compact and expressive representation of the probability
distribution without relying on a specific functional form. Modern machine
learning methods provide a perfect computational workbench to effectively
learn such flows from the available data (see, e.g., [10, 11] for other
applications in Bayesian analysis).

In this work, we exploit the flexibility of normalizing flows to devise a
framework of efficient compression and distribution of high-dimensional
probability distributions. Such a framework allows us to store and generate
samples from arbitrary probability distributions, all at a minimal cost in
terms of computing resources and development effort. We organize the rest of
the manuscript as follows. In Section II we give an overview of Bayesian
statistics, explaining how posterior distributions are constructed from
physical models and data. In this section we also review the theory behind
normalizing flows and outline the machine learning framework we use for
training such maps. In Section III we demonstrate the effectiveness of the
approach on analytic distributions as well as published datasets from nuclear
theory. We close in Section IV with our conclusions and outlooks.

##  II Formalism

###  II.1 Bayesian Statistics

Bayesian statistics offers a probabilistic framework for making inferences
based on observed data within the context of uncertainty. One of its main
advantages is the explicit and clear statement of all assumptions concerning
the various sources of uncertainties and how they interact. This explicitness
allows for the transparent scrutiny and challenging of these assumptions,
facilitating improvements in the developed statistical framework. Another
major advantage is its capacity to continuously adapt knowledge based on new
data, all within the language of probabilities. In the Bayesian framework,
prior beliefs about model parameters and predictions are combined with
observed data to produce a posterior distribution, effectively refining our
understanding as new information becomes available. This process is
constructed through Bayes’ theorem [12, 13]:

$P(\boldsymbol{\alpha}|\boldsymbol{D})=\frac{P(\boldsymbol{D}|\boldsymbol{\alpha})P(\boldsymbol{\alpha})}{P(\boldsymbol{D})}.$
(1)

The _posterior_ distribution $P(\boldsymbol{\alpha}|\boldsymbol{D})$ encodes
the new information on the model parameters $\boldsymbol{\alpha}$ -or any
other quantities of interest- given the new observed data $\boldsymbol{D}$.
This posterior is constructed from the combination of the _likelihood_
$P(\boldsymbol{D}|\boldsymbol{\alpha})$, which models how the observed data
$\boldsymbol{D}$ could have arisen for a given parameter configuration
$\boldsymbol{\alpha}$, and the prior $P(\boldsymbol{\alpha})$, which accounts
for all the previous information and expert knowledge on the parameters
$\boldsymbol{\alpha}$. The _evidence_ $P(\boldsymbol{D})$ gives a sense on how
well the current model, including the statistical treatment for the
uncertainties, explains overall the data $\boldsymbol{D}$. In many
applications, this evidence term is often regarded as a normalization constant
that makes the posterior a true normalized probability distribution.

The specific form of the posterior distribution (1) varies wildly across the
many domains and applications of Bayesian statistics, adapting for each
problem’s specific characteristics. For relatively simple scenarios these
posteriors could take explicit forms -such as multivariate Gaussians- that are
easy to sample from and to calculate moments such as means and covariances.
For more complicated scenarios where non-linear models or unknown error scales
are involved, for example, obtaining information from the posterior can become
a challenge. Markov-Chain Monte Carlo (MCMC) sampling techniques, such as the
Metropolis Hastings algorithm and variations such as Hamiltonian Monte Carlo
methods [12], are usually employed in these situations. Many of these sampling
methods could require several chains (walkers) to go through thousands to
millions of samples to fully explore the parameter space, with the situation
turning more problematic the more parameters are involved (the curse of
dimensionality), or if the distributions show multi-modality and high-
correlation structures. The correlation function:

$\displaystyle C_{i}(\tau)$ $\displaystyle=$
$\displaystyle\Big{\langle}(\alpha_{i}(0)-\langle\alpha_{i}(0)\rangle)(\alpha_{i}(\tau)-\langle\alpha_{i}(\tau)\rangle)\Big{\rangle}$
(2) $\displaystyle\approx$
$\displaystyle\frac{1}{N}\sum_{t=1}^{N}(\alpha_{i}(t)-\langle\alpha_{i}\rangle)(\alpha_{i}(t+\tau)-\langle\alpha_{i}\rangle)$
$\displaystyle\propto$ $\displaystyle e^{-\tau/l}$

can be used to estimate the auto-correlation length $l$ of a parameter
$\alpha_{i}$ in the chain. However a precise measurement of $l$ typically
requires a large sampling size.

The usual requirement of such a vast number of model evaluations for
estimating any quantity of interest is one of the main computational
bottlenecks for proper uncertainty quantification studies, and one of the main
motivations for the development and use of several model emulation strategies
[14, 15, 16]. On the other hand, even when the sampling rate can be made fast
enough, relying on storing the samples for later distribution and for future
refining when new data becomes available quickly can become an inefficient
process, specially when dealing with large sets of calculated (model) data.
Normalizing flows, as we discuss in the following section, can serve to
improve this process.

###  II.2 Normalizing Flow

A normalizing flow (NF) is a map $\mathbb{R}^{N}\rightarrow\mathbb{R}^{N}$
which induces a non-trivial distribution from a trivial distribution. In our
context, NFs induce a desired $N-$dimensional density $p(\boldsymbol{\alpha})$
from the $N$-dimensional Gaussian distribution $g_{N}(\boldsymbol{\omega})$.
For the rest of this work we will use the shorthand notation
$p(\boldsymbol{\alpha})$ to represent any general distribution, not
necessarily normalized to integrate 1, a specific case being the un-normalized
posterior $\propto
P(\boldsymbol{D}|\boldsymbol{\alpha})P(\boldsymbol{\alpha})$ from Eq. (1).

In mathematical terms, the NF $f$ with
$\boldsymbol{\alpha}=f(\boldsymbol{\omega})$ for the density
$p(\boldsymbol{\alpha})$ satisfies the equation:

$d\boldsymbol{\alpha}\;p(\boldsymbol{\alpha})=d\boldsymbol{\omega}\;\det\left(\frac{\partial\boldsymbol{\alpha}}{\partial\boldsymbol{\omega}}\right)p(f(\boldsymbol{\omega}))=d\boldsymbol{\omega}\;g_{N}(\boldsymbol{\omega}),$
(3)

where
$\det\left(\frac{\partial\boldsymbol{\alpha}}{\partial\boldsymbol{\omega}}\right)$
is the Jacobian matrix. Under reasonable assumptions, such NF exist for any
normalizable (real and non-negative) distributions (see, e.g., [9] and
references therein). As illustration, in 1-dimension, a NF
$\omega\rightarrow\alpha$ for a distribution $\pi(\alpha)$ of a single
variable:

$\frac{d\alpha}{d\omega}\pi(\alpha)=g_{1}(\omega)$ (4)

can be constructed via the cumulative distribution functions:

$\Pi(\alpha)=\int_{-\infty}^{\alpha}d\alpha^{\prime}\;\pi(\alpha^{\prime})\;\;\mathrm{and}\;\;G_{1}(\omega)=\int_{-\infty}^{\omega}d\omega^{\prime}\;g_{1}(\omega^{\prime})\;\text{.}$
(5)

The cumulative distribution functions can be seen as a map from a probability
distribution to the uniform distribution on a unit interval. Therefore a NF
for $\pi$ is constructed as $f=\Pi^{-1}\circ G_{1}$, where “$\circ$” denotes
composition . In higher dimensions, NFs are known to exist albeit they cease
to be unique. Finding NFs with desirable properties is an active area of
research (see [17] for a review).

An exact NF, once constructed, enables us to obtain samples
$\boldsymbol{\alpha}$ from the original probability distribution
$p(\boldsymbol{\alpha})$ by simply sampling $\boldsymbol{\omega}$ from the
$N$-dimensional Gaussian distribution and applying the NF, i.e.
$\boldsymbol{\alpha}=f(\boldsymbol{\omega})$. In the same way, the expectation
value of any function $\mathcal{O}(\boldsymbol{\alpha})$ with respect to the
distribution $p(\boldsymbol{\alpha})$ can be evaluated as:

$\langle\mathcal{O}\rangle=\frac{\int
d\boldsymbol{\alpha}\;p(\boldsymbol{\alpha})\mathcal{O}(\boldsymbol{\alpha})}{\int
d\boldsymbol{\alpha}\;p(\boldsymbol{\alpha})}=\frac{\int
d\boldsymbol{\omega}\;g_{N}(\boldsymbol{\omega})\mathcal{O}(f(\boldsymbol{\omega}))}{\int
d\boldsymbol{\omega}\;g_{N}(\boldsymbol{\omega})},$ (6)

effectively computing $\mathcal{O}(\boldsymbol{\alpha})$ for each sample
$\boldsymbol{\alpha}=f(\boldsymbol{\omega})$ and taking the average.
Performing such computations from the NF offers two main advantages when
compared to MCMC approaches, especially for distributions in high dimensions
or those with large correlations. First, we can take samples from the trivial
Gaussian distribution and apply the NF in parallel, a sharp contrast with MCMC
methods where configurations are generated sequentially. This, in turn,
enables the parallelization of possibly time-consuming computations, for
example when calculating expectation values as in Eq. (6) with expensive
models. Second, there is no auto-correlations between samples in the Gaussian
sampling process, in contrast with MCMC methods, where auto-correlation
between neighboring samples is unavoidable for any distributions in any
dimensions. As will be discussed in Sec. III, a sequence of samples generated
by a MCMC algorithm can possess large auto-correlation lengths even for a
several-dimensional distribution.

Although an exact NF should exist for reasonably well behaved distributions
[9], finding such a map is often not a feasible task especially in high
dimensions. Therefore, in practice, we only aim to construct an approximate NF
by relaxing the condition in Eq. (3) and look for a function
$f:\boldsymbol{\omega}\rightarrow\boldsymbol{\alpha}$ that satisfies

$d\boldsymbol{\omega}\;g_{N}(\boldsymbol{\omega})=d\boldsymbol{\alpha}\;p^{\prime}(\boldsymbol{\alpha})\approx
d\boldsymbol{\alpha}\;p(\boldsymbol{\alpha})\;\text{,}$ (7)

where the level of the approximation is measured by a distance between two
distributions, such as the Jeffreys divergence (see Eq. (11) below), which is
what we use in this study.

When using such an approximate NF to sample from the original distribution
$p(\boldsymbol{\alpha})$, a correction via reweighting could be used to
improve the approximation. Simply, we generate Gaussian samples, apply the NF
and correct their weight by the ratio
$p(\boldsymbol{\alpha})/p^{\prime}(\boldsymbol{\alpha})$. In the same way, to
compute the expectation values using reweighting, we take the ratio:

$\langle\mathcal{O}\rangle_{p}=\frac{\Big{\langle}\frac{p(\boldsymbol{\alpha})}{p^{\prime}(\boldsymbol{\alpha})}\mathcal{O}(\boldsymbol{\alpha})\Big{\rangle}_{p^{\prime}}}{\Big{\langle}\frac{p(\boldsymbol{\alpha})}{p^{\prime}(\boldsymbol{\alpha})}\Big{\rangle}_{p^{\prime}}}$
(8)

where the subscript $p^{\prime}$ means that the expectation value is computed
according to $p^{\prime}(\alpha)$. For example,

$\Big{\langle}\frac{p(\boldsymbol{\alpha})}{p^{\prime}(\boldsymbol{\alpha})}\mathcal{O}(\boldsymbol{\alpha})\Big{\rangle}_{p^{\prime}}=\frac{\int
d\boldsymbol{\alpha}\;p(\boldsymbol{\alpha})\mathcal{O}(\boldsymbol{\alpha})}{\int
d\boldsymbol{\alpha}\;p^{\prime}(\boldsymbol{\alpha})}\;\text{.}$ (9)

The effective sampling size of the estimator Eq. (8) is given by
$1/\langle(p/p^{\prime})^{2}\rangle_{p^{\prime}}\in\left[0,1\right]$, which
takes the value 1 for a perfectly trained NF [18]. In practice, as we
demonstrate in Sec. III, a well-trained NF can closely reproduce the true
posterior density, in which case the reweighing is not necessary when the NF
is deployed. For the training of the NF we do use reweighing to accurately
evaluate the loss function, as we explain next.

###  II.3 Machine learning

The construction of an approximate normalizing flow (NF) for the posterior
density can be mapped onto a machine learning problem straightforwardly (see,
e.g., [19, 20, 21] for applications in lattice field theory). A neural network
(NN) can be trained to represent the NF so that when it is applied to samples
from the Gaussian distribution, the mapped samples approximately reconstruct
the target density. The use of machine learning is motivated for the following
reasons. First, since a well-designed NN is able to represent a large family
of functions efficiently, it has potential for representing a good
approximation to exact NFs for various posterior densities with a rather small
number of network parameters. Second, once an approximate NF is obtained, it
can be stored cheaply — it costs much smaller memory to store the trained NN
than keeping millions of samples of the model parameters. Furthermore, a NF
enables us to access more number of samples from the posterior density. In the
rest of the section, we detail our framework for building and deploying the NF
through a machine learning framework.

We begin by preparing the training data for which a NF is trained. We take a
supervised learning approach and prepare the data by sampling $N_{t}$ numbers
of $\boldsymbol{\alpha}$ via MCMC from the posterior density
$p(\boldsymbol{\alpha})$ and storing the pairs
$\\{\boldsymbol{\alpha}_{i},p(\boldsymbol{\alpha}_{i})\\}$
($i=1,\cdots,N_{t}$). The quantity and quality of the training data directly
affects the quality of the resulting NF. In particular, it is crucial that the
training data covers the range of $\boldsymbol{\alpha}$ where
$p(\boldsymbol{\alpha})$ has appreciable support.

The training stage is where various choices can be made for the neural network
architecture, the loss function, and the training procedure such as the
optimizer, the learning rate and its schedule, and the batch size. For the
neural network, we use a 6-layer Real NVP (non-volume preserving) [22]
followed by a layer which scales and shifts each output of the Real-NVP.
Schematically, the entire NN is:

$f(\boldsymbol{\omega})=L\circ\left(A_{o}\circ
A_{e}\right)^{6}(\boldsymbol{\omega})$ (10)

where $A_{e}$ and $A_{o}$ denote the affine coupling layers with even and odd
masking respectively. The final layer $L$, parameterized by two
$N$-dimensional vectors $\boldsymbol{a}$ and $\boldsymbol{b}$, acts on each
component of the output of Real NVP separately as
$L(\omega_{j};\boldsymbol{a},\boldsymbol{b})=a_{j}\omega_{j}+b_{j}$. A great
advantage of this NN is that the inverse map $f^{-1}$, and the determinant of
the Jacobian of $f$ and $f^{-1}$ can be computed very efficiently. Prior to
the training, the NN is initialized so that it represents an approximate
Gaussian fit to the train data. For more details of the NN’s architecture and
its initialization, see Appendix A.

For the loss function, we use the symmetrized Kullback–Leibler divergence,
called Jeffrey’s divergence [23]:

$D_{J}(p,p^{\prime})=\int
d\boldsymbol{\alpha}\;\left(\tilde{p}(\boldsymbol{\alpha})\log\frac{p(\boldsymbol{\alpha})}{p^{\prime}(\boldsymbol{\alpha})}+\tilde{p}^{\prime}(\boldsymbol{\alpha})\log\frac{p^{\prime}(\boldsymbol{\alpha})}{p(\boldsymbol{\alpha})}\right)\;\text{.}$
(11)

In our context, $p(\boldsymbol{\alpha})$ and $p^{\prime}(\boldsymbol{\alpha})$
are the exact density and the density induced by the NF, respectively, while
$\tilde{p}(\boldsymbol{\alpha})$ and $\tilde{p}^{\prime}(\boldsymbol{\alpha})$
are those normalized to unity. Noting that the train data is sampled according
to the true probability density $p(\boldsymbol{\alpha})$, to compute the first
term in the loss function, we simply sample uniformly a subset $N_{b}$ (batch)
of the training data ($N_{b}\leq N_{t}$) and compute the expectation value:

$\bigg{\langle}\log\frac{p(\boldsymbol{\alpha})}{p^{\prime}(\boldsymbol{\alpha})}\bigg{\rangle}_{p}\approx\sum_{k=1}^{N_{b}}\log\frac{p(\boldsymbol{\alpha}_{k})}{p^{\prime}(\boldsymbol{\alpha}_{k})}\;\text{.}$
(12)

The subscript $p$ means that the expectation value is computed according to
$p(\boldsymbol{\alpha})$, therefore the samples $\boldsymbol{\alpha}_{k}$ on
the right hand side are sampled from $p$. The value of the NF-induced
distribution for a given sample $\boldsymbol{\alpha}_{k}$ is computed as:

$p^{\prime}(\boldsymbol{\alpha}_{k})=\det\left(\frac{\partial\boldsymbol{\omega}}{\partial\boldsymbol{\alpha}}\right)g_{N}(\boldsymbol{\omega}_{k}),$
(13)

by using the inverse of the NF, $f^{-1}$, and obtaining
$\boldsymbol{\omega}_{k}=f^{-1}(\boldsymbol{\alpha}_{k})$. Given that we take
the supervised learning approach, the computation of the second term in Eq.
(11), which is the expectation value with respect to the NF-induced
distribution, requires reweighting:

$\bigg{\langle}\log\frac{p^{\prime}(\boldsymbol{\alpha})}{p(\boldsymbol{\alpha})}\bigg{\rangle}_{p^{\prime}}=\frac{\big{\langle}\frac{p^{\prime}(\boldsymbol{\alpha})}{p(\boldsymbol{\alpha})}\log\frac{p^{\prime}(\boldsymbol{\alpha})}{p(\boldsymbol{\alpha})}\big{\rangle}_{p}}{\big{\langle}\frac{p^{\prime}(\boldsymbol{\alpha})}{p(\boldsymbol{\alpha})}\big{\rangle}_{p}}\;\text{.}$
(14)

By combining the two terms, at each training step, after $N_{b}$ samples are
uniformly taken from the train data, the loss function is approximated by:

$D_{J}(p,p^{\prime})\approx\sum^{N_{b}}_{k}\log\frac{p(\boldsymbol{\alpha}_{k})}{p^{\prime}(\boldsymbol{\alpha}_{k})}+\frac{\sum_{k}^{N_{b}}\frac{p^{\prime}(\boldsymbol{\alpha}_{k})}{p(\boldsymbol{\alpha}_{k})}\log\frac{p^{\prime}(\boldsymbol{\alpha}_{k})}{p(\boldsymbol{\alpha}_{k})}}{\sum_{k}^{N_{b}}\frac{p^{\prime}(\boldsymbol{\alpha}_{k})}{p(\boldsymbol{\alpha}_{k})}}\text{.}$
(15)

Unlike the commonly-used Kullback–Leibler divergence, Jeffrey’s divergence
takes its minimal value zero only when the exact NF is represented by the
neural networks, i.e. $p=p^{\prime}$, regardless of the normalization of the
posteriors. This feature of Jeffrey’s divergence is suitable for the Bayesian
framework since we typically do not know the normalization of the posterior,
i.e, the evidence.

For the training hyperparameters, in this study, we use Adam [24] optimizer,
and fix the learning rate to $10^{-3}$ throughout the training process. At
each training step, $N_{b}=10^{3}$ samples are randomly taken from the train
data to estimate Jeffrey’s loss function, and the NN parameters are optimized
accordingly. We do not use any further techniques for training deep neural
networks. With our algorithm, we observed a consistent improvement in the loss
value as we increased the number of layers in Real NVP up to 6 for the study
of the relativistic mean field model in Sec. III.2.

After we have trained the NN we can directly obtain uncorrelated parameter
samples from the approximated posterior density. For computing expectation
values as in Eq. (6), we could furthermore correct for the slight deviation
from the true posterior distribution, as discussed in the previous section, by
computing the ratio
$p(\boldsymbol{\alpha}_{k})/p^{\prime}(\boldsymbol{\alpha}_{k})$ at each
sample $k$.

##  III Applications

Figure 1: Two-dimensional distributions $p(\alpha_{1},\alpha_{2})$ we consider
in this work for testing the NF framework: Rosenbrock-banana [25] (left),
Himmelblau [26] (center), and Student-t [12] (right). The colors in each
legend indicate the numerical value of each normalized density along the
respective contours.

To demonstrate the algorithm described in the previous section, we analyze
three two-dimensional distributions as tests: the Rosenbrock’s banana function
[25], an exponentiated Himmelblau distribution [26], and the Student’s
t-distribution [12]. Afterwards, we focus on a realistic case by learning the
NF for the distribution obtained in [15], which was constructed by performing
a Bayesian calibration of a relativistic mean field nuclear model. We used the
JAX library Equinox [27] for the implementation of the neural network in this
work. All the data, the codes, and the trained NF used to produce these
results are publicly available here.

###  III.1 2-dimensional test distributions

The three distributions we selected for testing the NF highlight three
challenging features: the Rosenbrock’s banana function has a highly correlated
non-linear structure, the exponentiated Himmelblau distribution is multi-
modal, and the Student’s t-distribution has heavy tails. These test
distributions are shown in Fig. 1.

For each distribution, we generated a single MCMC chain for the training data.
In the chain, 1000 steps were taken first for burn-out and the acceptance rate
was adjusted to be in the range $\left[0.3,0.5\right]$. Then one sample per
100 MCMC steps was taken to partially reduce the auto-correlations. A total of
$N_{t}=10^{5}$ samples were collected in each case and used in the training of
the respective NFs.

The Rosenbrock-banana distribution, up to a normalization constant, can be
written as:

$p_{b}(\alpha_{1},\alpha_{2})=e^{-(\alpha_{1}-c_{1})^{2}-c_{2}(\alpha_{1}^{2}-\alpha_{2})^{2}}\;\text{.}$
(16)

We choose to work with the coefficients $c_{1}=1.0$ and $c_{2}=20.0$. The
auto-correlation length in the MCMC chain was estimated to be
$l\sim[500-1000]$ from the correlation function in Eq. (2) — it is large due
to the high correlation structure between the two parameters. This structure
is apparent in Fig. 1: the distribution is not localized within a convex
boundary but rather extends and curves across the parameter space.

The exponentiated Himmelblau distribution, up to a normalization constant, can
be written as:

$p_{h}(\alpha_{1},\alpha_{2})=e^{-\frac{1}{c_{3}}[(\alpha_{1}^{2}+\alpha_{2}-c_{1})^{2}+(\alpha_{1}+\alpha_{2}^{2}-c_{2})^{2}]}.$
(17)

We choose to work with the coefficients $c_{1}=11.0$, $c_{2}=7.0$, and
$c_{3}=100.0$. The distribution has, as is shown in Fig. 1, four peaks that
are connected by a non-vanishing density. Therefore a MCMC chain could
reasonably explore those four regions and offer sufficient training data. The
auto-correlation length in the MCMC sampling was estimated to be $l\sim 100$
from the correlation function in Eq. (2).

The bivariate Student-t distribution, up to a normalization constant, can be
written as:

$p_{s}(\alpha_{1},\alpha_{2})=\left[1+\frac{y^{T}My}{c}\right]^{-(2+c)/2}.$
(18)

We choose to work with the mean and scaled variance matrix:

$y=\begin{pmatrix}\alpha_{1}-1.0\\\
\alpha_{2}-1.0\end{pmatrix}\;\;\mathrm{and}\;\;M=\begin{pmatrix}4.0&4.8\\\
4.8&9.0\end{pmatrix}^{-1}.$ (19)

The parameter $c$ determines how quickly the density decays to zero, with
smaller values creating heavier tails. We choose $c=3.$ in our study to have
very heavy tails, creating a distribution with only two well defined moments
(means and variances). As can be seen in Fig. 1, the density extends over a
big region of the 2D parameter space, with an appreciable linear correlation
between them. The auto-correlation length in the MCMC sampling was estimated
to be $l\sim 100$ from the correlation function in Eq. (2).

Figure 2: Comparison between the samples obtained by MCMC (black) and NF (red)
for the three test distributions $p(\alpha_{1},\alpha_{2})$: Rosenbrock-banana
(left), Himmelblau (center), and Student-t (right). For comparison, the
original distributions defined in Eqs. (16), (17), and (19) are also plotted
(blue) and labeled as “Int.” to stand for numerical integration. The NF
samples reproduce all the features of the corner plot of the MCMC, overlapping
well enough to make it hard to distinguish between them.

The $N_{t}$ MCMC training samples together with the same number of samples
obtained from the trained NFs are compared in Fig. 2. In all cases, the NFs
successfully reproduced key features of the test cases: the non-linear high
correlation, the multi-modal structure and the heavy tails.

To quantitatively assess the performance of the NF Table 1 presents the
calculations of the first two moments (means and covariances) of the three
distributions, as well as the final value of the loss function (15) during
training. In all cases, the NFs obtained very small values for Jeffrey’s
divergence, implying that the effective sampling size be close to one.
Therefore all results, means and covariances, from the NFs were computed with
respect to the approximate density $p^{\prime}$, i.e., without reweighting. In
all cases the obtained means and covariances from the NF are in good agreement
with the estimate from the train data (MCMC in the table). Both calculations
are also in agreement with those obtained by numerical integration (Int.),
with the Student-t variances being moderately under-estimated by the sampling
methods.

|  | Banana | Himmelblau | Student-t
---|---|---|---|---
$\langle\alpha_{1}\rangle$ | Int. | 1.00 | 0.110 | 1.00
MC | 1.004(8) | 0.124(8) | 1.00(1)
NF | 1.002(2) | 0.11(1) | 0.99(1)
$\langle\alpha_{2}\rangle$ | Int. | 1.50 | 0.230 | 1.00
MC | 1.51(2) | 0.238(7) | 1.00(3)
NF | 1.503(5) | 0.211(8) | 1.00(2)
$\text{var}(\alpha_{1})$ | Int. | 0.500 | 8.96 | 11.9
MC | 0.505(9) | 8.92(2) | 11.3(4)
NF | 0.500(2) | 8.97(2) | 11.2(2)
$\text{var}(\alpha_{2})$ | Int. | 2.53 | 6.39 | 26.8
MC | 2.5(1) | 6.42(2) | 26(1)
NF | 2.54(2) | 6.39(2) | 25.8(7)
$\text{cov}(\alpha_{1}\alpha_{2})$ | Int. | 1.00 | 0.225 | 14.4
MC | 1.01(3) | 0.20(3) | 13.7(6)
NF | 1.005(6) | 0.22(3) | 13.4(3)
Jeff. Div. |  | $1.0\times 10^{-5}$ | $3.8\times 10^{-3}$ | $6.0\times 10^{-5}$
Table 1: Comparison between the first two moments for each of the three test
distributions we considered. These moments are calculated using numerical
integration in python (Int.), as well as samples from Markov-Chain Monte Carlo
(MC), and Normalizing Flow (NF) in each case. The last digit in parentheses
for MC and NF shows the associated standard deviation on the calculation
obtained by data blocking and statistical bootstrap (see, e.g. [28] for the
methods). Jeffreys divergence between MCMC and NF -the loss function used to
train the NF in Eq. (15)- is also shown in the last row.

###  III.2 RMF

Having tested the NF framework on the three known 2-dimensional distributions,
we proceed to study a realistic Bayesian posterior without an explicit
functional form [15]. This posterior constrains the 8 parameters
$\boldsymbol{\alpha}$ from a relativistic mean field model [29] used to
calculate properties of finite nuclei and neutron stars. The original data
consisted of $5\times 10^{6}$ MCMC samples that were obtained from 8 walkers,
each with a 100,000 burn-out period. Obtaining such a large number of samples
in the original work was possible thanks to the use of a reduced order model,
since each full evaluation of the original model could take several minutes.
In constructing the train data for the NF, we took one sample per 100 MCMC
steps in each chain, resulting in a total of $N_{t}=5\times 10^{4}$ samples.

Figure 3: Comparison between the samples obtained by MCMC (black) and NF (red)
for the relativistic mean field Bayesian calibration [15]. Similar to the
situation in Fig. 2, the NF samples reproduce all the features of the corner
plot of the MCMC, overlapping well enough to make it hard to distinguish
between them. The distinction is most noticeable in the $J$ and $L$ univariate
plots.

The NF for this model was built following the steps described in Sec. II.3. We
initialized the neural network in Eq. (10) to an approximate Gaussian fit with
means and covariances estimated from the train data. The NN was then trained
to minimize the loss function in Eq. (15) by estimating it using
$N_{b}=10^{3}$ randomly-drawn samples from the train data at each training
step. As for the optimizer, we used again Adam with the learning rate of
$10^{-3}$.

Jeffrey’s divergence was $\sim 10$ at the initialization (with an approximate
Gaussian fit) and reduced to $2\times 10^{-3}$ after training. Fig. 3 shows
the comparison between the $5\times 10^{4}$ training MCMC samples, and the
same number of samples obtained from the trained NF. Qualitatively, the NF
reproduces all the features from the distribution estimated from the MCMC
samples, in particular the non-linear correlation between the symmetry energy
$J$ and its slope $L$. For a quantitative assessment, we show on top of the
diagonal univariate plots of Fig .3 for both the MCMC and the NF samples the
calculated means and the $95\%$ confidence interval around them. In all cases
the NF calculation matches the MCMC one almost perfectly.

##  IV Conclusions and Outlook

We have presented an efficient method to accurately reconstruct Bayesian
posterior distributions, through the use of NFs approximated through neural
networks, for the purpose of further sampling and deployment. The NF map is
learned by a conventional machine learning procedure by minimizing a loss
function -Jeffreys divergence in Eq. (11)- over training data in the form of a
finite number of MCMC samples from the posterior. The NF framework was able to
successfully learn the three test distributions with known forms we explored,
even in the presence of high non-linear correlations, multi-modality, and
heavy tails. In the more realistic scenario involving the Bayesian posterior
of the parameters of a physical model [15], the NF reproduced the target
distribution with a similar level of accuracy.

Once trained, the deployed NF enables us to obtain uncorrelated samples from
the posterior by directly sampling the Gaussian distribution. This in turns
can speed-up the evaluation of Bayesian integrals -such as expectation values-
by parallelizing the calculation over several independent nodes. Furthermore,
the trained NF requires only a small amount of memory to be stored,
facilitating both the dissemination of results from calibration efforts, as
well as increasing the reproducibility of studies based on those results. The
NF could also provide an alternative to MCMC based methods for the
construction of the posterior from the beginning, bypassing the need to
compute long correlated chains. Such self-learning approaches will also enable
us to overcome issues with supervised learning and improve the quality of NFs.

We believe that the implementation of NFs could become a standard procedure in
the workflow of distributing and updating results derived from Bayesian
calibrations. By using trained NFs from previous calibrations as the starting
point for updating posteriors distributions when new data becomes available,
this framework could enable a continuous calibration of complex models. In
this way, the full Bayesian philosophy is embraced by using all previous
information as prior knowledge, and continuously updating it as new vetted
data is made available. By deploying the trained NFs in turn, researchers can
directly draw samples from the distributions for their own purposes.
Furthermore, each trained NF can be tagged and distributed alongside metadata
defining the training process to provide a traceable and reproducible
workflow.

###### Acknowledgements.

Y.Y. is grateful to Scott Lawrence for many useful discussions and suggestions
for the machine learning tools and techniques used in this work. P.G. and K.G.
are grateful to Edgard Bonilla and Frederi Viens for useful discussions on the
selection of distributions and on the interpretation of the results. Y.Y. is
supported by the INT’s U.S. Department of Energy grant No. DE-FG02-00ER41132.
This work has also been supported by the National Science Foundation CSSI
program under award number 2004601 (BAND collaboration) and the U.S.
Department of Energy under award number DE-SC0023175 (NUCLEI SciDAC-5
collaboration). This work was supported in part through computational
resources and services provided by the Institute for Cyber-Enabled Research at
Michigan State University.

## References

  * Chen and Piekarewicz [2014] W.-C. Chen and J. Piekarewicz, Physical Review C 90, 044305 (2014).
  * Klüpfel  _et al._ [2009] P. Klüpfel, P.-G. Reinhard, T. Bürvenich, and J. Maruhn, Physical Review C 79, 034310 (2009).
  * Erler  _et al._ [2013] J. Erler, C. Horowitz, W. Nazarewicz, M. Rafalski, and P.-G. Reinhard, Physical Review C 87, 044320 (2013).
  * Jordan  _et al._ [1999] M. I. Jordan, Z. Ghahramani, T. S. Jaakkola, and L. K. Saul, Mach. Learn. 37, 183–233 (1999).
  * Titsias and Lázaro-Gredilla [2014] M. Titsias and M. Lázaro-Gredilla, in  _Proceedings of the 31st International Conference on Machine Learning_, Proceedings of Machine Learning Research, Vol. 32, edited by E. P. Xing and T. Jebara (PMLR, Bejing, China, 2014) pp. 1971–1979.
  * Blei  _et al._ [2017] D. M. Blei, A. Kucukelbir, and J. D. McAuliffe, Journal of the American statistical Association 112, 859 (2017).
  * Kejzlar  _et al._ [2023] V. Kejzlar, S. Bhattacharya, M. Son, and T. Maiti, The American Statistician 77, 85 (2023).
  * Kejzlar and Maiti [2023] V. Kejzlar and T. Maiti, Statistics and Computing 33, 18 (2023).
  * Kobyzev  _et al._ [2021] I. Kobyzev, S. J. Prince, and M. A. Brubaker, IEEE Transactions on Pattern Analysis and Machine Intelligence 43, 3964 (2021).
  * Gabrié  _et al._ [2021] M. Gabrié, G. M. Rotskoff, and E. Vanden-Eijnden, “Efficient bayesian sampling using normalizing flows to assist markov chain monte carlo methods,” (2021), arXiv:2107.08001 [stat.ML] .
  * Rizvi  _et al._ [2023] S. Rizvi, M. Pettee, and B. Nachman, arXiv preprint arXiv:2305.10500 (2023).
  * Gelman  _et al._ [1995] A. Gelman, J. B. Carlin, H. S. Stern, and D. B. Rubin,  _Bayesian Data Analysis_ (Chapman and Hall/CRC, 1995).
  * Phillips  _et al._ [2021] D. Phillips, R. Furnstahl, U. Heinz, T. Maiti, W. Nazarewicz, F. Nunes, M. Plumlee, M. Pratola, S. Pratt, F. Viens, _et al._ , Journal of Physics G: Nuclear and Particle Physics 48, 072001 (2021).
  * McDonnell  _et al._ [2015] J. D. McDonnell, N. Schunck, D. Higdon, J. Sarich, S. M. Wild, and W. Nazarewicz, Phys. Rev. Lett. 114, 122501 (2015).
  * Giuliani  _et al._ [2023] P. Giuliani, K. Godbey, E. Bonilla, F. Viens, and J. Piekarewicz, Frontiers in Physics 10, 1212 (2023).
  * Özge Sürer  _et al._ [2023] Özge Sürer, M. Plumlee, and S. M. Wild, Technometrics 0, 1 (2023), https://doi.org/10.1080/00401706.2023.2246157 .
  * Villani [2003] C. Villani,  _Topics in optimal transportation_ , 58 (American Mathematical Soc., 2003).
  * Vaitl  _et al._ [2022] L. Vaitl, K. A. Nicoli, S. Nakajima, and P. Kessel, Machine Learning: Science and Technology 3, 045006 (2022).
  * Albergo  _et al._ [2019] M. S. Albergo, G. Kanwar, and P. E. Shanahan, Phys. Rev. D 100, 034515 (2019), arXiv:1904.12072 [hep-lat] .
  * Nicoli  _et al._ [2021] K. A. Nicoli, C. J. Anders, L. Funcke, T. Hartung, K. Jansen, P. Kessel, S. Nakajima, and P. Stornati, Phys. Rev. Lett. 126, 032001 (2021), arXiv:2007.07115 [hep-lat] .
  * Lawrence and Yamauchi [2021] S. Lawrence and Y. Yamauchi, Phys. Rev. D 103, 114509 (2021), arXiv:2101.05755 [hep-lat] .
  * Dinh  _et al._ [2016] L. Dinh, J. Sohl-Dickstein, and S. Bengio, CoRR abs/1605.08803 (2016), 1605.08803 .
  * Jeffreys [1946] H. Jeffreys, Proceedings of the Royal Society of London. Series A. Mathematical and Physical Sciences 186, 453 (1946).
  * Kingma and Ba [2017] D. P. Kingma and J. Ba, “Adam: A method for stochastic optimization,” (2017), arXiv:1412.6980 [cs.LG] .
  * Rosenbrock [1960] H. Rosenbrock, The computer journal 3, 175 (1960).
  * Himmelblau [1972] D. M. Himmelblau,  _Applied Nonlinear Programming_ (McGraw-Hill, 1972).
  * Kidger and Garcia [2021] P. Kidger and C. Garcia, Differentiable Programming workshop at Neural Information Processing Systems 2021 (2021).
  * Gattringer and Lang [2010] C. Gattringer and C. B. Lang,  _Quantum chromodynamics on the lattice_, Vol. 788 (Springer, Berlin, 2010).
  * Todd and Piekarewicz [2003] B. G. Todd and J. Piekarewicz, Phys. Rev. C 67, 044317 (2003).

##  Appendix A Neural Network

The structure of the neural network used in this study is, as in Eq. (10),

$f(\boldsymbol{\omega})=L\circ\left(A_{o}\circ
A_{e}\right)^{6}(\boldsymbol{\omega})\;\text{.}$ (20)

The first block that the input $\boldsymbol{\omega}$ goes through is a 6-layer
Real-NVP (non-volume preserving) [22] neural network. Each layer consists of
two affine coupling layers — $A_{e}$ with the even masking $b_{e}$ followed by
$A_{o}$ with the odd masking $b_{o}$. The affine coupling layer is constructed
as

$A(\boldsymbol{\omega};b)=b\odot\boldsymbol{\omega}+(1-b)\odot\left[\boldsymbol{\omega}\odot\exp(s(b\odot\boldsymbol{\omega}))+t(b\odot\boldsymbol{\omega})\right]$
(21)

where $\odot$ is the element-wise product. The masking patterns are

$\displaystyle
b_{e}=\left[0,1,0,1,\cdots\right]\;\;\;\mathrm{and}\;\;\;b_{o}=\left[1,0,1,0,\cdots\right]$
(22)

We used a single linear layer for the network $t$ and an MLP with one hidden
layer and the hyperbolic tangent activation function for the network $s$. Note
that when the networks $s$ and $t$ output zero, the affine coupling layer acts
as the identity function. This can be arranged by setting all parameters in
$t$ and the final linear layer in $s$ to zero.

The second block of the neural network is a simple linear transformation that
acts on each output of the Real-NVP separately:

$L(\boldsymbol{\omega};\boldsymbol{a},\boldsymbol{b})=\boldsymbol{a}\odot\boldsymbol{\omega}+\boldsymbol{b}\;\text{.}$
(23)

The purpose of this layer is to encode a rough estimate of the mean and
variance of the parameters $\boldsymbol{\alpha}$, so the Real-NVP deals only
with numbers of order unity. This implementation not only makes the training
efficient but also was necessary for stabilizing the evaluation of Jeffreys’
divergence.

Prior to the training of the neural network, is it initialized so that it
yields an approximate Gaussian fit to the train data. This is arranged by
initializing the Real-NVP as the identity, and setting the coefficients in $L$
by the mean and variance approximated from the train data. The parameters in
the first layer of $s$ were set randomly from the Gaussian distribution.

In evaluating the loss function Eq. (11), it is important that the inverse of
the NF and the determinant of the Jacobian is computed efficiently for this
neural network via convolution. The inverse function of an affine coupling
layer is

$A^{-1}(\boldsymbol{\alpha};b)=b\odot\boldsymbol{\alpha}+(1-b)\odot(\boldsymbol{\alpha}-t(b\odot\boldsymbol{\alpha}))\odot\exp(-s(b\odot\boldsymbol{\alpha}))$
(24)

Therefore the inverse of the normalizing flow is

$f^{-1}(\boldsymbol{\alpha})=\left(A^{-1}_{e}\circ A^{-1}_{e}\right)^{6}\circ
L^{-1}(\boldsymbol{a},\boldsymbol{b})(\boldsymbol{\alpha})\;\text{.}$ (25)

The determinant of the Jacobian of the map $f$ Eq. (10) can be computed very
efficiently by combining that of the single affine coupling layer

$\displaystyle\log\det\left(\frac{\partial
A(\boldsymbol{\omega};b)}{\partial\boldsymbol{\omega}}\right)=(1-b)\cdot
s(b\odot\boldsymbol{\omega})$ (26)

and likewise for the inverse map (see [22] for further details).


# An intrinsic connection of space-time points

Ty Shedleski<EMAIL_ADDRESS>Alumni, University of Pittsburgh, USA
Muhammad Usman<EMAIL_ADDRESS><EMAIL_ADDRESS>Department of Space Science, Institute of Space Technology (IST), Pakistan
Space and Astrophysics Research Lab (SARL), National Center for GIS and Space
Applications (NCGSA), Pakistan

###### Abstract

Quantum field theory (QFT) describes the dynamics of quantum particles in the
quantum realm in the Minkowski space-time, whereas the General Relativity (GR)
is a classical theory describing the nature of dynamical behavior of large
bodies in different space-times. This research is a proposal to the proof of
concept that through the Einstein-Rosen bridge (also known as wormhole) the
information can travel between two points thus proving the classical
entanglement connection between two spatially distant points which are not
causally connected. These results introduce the classical entanglement between
the galactic black hole with its surrounding.

##  I Introduction

Quantum Field Theory (QFT) is the concept of explaining the nature and
dynamics of the particles at the very basic level depending upon their spin.
Spin $0\hbar$ particles are described by the Klein-Gordon equation, spin
$\nicefrac{{1}}{{2}}\hbar$ particles are described by the Dirac equation and
spin $1\hbar$ are described by the Proca equation. The hypothetical spin
$2\hbar$ particles i.e. graviton is yet to be completely understood. On the
other side (of relativity), the dynamics of any object or system of objects
can be described using the Einstein Field Equations (EFEs) but its
consequences and implications on the level of black holes, galaxies and
cosmology are yet to be completely understood as well as checked
observationally.

In this article, we demonstrate the energy flow from a harmonic oscillator the
wormhole. The results enable us to understand the concept of the wormhole
which is the prediction of the general relativity. Wormholes mathematically
come as the solution of the EFEs that are the (hypothetical) structures of
space-time that connect the causally disconnected space-time points. We find
that the energy flow will happen when a harmonic oscillator is placed in the
space-time of a wormhole and the phase-space is not restricted. The article’s
structure is as follows: in section (II), we describe the Klein-Gordon
equation in the curved space-time, section (III) describes the Dirac equation
in the curved space-time. In section (IV), the dynamics of wormhole, harmonic
oscillator when placed in the space-time of wormhole and energy flow is
presented. Section (V) concludes the article with summery of the work.

##  II The Klein-Gordon Equation in Curved Space-Time

We start this section by writing the Klein-Gordon Lagrangian density in flat
space-time

$\mathcal{L}=\frac{1}{2}\eta^{\mu\nu}\left(\partial_{\nu}\phi\right)\left(\partial_{\mu}\phi\right)-\frac{1}{2}m^{2}\phi^{2}$
(1)

$\phi$ is a scalar, therefore its covariant derivative will be equivalent to
the partial derivative $\partial_{\nu}\phi$. The space-time curvature will
lead to the addition of another term in eq. (1), proportional to the Ricci
curvature $R.$ Therefore, the Klein-Gordon Lagrangian density in curved space-
time can be written in a general form

$\displaystyle\mathcal{L}$
$\displaystyle=\mathcal{L}_{flat}+\mathcal{L}_{curved}$
$\displaystyle=\frac{1}{2}\eta^{\mu\nu}\left(\partial_{\nu}\phi\right)\left(\partial_{\mu}\phi\right)-\frac{1}{2}m^{2}\phi^{2}-\frac{1}{2}\xi
R\phi^{2}$ (2)

where $\xi$ is the gravitational coupling. Let us recall the action in curved
space-time and the d’Alembert rule

$\displaystyle S$ $\displaystyle=\int d^{4}x\sqrt{-g}\mathcal{L}$
$\displaystyle\delta S$ $\displaystyle=0$

the solution of the d’Alembert rule is the Euler-Lagrange equations

$\displaystyle\delta S$
$\displaystyle=\partial_{\mu}\left(\dfrac{\partial\left(\sqrt{-g}\mathcal{L}\right)}{\partial\left(\partial_{\mu}\phi\right)}\right)-\dfrac{\partial\left(\sqrt{-g}\mathcal{L}\right)}{\partial\phi}=0$

Thus

$\partial_{\mu}\left(\dfrac{\partial\left(\sqrt{-g}\mathcal{L}\right)}{\partial\left(\partial_{\mu}\phi\right)}\right)-\dfrac{\partial\left(\sqrt{-g}\mathcal{L}\right)}{\partial\phi}=0$
(3)

by inserting eq. (2) into Euler-Lagrange equations, we obtain the Klein-Gordon
equation in curved space-time

$\frac{1}{\sqrt{-g}}\partial_{\mu}\left(\sqrt{-g}g^{\mu\nu}\partial_{\nu}\phi\right)+\left(m^{2}+\xi
R\right)\phi=0$ (4)

As an example, for the FRW metric $g$, we write

$\sqrt{-g}=\frac{a^{3}r^{2}\sin\theta}{\sqrt{1-kr^{2}}}$ (5)

the eq. (4) can be rearranged

$\displaystyle\frac{1}{\sqrt{-g}}\partial_{\mu}\left(\sqrt{-g}g^{\mu\nu}\partial_{\nu}\phi\right)+\left(m^{2}+\xi
R\right)\phi$ $\displaystyle=0$
$\displaystyle\frac{1}{\sqrt{-g}}\partial_{0}\left(\sqrt{-g}g^{00}\partial_{0}\phi\right)+\frac{1}{\sqrt{-g}}\partial_{i}\left(\sqrt{-g}g^{ij}\partial_{j}\phi\right)+\left(m^{2}+\xi
R\right)\phi$ $\displaystyle=0$
$\displaystyle\partial_{0}^{2}\phi+\frac{1}{\sqrt{-g}}\partial_{i}\left(\sqrt{-g}g^{ij}\partial_{j}\phi\right)+\left(m^{2}+\xi
R\right)\phi$ $\displaystyle=0$ (6)

where the 0 indice is for the temporal component of the metric and $i,j=1,3$
are the spatial components. Note that there are only three equations instead
of four because $\theta\theta$ and $\phi\phi$ equations are the same. The
Laplacian for FRW space-time is

$\nabla=\frac{\sqrt{1-kr^{2}}}{r^{2}}\frac{\partial}{\partial
r}\left(\sqrt{1-kr^{2}}r^{2}\frac{\partial}{\partial
r}\right)+\frac{1}{r^{2}\sin^{2}\theta}\frac{\partial}{\partial\theta}\left(\sin\theta\frac{\partial}{\partial\theta}\right)+\frac{1}{r^{2}\sin^{2}\theta}\left(\frac{\partial^{2}}{\partial\varphi^{2}}\right)$
(7)

##  III Dirac Equation in Curved Space-Time

Let us also discuss the Dirac equation in curved space-time, i.e. the
covariant derivative and the coordinates are computed in a Riemannian
Manifold. Recall that the general Dirac equation is given by

$\left(i\gamma^{\mu}\nabla_{\mu}-m\right)\psi=0\leavevmode\nobreak\ ,$ (8)

the covariant derivative in the above equation $\nabla_{\mu}$ is equivalent to
$\partial_{\mu}$ in a Minkowskian manifold (flat space-time). The covariant
derivative is defined as

$\nabla_{\mu}\equiv\partial_{\mu}+\Gamma_{\mu}\leavevmode\nobreak\ ,$ (9)

where $\Gamma_{\mu}$ are the spin connections also known as the Christoffel’s
symbols (for torsion free space). They are given as

$\Gamma_{\mu}=\frac{1}{2}\sigma^{bc}e_{b}^{\nu}e_{c\nu;\mu}\leavevmode\nobreak\
.$ (10)

It is important to discuss the appropriate Dirac matrices, $\gamma^{\mu}$, in
the given curved space-time. The metric tensor of the curved space-time is
defined in a general coordinate system by

$g_{\mu\nu}\left(x\right)=e_{\mu}^{a}\left(x\right)e_{\nu}^{b}\left(x\right)\eta_{ab}$
(11)

where $e_{\mu}^{a}$, known as the dyads or the vielbein or simply the vector
directions, are

$e_{\mu}^{a}\left(x\right)=\left(\frac{\partial Y_{X}^{a}}{\partial
x^{\mu}}\right)_{x=X},\qquad a=0,\,1$ (12)

therefore, the $\gamma$ matrices in curved space-time are defined over the
dyads by

$\gamma^{\mu}=e_{\mu}^{a}\gamma^{a}$ (13)

in the eq. (13), the $\gamma^{a}$ are the flat Mikowskian $\gamma$ matrices
which satisfy

$\left\\{\gamma^{a},\gamma^{b}\right\\}=2\eta^{ab}\leavevmode\nobreak\
,\qquad\frac{1}{4}\left[\gamma^{a},\gamma^{b}\right]=\sigma^{ab}\leavevmode\nobreak\
.$ (14)

Using eq. (11) and eq. (13), a straightforward calculation gives us the
definition of the $\gamma$ matrices in curved space-time

$\left\\{\gamma^{\mu},\gamma^{\nu}\right\\}=2g^{\mu\nu}$ (15)

The importance of equation 14 is in the definition of the commutation rule
when used in the covariant derivative commutation relation with the gamma
matrices

$\left[\Gamma_{\mu},\gamma^{\mu}\right]=\frac{\partial\gamma^{\mu}}{\partial
x^{\mu}}+\Gamma^{\mu}_{\mu\rho}\gamma^{\rho},$ (16)

that leads us to introduce the Christoffel symbols

$\displaystyle\Gamma_{\nu\rho}^{\mu}$
$\displaystyle=\frac{1}{2}g^{\lambda\mu}\left[g_{\lambda\nu,\rho}+g_{\lambda\rho,\nu}-g_{\nu\rho,\lambda}\right]$
$\displaystyle\Gamma_{\mu\rho}^{\mu}$
$\displaystyle=\frac{1}{\sqrt{-g}}\partial_{\rho}\left(\sqrt{-g}\right)\leavevmode\nobreak\
.$ (17)

For simplicity, considering the simple (1+1) space-time case, eq. (16) gives

$\gamma^{\mu}\Gamma_{\mu}=\frac{1}{2}\gamma^{a}\frac{1}{\sqrt{-g}}\partial_{\mu}\left(\sqrt{-g}e_{a}^{\mu}\right)$
(18)

Then the Dirac equation given by eq. (8) now can be written

$\displaystyle\left(i\gamma^{\mu}\nabla_{\mu}-m\right)\psi=0$
$\displaystyle\left[i\gamma^{\mu}\left(\partial_{\mu}+\Gamma_{\mu}\right)-m\right]\psi=0$
$\displaystyle\left[i\gamma^{\mu}\partial_{\mu}+i\gamma^{\mu}\Gamma_{\mu}-m\right]\psi=0$
$\displaystyle\left[ie_{a}^{\mu}\gamma^{a}\partial_{\mu}+i\frac{1}{2}\gamma^{a}\frac{1}{\sqrt{-g}}\partial_{\mu}\left(\sqrt{-g}e_{a}^{\mu}\right)-m\right]\psi=0$

The Dirac equation, in natural units, in curved space-time is

$\left[i\gamma^{a}e_{a}^{\mu}\partial_{\mu}+\frac{i}{2}\gamma^{a}\frac{1}{\sqrt{-g}}\partial_{\mu}\left(\sqrt{-g}e_{a}^{\mu}\right)-m\right]\psi=0$
(19)

Let us now focus on the space-time part, since we are considering the work in
(1+1) space-time, a special feature of this space-time is that the metric can
always be reduced to the conformally flat form

$ds^{2}=\Omega^{2}\left(dt^{2}-dx^{2}\right)$ (20)

the function $\Omega\left(x,t\right)$ is pre-defined folowing the chosen
metric. We here are more interested at this moment in the general behavior of
this function to derive the Dirac equation. The Christoffel symbols are given
by

$\begin{cases}\Gamma_{00}^{0}=\Gamma_{11}^{0}=\Gamma_{10}^{1}=\Gamma_{01}^{1}=\frac{\dot{\Omega}}{\Omega}\\\
\Gamma_{01}^{0}=\Gamma_{10}^{0}=\Gamma_{00}^{1}=\Gamma_{11}^{1}=\frac{\Omega^{\prime}}{\Omega}\end{cases}$
(21)

The 0 and 1 stands for the time and space coordinates respectively, while
$\dot{\Omega}=\frac{\partial\Omega}{\partial
x^{0}}=\frac{\partial\Omega}{\partial t}$ and
$\Omega^{\prime}=\frac{\partial\Omega}{\partial
x^{1}}=\frac{\partial\Omega}{\partial x}$. The vielbein
$e_{a}^{\mu}\left(x\right)$ are given by

$e_{a}^{\mu}\left(x\right)=\left.\frac{\partial x^{\mu}}{\partial
y^{a}}\right|_{x^{\mu}=X^{\mu}}\rightarrow\left\\{e_{a}^{\mu}=\frac{1}{\Omega}\right\\}$
(22)

Another way to write the spin connections is

$\Gamma_{\mu}\left(x\right)=-\frac{i}{4}\omega_{ab\nu}\left(x\right)\sigma^{ab}$
(23)

where $\sigma^{ab}=\frac{i}{2}\left[\gamma^{a},\gamma^{b}\right]$. It is now
straightforward to show
$\left[\gamma^{a},\Gamma_{\mu}\right]=\omega_{b\nu}^{a}\gamma^{b}$. Thus we
can easily obtain the spin connection in terms of the affine connection

$\omega_{b\nu}^{a}=e_{\mu}^{a}\partial_{\nu}\left(e_{b}^{\mu}\right)+e_{\mu}^{a}e_{b}^{\sigma}\Gamma_{\sigma\nu}^{\mu}$
(24)

Now we compute the spin connection on our (1+1) spacetime

$\displaystyle\omega_{10}^{0}$
$\displaystyle=e_{\mu}^{0}\partial_{0}\left(e_{1}^{\mu}\right)+e_{\mu}^{0}e_{1}^{\sigma}\Gamma_{\sigma
0}^{\mu}$
$\displaystyle=e_{1}^{0}\partial_{0}\left(e_{1}^{1}\right)+e_{0}^{0}\partial_{0}\left(e_{1}^{0}\right)e_{1}^{0}e_{1}^{0}\Gamma_{00}^{1}+e_{0}^{0}e_{1}^{0}\Gamma_{00}^{0}+e_{1}^{0}e_{1}^{1}\Gamma_{10}^{1}+e_{0}^{0}e_{1}^{1}\Gamma_{10}^{0}$
$\displaystyle=\frac{\Omega^{\prime}}{\Omega}$ $\displaystyle\omega_{11}^{0}$
$\displaystyle=e_{\mu}^{0}\partial_{0}\left(e_{1}^{\mu}\right)+e_{\mu}^{0}e_{1}^{\sigma}\Gamma_{\sigma
1}^{\mu}$
$\displaystyle=e_{1}^{0}\partial_{0}\left(e_{1}^{1}\right)+e_{0}^{0}\partial_{0}\left(e_{1}^{0}\right)+e_{1}^{0}e_{1}^{0}\Gamma_{01}^{1}+e_{0}^{0}e_{1}^{0}\Gamma_{01}^{0}+e_{1}^{0}e_{1}^{1}\Gamma_{11}^{1}+e_{0}^{0}e_{1}^{1}\Gamma_{11}^{0}$
$\displaystyle=\frac{\dot{\Omega}}{\Omega}$

Thus using symmetry

$\begin{cases}\omega_{10}^{0}=\omega_{00}^{1}=\frac{\Omega^{\prime}}{\Omega}\\\
\omega_{11}^{0}=\omega_{01}^{1}=\frac{\dot{\Omega}}{\Omega}\end{cases}$ (25)

The equation (23) now leads to a general definition of $\Gamma_{\mu}$ for
$\mu=0,\,1$

$\displaystyle\Gamma_{0}$
$\displaystyle=-\frac{i}{4}\omega_{010}\frac{i}{2}\left[\gamma^{0},\gamma^{1}\right]$
$\displaystyle=\frac{\Omega^{\prime}}{8\Omega}\left[\gamma^{0},\gamma^{1}\right]$
(26) $\displaystyle\Gamma_{1}$
$\displaystyle=-\frac{i}{4}\omega_{101}\frac{i}{2}\left[\gamma^{0},\gamma^{1}\right]$
$\displaystyle=\frac{\dot{\Omega}}{8\Omega}\left[\gamma^{0},\gamma^{1}\right]$
(27)

Remarks: Recall the Dirac equation given by eq. (8)

$\displaystyle\left(i\gamma^{\mu}\nabla_{\mu}-m\right)\psi\left(x\right)=0$
$\displaystyle\left[i\gamma^{\mu}\partial_{\mu}+i\gamma^{\mu}\Gamma_{\mu}-m\right]\psi\left(x\right)=0$
$\displaystyle\left[i\gamma^{0}\partial_{0}+i\gamma^{1}\partial_{1}+i\gamma^{0}\Gamma_{0}+i\gamma^{1}\Gamma_{1}-m\right]\psi\left(x\right)=0$
$\displaystyle\left[i\gamma^{0}\partial_{0}+i\gamma^{1}\partial_{1}+i\gamma^{0}\frac{\Omega^{\prime}}{8\Omega}\left[\gamma^{0},\gamma^{1}\right]+i\gamma^{1}\frac{\dot{\Omega}}{8\Omega}\left[\gamma^{0},\gamma^{1}\right]-m\right]\psi\left(x\right)=0$
(28)

multiplying by $\Omega\gamma^{0}$ and using
$\gamma^{\mu}=e_{\mu}^{a}\gamma^{a}$

$\displaystyle\left[i\Omega\gamma^{0}\gamma^{0}\partial_{0}+i\Omega\gamma^{0}\gamma^{1}\partial_{1}+i\Omega\gamma^{0}\gamma^{0}\frac{\Omega^{\prime}}{8\Omega}\left[\gamma^{0},\gamma^{1}\right]+i\Omega\gamma^{0}\gamma^{1}\frac{\dot{\Omega}}{8\Omega}\left[\gamma^{0},\gamma^{1}\right]-\Omega\gamma^{0}m\right]\psi\left(x\right)=0$
$\displaystyle\left[i\partial_{0}+i\gamma^{0}\gamma^{1}\partial_{1}+i\frac{\Omega^{\prime}}{8\Omega}\left[\gamma^{0},\gamma^{1}\right]+i\gamma^{0}\gamma^{1}\frac{\dot{\Omega}}{8\Omega}\left[\gamma^{0},\gamma^{1}\right]-\Omega\gamma^{0}m\right]\psi\left(x\right)=0$

With some trivial rearrangement, we obtain

$i\left[\partial_{0}+\gamma^{0}\gamma^{1}\frac{\left[\gamma^{0},\gamma^{1}\right]}{8}\frac{\dot{\Omega}}{\Omega}\right]\psi\left(x\right)=i\left[\gamma^{0}\gamma^{1}\partial_{1}+\frac{\left[\gamma^{0},\gamma^{1}\right]}{8}\frac{\Omega^{\prime}}{\Omega}\right]\psi\left(x\right)+\Omega\gamma^{0}m\,\psi\left(x\right)$
(29)

If we choose $\gamma^{0}=\sigma_{z}$ and $\gamma^{1}=i\sigma_{y},$where
$\sigma_{i}$ are the Pauli matrices. The eq. (29) can be rewritten as

$i\left(\partial_{0}+\frac{\dot{\Omega}}{4\Omega}\right)\psi\left(x\right)=\left[-i\sigma_{x}\left(\partial_{1}+\frac{\Omega^{\prime}}{4\Omega}\right)+\Omega\sigma_{z}m\right]\psi\left(x\right)$
(30)

which completely is in agreement with the Dirac equation in flat space-time,
i.e $\Omega=1$,

$i\partial_{0}\psi\left(x\right)=-i\sigma_{x}\partial_{1}\psi\left(x\right)+\sigma_{z}m\psi\left(x\right)$
(31)

##  IV The Einstein-Rosen Bridge, Wormhole and Energy Flow

Time travel in general relativity is described by the theories of wormhole.
Consider the charge free Schwarzschild space-time

$\mathrm{d}s^{2}=-c^{2}\left(1-\dfrac{2Gm}{rc^{2}}\right)\mathrm{d}t^{2}+\dfrac{1}{1-\dfrac{2Gm}{rc^{2}}}\mathrm{d}r^{2}+r^{2}\left(\mathrm{d\theta^{2}+\sin^{2}\theta\mathrm{d}\phi^{2}}\right)\leavevmode\nobreak\
.$ (32)

Using the transformation

$u^{2}=c^{2}r-2Gm\leavevmode\nobreak\ ,$

the Schwarzschild space-time becomes

$\mathrm{d}s^{2}=-c^{2}\dfrac{u^{2}}{u^{2}+2Gm}\mathrm{d}t^{2}+4\dfrac{u^{2}+2Gm}{c^{4}}\mathrm{d}u^{2}+\left(\dfrac{u^{2}+2Gm}{c^{2}}\right)^{2}\left(\mathrm{d\theta^{2}+\sin^{2}\theta\mathrm{d}\phi^{2}}\right)\leavevmode\nobreak\
.$ (33)

Evidently, the above space-time can be described by two different congruent
parts or sheets, one corresponding to $u<0$ while the other corresponding to
$u>0$. These two sheets of planes are connected via a hyperplane $u=0$
(corresponding to $r=2Gm/c^{2}$) in which $g$ vanishes. This connection
between the two sheets is the Einstein-Rosen bridge or wormhole which is not
restricted by the gravity on either side.

The inclusion of electric charge in the Schwarzschild space-time leads to the
following Schwarzschild static spherically symmetric space-time

$\mathrm{d}s^{2}=-c^{2}\left(1-\dfrac{2Gm}{rc^{2}}-\dfrac{G}{4\pi\epsilon_{0}c^{4}}\dfrac{\varepsilon^{2}}{2r^{2}}\right)\mathrm{d}t^{2}+\dfrac{1}{1-\dfrac{2Gm}{rc^{2}}-\dfrac{G}{4\pi\epsilon_{0}c^{4}}\dfrac{\varepsilon^{2}}{2r^{2}}}\mathrm{d}r^{2}+r^{2}\left(\mathrm{d\theta^{2}+\sin^{2}\theta\mathrm{d}\phi^{2}}\right)\leavevmode\nobreak\
.$ (34)

where $\varepsilon$ is the electric charge. Using the transformation with
$m=0$

$u^{2}=\dfrac{1}{2}\left(8\pi\epsilon_{0}c^{4}r^{2}-G\varepsilon^{2}\right)\leavevmode\nobreak\
,$

the eq. (34) becomes

$\mathrm{d}s^{2}=-c^{2}\dfrac{2u^{2}}{\left(2u^{2}+G\varepsilon^{2}\right)}\mathrm{d}t^{2}+\mathrm{d}u^{2}+\left(\dfrac{2u^{2}+G\varepsilon^{2}}{8\pi\epsilon_{0}c^{4}}\right)\left(\mathrm{d\theta^{2}+\sin^{2}\theta\mathrm{d}\phi^{2}}\right)\leavevmode\nobreak\
.$ (35)

Once again, quite clear that the above space-time can be described by two
different congruent parts or sheets, one corresponding to $u<0$ while the
other corresponding to $u>0$. These two sheets of planes are connected via a
hyperplane $u=0$ (corresponding to
$r^{2}=G\varepsilon^{2}/8\pi\epsilon_{0}c^{4}$) in which $g$ vanishes. This
connection between the two sheets is the Einstein-Rosen bridge or wormhole
which is not restricted by the gravity on either side.

The neutral space-time interval corresponding to the non-traversable wormhole
is given by eq. (33). The phase-space of the wormhole allows for the energy
flow through the wormhole neck. Imagine a situation in which an harmonic
oscillator is oscillating. In this scenario, we want to see how the
oscillation in the oscillator will be effected by the presence of the
wormhole. In this case, the Einstein File Equations (EFEs) are needed to be
solved in the corresponding space-time with matter component taken as the
harmonic oscillator.

The Einstein tensor corresponding to the space-time given by eq. (33) is

$\displaystyle G_{tt}$ $\displaystyle=$
$\displaystyle\dfrac{u^{2}\left(c^{4}+16Gmu^{\prime 2}+8u^{2}u^{\prime
2}\right)}{\left(2Gm+u^{2}\right)^{3}}$ (36) $\displaystyle G_{rr}$
$\displaystyle=$
$\displaystyle\dfrac{4\left(-c^{6}r+8\left(-2Gmu-u^{3}\right)u^{\prime
2}-8\left(2Gm+u^{2}\right)u^{\prime\prime}\right)}{c^{6}\left(2Gmu+u^{3}\right)}$
(37) $\displaystyle G_{\theta\theta}$ $\displaystyle=$
$\displaystyle\dfrac{-4u\left(2Gm+u^{2}\right)u^{\prime
2}-6\left(2Gm+u^{2}\right)u^{\prime\prime}}{c^{6}u}$ (38) $\displaystyle
G_{\phi\phi}$ $\displaystyle=$
$\displaystyle\dfrac{\left(-4u\left(2Gm+u^{2}\right)u^{\prime
2}-6\left(2Gm+u^{2}\right)^{2}u^{\prime\prime}\right)\sin(\theta)^{2}}{c^{6}u}$
(39) $\displaystyle G_{\mu\nu}$ $\displaystyle=$ $\displaystyle
0;\qquad\qquad\qquad\qquad\qquad\qquad\text{when }\mu\neq\nu$ (40)

The stress-energy tensor for the harmonic oscillator is given by

$\displaystyle T_{tt}$ $\displaystyle=$
$\displaystyle-\left(\dfrac{1}{2}Mr^{\prime
2}+\dfrac{1}{2}M\omega^{2}r^{2}\right)$ (41) $\displaystyle T_{rr}$
$\displaystyle=$ $\displaystyle\dfrac{1}{2}Mr^{\prime
2}-\dfrac{1}{2}M\omega^{2}r^{2}$ (42) $\displaystyle T_{\theta\theta}$
$\displaystyle=$ $\displaystyle\dfrac{1}{2}Mr^{\prime
2}-\dfrac{1}{2}M\omega^{2}r^{2}$ (43) $\displaystyle T_{\phi\phi}$
$\displaystyle=$ $\displaystyle\dfrac{1}{2}Mr^{\prime
2}-\dfrac{1}{2}M\omega^{2}r^{2}$ (44) $\displaystyle T_{\mu\nu}$
$\displaystyle=$ $\displaystyle 0;\qquad\qquad\qquad\qquad\qquad\text{when
}\mu\neq\nu$ (45)

Equating both sides gives the Einstein Field Equations. Also the
${\theta\theta}$ and ${\phi\phi}$ components are the same so only one will be
kept. The solution of the Einstein field equations gives the dynamics of
harmonic oscillator placed in the wormhole. These are

$\displaystyle\dfrac{u^{2}\left(c^{4}+16Gmu^{\prime 2}+8u^{2}u^{\prime
2}\right)}{\left(2Gm+u^{2}\right)^{3}}$ $\displaystyle=$
$\displaystyle-\dfrac{8\pi G}{c^{4}}\left(\dfrac{1}{2}Mr^{\prime
2}+\dfrac{1}{2}M\omega^{2}r^{2}\right)$ (46)
$\displaystyle\dfrac{4\left(-c^{6}r+8\left(-2Gmu-u^{3}\right)u^{\prime
2}-8\left(2Gm+u^{2}\right)u^{\prime\prime}\right)}{c^{6}\left(2Gmu+u^{3}\right)}$
$\displaystyle=$ $\displaystyle\dfrac{8\pi
G}{c^{4}}\left(\dfrac{1}{2}Mr^{\prime 2}-\dfrac{1}{2}M\omega^{2}r^{2}\right)$
(47) $\displaystyle\dfrac{-4u\left(2Gm+u^{2}\right)u^{\prime
2}-6\left(2Gm+u^{2}\right)u^{\prime\prime}}{c^{6}u}$ $\displaystyle=$
$\displaystyle\dfrac{8\pi G}{c^{4}}\left(\dfrac{1}{2}Mr^{\prime
2}-\dfrac{1}{2}M\omega^{2}r^{2}\right)$ (48)

where $m$ is the mass of the black hole describing the geometry given by eq.
(33), $M$ is the mass of the harmonic oscillator placed in the space-time
given by eq. (33), $u(t)$ is the transformed spatial parameter of the black
hole geometry, $r(t)$ is the radial parameter of the harmonic oscillator and
$\omega$ is the angular frequency of the harmonic oscillator. From Eq. (46),
$M$ can be evaluated as

$M=-\dfrac{c^{4}}{4\pi G}\left(\dfrac{1}{r^{\prime
2}+\omega^{2}r^{2}}\right)\left(\dfrac{u^{2}\left(c^{4}+16Gmu^{\prime
2}+8u^{2}u^{\prime 2}\right)}{\left(2Gm+u^{2}\right)^{3}}\right)$ (49)

Using the above equation in the remaining two equations (eq. (47)and (48)) and
then solving the obtained equations numerically to get the dynamics of the
harmonic oscillator placed in the wormhole. The initial conditions used are
$u(0)=1$, $u^{\prime}(0)=0$ and $r(0)=1$. We have also used the natural units
so that the solution is simplified. The initial condition of harmonic
oscillator radial parameter $r(0)=1$ represents the perturbation at the
starting time $t=0$ while the initial conditions of the radial parameter of
the wormhole represent placement of oscillator is happening and due to the
perturbation the space-time of the wormhole is being disturbed equally. We
obtained two sets of solutions of the eqs. (46)-(48) which are

###  IV.1 Solution-I

In this outcome scenario, the dynamics of the oscillators happens such that
the velocity of the oscillator asymptotically increases such that the momentum
shoots to infinity as the oscillator approaches the event horizon. This is
shown in fig. (3).

(a) Evolution of radial parameter of wormhole. Figure 1: wormhole dynamics.

From the plot of the fig. (1), in which the radial parameter of the wormhole
is plotted as the function of time, we see that the wormhole does allow for
the oscillatory behavior to occur in an ideal scenario. There is once again no
damping effect in the ideal oscillatory behavior suggesting that the energy
can ideally flow without any loss (the case when there is some object placed
in the space-time is discussed in the next paragraph).

(a) Evolution of oscillator’s radial parameter.

(b) Evolution of oscillator’s radial parameter’s derivative.

Figure 2: Harmonic oscillator’s dynamics.

The plots shown in fig. (2) are the dynamics of harmonic oscillator when
placed in the wormhole space-time. Here we see that as we place an oscillator
in the region, the gravity of the massive object (Schwarzschild black hole
here) causes the oscillator to lose energy. The damping of the oscillations
are ultra high and the oscillations completely vanish in the first damping.
The initial energy $E_{0}$ of the oscillator is trasferred to the space-time
which would move through the wormhole.

Figure 3: Phase-space diagram of the oscillator.

The phase-space diagram of this scenario is presented in fig. (3). As
explained earlier, the dynamics of the oscillators happens such that the
velocity of the oscillator asymptotically increases such that the momentum
shoots to infinity as the oscillator approaches the event horizon. $0$ here
represents the location of the event horizon. This gain in the momentum would
come from the transfer of energy from the other side of the wormhole.

###  IV.2 Solution-II

In the second outcome scenario, the dynamics of the oscillators happens in
such a way that the velocity of the oscillator asymptotically increases such
that the momentum in this scenario does not shoot to infinity as the
oscillator approaches the event horizon. This is shown in fig. (6).

(a) Evolution of radial parameter of wormhole. Figure 4: wormhole dynamics.

From the plot presented in fig. (4), in which the radial parameter of the
wormhole is plotted as the function of time, we see that the wormhole here
allows for the oscillatory behavior to occur in the ideal scenario. This
result is as intuitively expected. There is no damping effect in the ideal
oscillatory behavior thus the energy flow is possible without any loss. The
situation when the oscillator is placed in the space-time is presented in fig.
(5) and fig. (6).

(a) Evolution of oscillator’s radial parameter.

(b) Evolution of oscillator’s radial parameter’s derivative.

Figure 5: Harmonic oscillator’s dynamics.

The plots shown in fig. (5) are the dynamics of the harmonic oscillator when
placed in the wormhole space-time. As we place an oscillator in the region,
the gravity of the massive object (Schwarzschild black hole here) causes the
oscillator to move towards the event horizon.

Figure 6: Phase-space diagram of the oscillator.

The phase-space plot of the oscillator is shown in fig. (6). Here we see that
the transfer of the momentum (and hence the energy transfer) is smooth. As
explained earlier, the dynamics of the oscillators happens in such a way that
the velocity of the oscillator asymptotically increases such that the momentum
in this scenario decreases as the event horizon is approached, therfore the
energy is being transferred from the oscillator to the other side of the
wormhole.

##  V Conclusion

We here present the wormhole formation in the very familiar Schwarzschild
space-time. The wormhole is typically described as the structure that
physically connects the two causally disconnected region of the space-time.
This connection is assisted/formed by the presence of some massive object.

In this article, we showed the proof of concept of the energy flow through the
wormhole by placing a harmonic oscillator in the wormhole space-time. The
harmonic oscillator is the natural choice of object since the objects are
expressible in terms of harmonic oscillators or their combination.

The Einstein Field Equations for the Schwarzschild geometry when a harmonic
oscillator is place in it gave us two sets of solutions. In one outcome
scenario, the dynamics of the oscillators occurred such that the velocity of
the oscillator increased asymptotically. The momentum of the oscillator
increased very rapidly as the consequence, we observe that as the oscillator
approaches event horizon the momentum of the oscillator approaches infinity.
This is shown in fig. (3). We noticed an ideal oscillatory behavior suggesting
that the energy can ideally flow without any loss which is presented in fig.
(1).

Another scenario which we obtained through EFEs solution showed oscillatory
behavior of the oscillator without damping ideally as in the earlier solution.
This once again suggested the energy flow without any loss in the ideal
scenario which is shown in fig. (4). The dynamics of the oscillators in this
situation occurred such that the velocity of the oscillator increased
asymptotically. The momentum of the oscillator however in this situation
decreases as the event horizon is approached, therefore the energy is being
transferred from the oscillator to the other side of the wormhole. This is
shown in fig. (6).

The results show that the harmonic oscillators on the opposite sides of the
wormhole transfer energy in between them through the wormhole. The oscillators
are causally disconnected yet the gain in the energy of the oscillator in the
first solution is occurring while the oscillator in the second solution is
losing energy. The energy conservation therefore would follow

$\Delta E_{osc-I}\leq\Delta E_{osc-II}\leavevmode\nobreak\ ,$

where $\Delta E_{osc-I}$ and $\Delta E_{osc-II}$ are the change in the energy
of the harmonic oscillator of solution-I and the change in the energy of the
harmonic oscillator of solution-II respectively.

## References

  * (1) Sinha A. Roychoudhury R. “Dirac Equation in (1+1)-Dimensional Curved Space-Time”. International Journal of Theoretical Physics, Vol. 33, No. 7, 1994.
  * (2) Syed Alwi B. Ahmad. “Fermion Quantum Field Theory in Black-Hole Spacetimes”. PhD thesis, Faculty of Virginia Polytechnic Institute and State University, 1997.
  * (3) Huang X. B. “Exact solutions of the Dirac equation in Robertson-Walker space-time”. http://arxiv.org/abs/gr-qc/0501077v1
  * (4) Yagdjian K., Galstian A. . “Fundamental solutions for the Dirac equation in curved spacetime and generalized Euler-Poisson-Darboux equation”. http://arxiv.org/abs/2103.01405v3
  * (5) J. Maldacena and L. Susskind. “Cool Horizons for Entangled Black Holes,” arXiv:1306.0533 [hep-th]
  * (6) L. Susskind. “ER=EPR, GHZ, and the Consistency of Quantum Measurements.”
  * (7) S. Ryu and T. Takayanagi, Phys. Rev. Lett. 96, 181602 (2006), [hep-th/0603001]
  * (8) K. Hori & al. “Mirror Symmetry”, 2003, American Mathematical Society
  * (9) Bom Soo Kim, “QFT Entanglement Entropy, 2D Fermion and Gauge Fields”, [hep-th/2201.02722v1]
  * (10) Landau, L.D.; Lifshitz, E.M. “The Classical Theory of Fields”, 4th ed.; Butterworth-Heinemann: Oxford, UK, 1975; Volume 2, ISBN 978-0-7506-2768-9.
  * (11) Tamburini F., Licata I., “General Relativistic Wormhole Connections from Planck-Scales and the ER = EPR Conjecture”, http://arxiv.org/abs/1912.12424v1.
  * (12) Kauffman L. H, “ER=EPR, Entanglement Topology and Tensor Networks”, http://arxiv.org/abs/2203.09797v2.
  * (13) Patrascu A. I, “Entanglement, space-time and the Mayer-Vietoris theorem”, http://arxiv.org/abs/1703.09610v2.
  * (14) Kauffman L. H, Mehrotra E., “Topological Aspects of Quantum Entanglement”, http://arxiv.org/abs/1611.08047v6 .
  * (15) Baez J. C., Vicary J., “Wormholes and Entanglement”, http://arxiv.org/abs/1401.3416v2.
  * (16) Numasawa T., Shiba N., Takayanagi T., Watanabe K., “EPR Pairs, Local Projections and Quantum Teleportation in Holography”, https://arxiv.org/abs/1604.01772.
  * (17) Sonner J., “Holographic Schwinger effect and the geometry of entanglement”, https://arxiv.org/abs/1307.6850