# Hamiltonian chaos and differential geometry of configuration space-time

Loris Di Cairano<EMAIL_ADDRESS>Institute of Neuroscience and
Medicine INM-9, and Institute for Advanced Simulation IAS-5, Forschungszentrum
Jülich, 52428 Jülich, Germany  Department of Physics, Faculty of Mathematics,
Computer Science and Natural Sciences, Aachen University, 52062 Aachen,
Germany  Matteo Gori<EMAIL_ADDRESS>Physics and Materials Science
Research Unit, University of Luxembourg, L-1511 Luxembourg, Luxembourg  Giulio
Pettini<EMAIL_ADDRESS>Dipartimento di Fisica Università di Firenze, and
I.N.F.N., Sezione di Firenze, via G. Sansone 1, I-50019 Sesto Fiorentino,
Italy  Marco Pettini<EMAIL_ADDRESS>Aix-Marseille University,
Marseille, France  CNRS Centre de Physique Théorique UMR7332, 13288 Marseille,
France

###### Abstract

This paper tackles Hamiltonian chaos by means of elementary tools of
Riemannian geometry. More precisely, a Hamiltonian flow is identified with a
geodesic flow on configuration space-time endowed with a suitable metric due
to Eisenhart. Until now, this framework has never been given attention to
describe chaotic dynamics. A gap that is filled in the present work. In a
Riemannian-geometric context, the stability/instability of the dynamics
depends on the curvature properties of the ambient manifold and is
investigated by means of the Jacobi–Levi-Civita (JLC) equation for geodesic
spread. It is confirmed that the dominant mechanism at the ground of chaotic
dynamics is parametric instability due to curvature variations along the
geodesics. A comparison is reported of the outcomes of the JLC equation
written also for the Jacobi metric on configuration space and for another
metric due to Eisenhart on an extended configuration space-time. This has been
applied to the Hénon-Heiles model, a two-degrees of freedom system. Then the
study has been extended to the 1D classical Heisenberg XY model at a large
number of degrees of freedom. Both the advantages and drawbacks of this
geometrization of Hamiltonian dynamics are discussed. Finally, a quick hint is
put forward concerning the possible extension of the differential-geometric
investigation of chaos in generic dynamical systems, including dissipative
ones, by resorting to Finsler manifolds.

Hamiltonian Chaos, Differential Geometry, Eisenhart metric

###### pacs:

05.20.Gg, 02.40.Vh, 05.20.- y, 05.70.- a

##  I Introduction

As is well known, a generic property of nonlinear dynamical systems, described
by a system of differential equations, is the presence of deterministic chaos.
This means that despite the deterministic nature of a dynamical system of this
kind, that is, despite the Cauchy’s theorem of existence and unicity of the
solutions of a system of differential equations, the property of
predictability of the dynamics for arbitrary times is lost in the absence of
stability of the dynamics chaos ; wiggins ; chaos1 . Such a dramatic
consequence of the breaking of integrability of a three body problem was
already pointed out by Poincaré while describing the complexity of the
homoclinic tangles in the proximity of hyperbolic points in phase space
poincare . It was at the beginning of the 60’s of the last century that for
the first time the consequences of homoclinic tangles in phase space of a
nonlinear Hamiltonian system became visually evident. This was thanks to the
numerical integration of the equations of motion of the celebrated Hénon-
Heiles model henon . The numerically worked out surfaces of section in phase
space displayed what Poincaré declared to be unable even to dare to attempt
drawing poincare . For many decades now, a huge amount of work has been done,
both numerical and mathematical, on deterministic chaos. However, especially
for many degrees of freedom systems, a theoretical explanation of the origin
of chaos has been lacking. Homoclinic intersections certainly provide an
elegant explanation of the origin of chaos in both dissipative and Hamiltonian
systems, but apply to 1.5 or two degrees of freedom systems. Beautiful
theorems on Axiom A systems chaos  and Anosov flows anosov  cannot account for
the emergence of chaos in dynamical systems of physical relevance. An
independent attempt to explain the origin of chaos in Hamiltonian systems was
put forward by N.S.Krylov who resorted to the possibility of identifying a
Hamiltonian flow with a geodesic flow in configuration space to try to explain
the origin of the dynamical instability (which we nowadays call deterministic
chaos) that could explain the spontaneous tendency to thermalization of many
body systems. Krylov’s pioneering approach focused on the search for negative
curvatures in configuration space equipped with a suitable metric krylov .
Krylov’s work inspired abstract ergodic theory but did not go too far to
explain the origin of chaos in Hamiltonian dynamical systems. For instance, in
the case of the already mentioned Hénon-Heiles model, it turns out that no
region of negative curvature can be found in configuration space, therefore
Krylov’s intuition has been discarded for a long time. However, more recently,
on the basis of numerical ”experiments” it has been shown that chaos in
Hamiltonian flows of physical relevance stems from another mechanism,
parametric instability, which will be discussed throughout this paper. The
Riemannian-geometric approach to explaining the origin of chaos in Hamiltonian
flows is based on two fundamental elements marco : i) the identification of a
Hamiltonian flow with a geodesic flow of a Riemannian manifold equipped with a
suitable metric, so that the geodesic equations

$\frac{d^{2}q^{i}}{ds^{2}}+\Gamma^{i}_{jk}\frac{dq^{j}}{ds}\frac{dq^{k}}{ds}=0~{}.$
(1)

coincide with Newton’s equations

$\frac{d^{2}q^{i}}{dt^{2}}=-\frac{\partial V(q)}{\partial q^{i}}~{}.$ (2)

a Hamiltonian flow - of which the kinetic energy is a quadratic form in the
velocities, that is, $\displaystyle
H=\frac{1}{2}a_{ik}p^{i}p^{k}+V(q_{1},\ldots,q_{N})$ \- is equivalent to the
solutions of Newton’s equations of motion stemming from a Lagrangian function
$\displaystyle
L=\frac{1}{2}a_{ik}\dot{q}^{i}\dot{q}^{k}-V(q_{1},\ldots,q_{N})$;

ii) the description of the stability/instability of the dynamics by means of
the Jacobi–Levi-Civita (JLC) equation for the geodesic spread measured by the
geodesic deviation vector field $\displaystyle J$ (which locally measures the
distance between nearby geodesics), which in a parallel-transported frame
reads

$\frac{d^{2}J^{k}}{ds^{2}}+R^{k}_{~{}ijr}\frac{dq^{i}}{ds}{J^{j}}\frac{dq^{r}}{ds}=0~{}.$
(3)

where $\displaystyle R^{k}_{~{}ijr}$ are the components of the Riemann-
Christoffel curvature tensor.

The most natural geometrization of Hamiltonian dynamics in a Riemannian
framework 111The natural and elegant geometric setting of Hamiltonian dynamics
is provided by symplectic geometry. This geometrical framework is very
powerful to study, for example, symmetries. However, symplectic manifolds are
not endowed with a metric, and without a metric we do not know how to measure
the distance between two nearby phase space trajectories and thus to study
their stability/instability through the time evolution of such a distance. is
a consequence of Maupertuis least action principle for isoenergetic paths

$\delta\ \int_{q(t_{0})}^{q(t_{1})}\ dt\ W(q,\dot{q})=0\ ,$ (4)

where $\displaystyle
W(q,\dot{q})=\\{[E-V(q)]a_{ik}{\dot{q}}^{i}{\dot{q}}^{i}\\}^{1/2}$, which is
equivalent to the variational definition of a geodesic line on a Riemannian
manifold, a line of stationary or minimum length joining the points
$\displaystyle A$ and $\displaystyle B$:

$\delta\ \int_{A}^{B}\ ds=0\ .$ (5)

If the subset of configuration space $\displaystyle
M_{E}=\\{(q_{1},\ldots,q_{N})\in{\mathbb{R}}^{N}|V(q_{1},\ldots,q_{N})<E\\}$
is given the non-Euclidean metric of components

$g_{ij}=2[E-V(q)]a_{ik}\ ,$ (6)

whence the infinitesimal arc element $\displaystyle
ds^{2}=4[E-V(q)]^{2}dq_{i}\ dq^{i}$, then Newton’s equations (2) are retrieved
from the geodesic equations (1).

The JLC equation for the geodesic spread can be rewritten as book

$\frac{d^{2}J^{k}}{ds^{2}}+2\Gamma^{k}_{ij}\frac{dq^{i}}{ds}\frac{dJ^{j}}{ds}+\left(\frac{\partial\Gamma^{k}_{ri}}{\partial
q^{j}}\right)\,\frac{dq^{r}}{ds}\frac{dq^{i}}{ds}\,{J^{j}}=0\ ,$ (7)

which has general validity independently of the metric of the ambient
manifold.

Importantly, there are other Riemannian manifolds, endowed with different
metric tensors, to geometrize Hamiltonian dynamics book . Two of these
alternatives are concisely described in the following. One brings about the
standard tangent dynamics equation as geodesic spread (JLC) equation, whereas
the second one has never been investigated hitherto to describe chaos in
Hamiltonian flows. This gap is filled in the present work. The choice among
these manifolds is driven by practical computational reasons as will be
discussed in what follows.

##  II Eisenhart Geometrization of Hamiltonian dynamics

It is worth summarizing some basic facts of a geometrization of Hamiltonian
dynamics which makes a direct and unexpected link between the standard tangent
dynamics equations, used to numerically compute Lyapunov exponents, and the
JLC equation for the geodesic spread book .

###  II.1 Eisenhart Metric on Enlarged Configuration Space-Time $\displaystyle
M\times\mathbb{R}^{2}$

L.P.Eisenhart proposed a geometric formulation of Newtonian dynamics that
makes use, as ambient space, of an enlarged configuration space-time
$\displaystyle M\times\mathbb{R}^{2}$ of local coordinates
$\displaystyle(q^{0},q^{1},\ldots,q^{i},\ldots,q^{N},q^{N+1})$. This space can
be endowed with a nondegenerate pseudo-Riemannian metric Eisenhart  whose arc
length is

$ds^{2}=\left(g_{e}\right)_{\mu\nu}\,dq^{\mu}dq^{\nu}=a_{ij}\,dq^{i}dq^{j}-2V(q)(dq^{0})^{2}+2\,dq^{0}dq^{N+1}~{},$
(8)

where $\displaystyle\mu$ and $\displaystyle\nu$ run from $\displaystyle 0$ to
$\displaystyle N+1$ and $\displaystyle i$ and $\displaystyle j$ run from 1 to
$\displaystyle N$. The relation between the geodesics of this manifold and the
natural motions of the dynamical system is contained in the following theorem
lichnerowicz :

Theorem. The natural motions of a Hamiltonian dynamical system are obtained as
the canonical projection of the geodesics of
$\displaystyle(M\times\mathbb{R}^{2},g_{e})$ on the configuration space-time,
$\displaystyle\pi:M\times\mathbb{R}^{2}\mapsto M\times\mathbb{R}$. Among the
totality of geodesics, only those whose arc lengths are positive definite and
are given by

$ds^{2}=c_{1}^{2}dt^{2}$ (9)

correspond to natural motions; the condition (9) can be equivalently cast in
the following integral form as a condition on the extra coordinate
$\displaystyle q^{N+1}$:

$q^{N+1}=\frac{c_{1}^{2}}{2}t+c^{2}_{2}-\int_{0}^{t}{L}\,d\tau~{},$ (10)

where $\displaystyle c_{1}$ and $\displaystyle c_{2}$ are given real
constants. Conversely, given a point $\displaystyle P\in M\times\mathbb{R}$
belonging to a trajectory of the system, and given two constants
$\displaystyle c_{1}$ and $\displaystyle c_{2}$, the point $\displaystyle
P^{\prime}=\pi^{-1}(P)\in M\times\mathbb{R}^{2}$, with $\displaystyle q^{N+1}$
given by (10), describes a geodesic curve in
$\displaystyle(M\times\mathbb{R}^{2},g_{e})$ such that $\displaystyle
ds^{2}=c_{1}^{2}dt^{2}$.

For the full proof, see lichnerowicz . Since the constant $\displaystyle
c_{1}$ is arbitrary, we will always set $\displaystyle c_{1}^{2}=1$ in order
that $\displaystyle ds^{2}=dt^{2}$ on the physical geodesics.

From (8) it follows that the explicit table of the components of the Eisenhart
metric is given by

$g_{e}=\left(\begin{array}[]{ccccc}-2V(q)&0&\cdots&0&1\\\
0&a_{11}&\cdots&a_{1N}&0\\\ \vdots&\vdots&\ddots&\vdots&\vdots\\\
0&a_{N1}&\cdots&a_{NN}&0\\\ 1&0&\cdots&0&0\\\ \end{array}\right)\ ,$ (11)

where $\displaystyle a_{ij}$ is the kinetic energy metric. The Christoffel
coefficients

$\Gamma^{i}_{jk}=\frac{1}{2}g^{im}\left(\frac{\partial g_{mk}}{\partial
q^{j}}+\frac{\partial g_{mj}}{\partial q^{k}}-\frac{\partial g_{jk}}{\partial
q^{m}}\right)$ (12)

for $\displaystyle g_{e}$ and with $\displaystyle a_{ij}=\delta_{ij}$ are
found to be non-vanishing only in the following cases

$\Gamma^{i}_{00}=-\Gamma^{N+1}_{0i}=\partial_{i}V~{},$ (13)

where $\displaystyle\partial_{i}=\partial/\partial q^{i}$ so that the geodesic
equations read

$\displaystyle\displaystyle\frac{d^{2}q^{0}}{ds^{2}}$
$\displaystyle\displaystyle=$ $\displaystyle\displaystyle 0~{},$ (14)
$\displaystyle\displaystyle\frac{d^{2}q^{i}}{ds^{2}}+\Gamma^{i}_{00}\frac{dq^{0}}{ds}\frac{dq^{0}}{ds}$
$\displaystyle\displaystyle=$ $\displaystyle\displaystyle 0\ ,$ (15)
$\displaystyle\displaystyle\frac{d^{2}q^{N+1}}{ds^{2}}+\Gamma^{N+1}_{0i}\frac{dq^{0}}{ds}\frac{dq^{i}}{ds}$
$\displaystyle\displaystyle=$ $\displaystyle\displaystyle 0\ ;$ (16)

using $\displaystyle ds=dt$ one obtains

$\displaystyle\displaystyle\frac{d^{2}q^{0}}{dt^{2}}$
$\displaystyle\displaystyle=$ $\displaystyle\displaystyle 0\ ,$ (17)
$\displaystyle\displaystyle\frac{d^{2}q^{i}}{dt^{2}}$
$\displaystyle\displaystyle=$ $\displaystyle\displaystyle-\frac{\partial
V}{\partial q_{i}}~{},$ (18)
$\displaystyle\displaystyle\frac{d^{2}q^{N+1}}{dt^{2}}$
$\displaystyle\displaystyle=$ $\displaystyle\displaystyle-\frac{d{L}}{dt}~{}.$
(19)

Equation (17) states only that $\displaystyle q^{0}=t$. The $\displaystyle N$
equations (18) are Newton’s equations, and (19) is the differential version of
(10).

The fact that in the framework of the Eisenhart metric the dynamics can be
geometrized with an affine parametrization of the arc length, i.e.,
$\displaystyle ds=dt$, will be extremely useful in the following, together
with the remarkably simple curvature properties of the Eisenhart metric.

####  II.1.1 Curvature of $\displaystyle(M\times\mathbb{R}^{2},g_{e})$

The curvature properties of the Eisenhart metric $\displaystyle g_{e}$ are
much simpler than those of the Jacobi metric, and this is obviously a great
advantage from a computational point of view. The components of the
Riemann–Christoffel curvature tensor are

$R^{k}_{~{}ijr}=\left(\Gamma^{t}_{ri}\Gamma^{k}_{jt}-\Gamma^{t}_{ji}\Gamma^{k}_{rt}+\partial_{j}\Gamma^{k}_{ri}-\partial_{r}\Gamma^{k}_{ji}\right)\
.$ (20)

Hence, and after Eq.(13), the only non-vanishing components of the curvature
tensor are

$R_{0i0j}=\partial_{i}\partial_{j}V$ (21)

hence the Ricci tensor has only one nonzero component

$R_{00}=\triangle V$ (22)

so that the Ricci curvature is

$K_{R}(q,\dot{q})=R_{00}\dot{q}^{0}\dot{q}^{0}\equiv\triangle V\ ,$ (23)

and the scalar curvature is identically vanishing
$\displaystyle{\mathscr{R}}(q)=0~{}.$

####  II.1.2 Geodesic Spread Equation for the Eisenhart Metric $\displaystyle
g_{e}$

The Jacobi equation (3) for $\displaystyle(M\times\mathbb{R}^{2},g_{e})$ takes
the form

$\displaystyle\displaystyle\frac{\nabla^{2}J^{0}}{ds^{2}}+R^{0}_{i0j}\frac{dq^{i}}{ds}J^{0}\frac{dq^{j}}{ds}+R^{0}_{0ij}\frac{dq^{0}}{ds}J^{i}\frac{dq^{j}}{ds}$
$\displaystyle\displaystyle=$ $\displaystyle\displaystyle 0\ ,~{}~{}~{}$ (24)
$\displaystyle\displaystyle\frac{\nabla^{2}J^{i}}{ds^{2}}+R^{i}_{0j0}\left(\frac{dq^{0}}{ds}\right)^{2}J^{j}+R^{i}_{00j}\frac{dq^{0}}{ds}J^{0}\frac{dq^{j}}{ds}+R^{i}_{j00}\frac{dq^{j}}{ds}J^{0}\frac{dq^{0}}{ds}$
$\displaystyle\displaystyle=$ $\displaystyle\displaystyle 0\ ,~{}~{}~{}$ (25)
$\displaystyle\displaystyle\frac{\nabla^{2}J^{N+1}}{ds^{2}}+R^{N+1}_{i0j}\frac{dq^{i}}{ds}J^{0}\frac{dq^{j}}{ds}+R^{N+1}_{ij0}\frac{dq^{i}}{ds}J^{j}\frac{dq^{0}}{ds}$
$\displaystyle\displaystyle=$ $\displaystyle\displaystyle 0\ ,~{}~{}~{}$ (26)

and since $\displaystyle\Gamma^{0}_{ij}=0$ and
$\displaystyle\Gamma^{i}_{0k}=0$ it is $\displaystyle\nabla
J^{0}/ds=dJ^{0}/ds$, $\displaystyle R^{0}_{~{}ijk}=0$, and
$\displaystyle{\nabla J^{i}}/{ds}={dJ^{i}}/{ds}$, the only accelerating
components of the vector field $\displaystyle J$ are found to obey the
equations

$\frac{d^{2}J^{i}}{ds^{2}}+\frac{\partial^{2}V}{\partial q_{i}\partial
q^{k}}\left(\frac{dq^{0}}{ds}\right)^{2}J^{k}=0\ .$ (27)

and using $\displaystyle dq^{0}/ds=1$ one is left with

$\frac{d^{2}J^{i}}{dt^{2}}+\frac{\partial^{2}V}{\partial q_{i}\partial q^{k}}\
J^{k}=0\ ,$ (28)

the usual tangent dynamics equations. This fact is a crucial point in the
development of a geometric theory of Hamiltonian chaos because there is no new
definition of chaos in the geometric context. In fact, the numerical Lyapunov
exponents computed by means of Eqs.(28) already belong to geometric treatment
of chaotic geodesic flows.

###  II.2 Eisenhart Metric on Configuration Space-Time $\displaystyle
M\times\mathbb{R}$

Another interesting choice of the ambient space and Riemannian metric to
reformulate Newtonian dynamics in a geometric language was also proposed by
Eisenhart Eisenhart . If and how the description of Hamiltonian chaos in this
framework is coherent with the results obtained by standard treatment based on
the tangent-dynamics/JLC equations discussed in the preceding section has
never been investigated before.

This geometric formulation makes use of an enlarged configuration space
$\displaystyle M\times\mathbb{R}$, with local coordinates
$\displaystyle(q^{0},q^{1},\ldots,q^{N})$, where a proper Riemannian metric
$\displaystyle G_{e}$ is defined to give

$ds^{2}=\left(G_{e}\right)_{\mu\nu}\,dq^{\mu}dq^{\nu}=a_{ij}\,dq^{i}dq^{j}+A(q)\,(dq^{0})^{2}~{},$
(29)

where $\displaystyle\mu$ and $\displaystyle\nu$ run from $\displaystyle 0$ to
$\displaystyle N$ and $\displaystyle i$ and $\displaystyle j$ run from 1 to
$\displaystyle N$, and the function $\displaystyle A(q)$ does not explicitly
depend on time. With the choice $\displaystyle 1/[2A(q)]=V(q)+\eta$ and under
the condition

$q^{0}=2\int_{0}^{t}V(q)\,d\tau+2\eta t\ ,$ (30)

for the extra variable it can easily be seen that the geodesics of the
manifold $\displaystyle(M\times\mathbb{R},G_{e})$ are the natural motions of
standard autonomous Hamiltonian systems. Since
$\displaystyle\frac{1}{2}a_{ij}\dot{q}^{i}\dot{q}^{j}+V(q)=E$, where
$\displaystyle E$ is the energy constant along a geodesic, we can see that the
following relation exists between $\displaystyle q^{0}$ and the action:

$q^{0}=-2\int_{0}^{t}T\,d\tau+2(E+\eta)t\ .$ (31)

Explicitly, the metric $\displaystyle G_{e}$ reads as

$G_{e}=\left(\begin{array}[]{cccc}[2V(q)+2\eta]^{-1}&0&\cdots&0\\\
0&a_{11}&\cdots&a_{1N}\\\ \vdots&\vdots&\ddots&\vdots\\\
0&a_{N1}&\cdots&a_{NN}\\\ \end{array}\right)\ ,$ (32)

and together with the condition (31), this gives an affine parametrization of
the arc length with the physical time, i.e., $\displaystyle
ds^{2}=2(E+\eta)dt^{2}$, along the geodesics that coincide with natural
motions. The constant $\displaystyle\eta$ can be set equal to an arbitrary
value greater than the largest value of $\displaystyle|E|$ so that the metric
$\displaystyle G_{e}$ is nonsingular. This metric is a priori very interesting
because it seems to have some better property than the Jacobi metric and than
the previous metric $\displaystyle g_{e}$. In fact, at variance with the
Jacobi metric $\displaystyle g_{J}$ in Eq.(6), the metric $\displaystyle
G_{e}$ is nonsingular on the boundary $\displaystyle V(q)=E$; moreover, by
varying the total energy $\displaystyle E$ we get a family of different
metrics $\displaystyle g_{J}$, whereas by choosing a convenient value of
$\displaystyle\eta$, at different values of the energy the metric
$\displaystyle G_{e}$ remains the same. The consequence is that a comparison
among the geometries of the submanifolds of
$\displaystyle(M\times\mathbb{R},G_{e})$—where the geodesic flows of different
energies “live”—is meaningful. To the contrary, this is not true with
$\displaystyle(M_{E},g_{J})$. In some cases, the possibility of making this
kind of comparison can be important. With respect to the Eisenhart metric
$\displaystyle g_{e}$ on $\displaystyle M\times\mathbb{R}^{2}$ in the previous
section, the metric $\displaystyle G_{e}$ on $\displaystyle M\times\mathbb{R}$
defines a somewhat richer geometry, for example the scalar curvature of
$\displaystyle g_{e}$ is identically vanishing, which is not the case of
$\displaystyle G_{e}$.

In the case of a diagonal kinetic-energy metric, i.e. $\displaystyle
a_{ij}\equiv\delta_{ij}$, the only non vanishing Christoffel symbols are

$\Gamma_{00}^{i}=\frac{(\partial V/\partial
q^{i})}{[2V(q)+2\eta]^{2}},~{}~{}~{}~{}~{}\Gamma_{i0}^{0}=-\frac{(\partial
V/\partial q^{i})}{[2V(q)+2\eta]}\ ,$ (33)

whence the geodesic equations

$\displaystyle\displaystyle\frac{d^{2}q^{0}}{ds^{2}}+\Gamma^{0}_{i0}\frac{dq^{i}}{ds}\frac{dq^{0}}{ds}+\Gamma^{0}_{0i}\frac{dq^{0}}{ds}\frac{dq^{i}}{ds}$
$\displaystyle\displaystyle=$ $\displaystyle\displaystyle 0\ ,$ (34)
$\displaystyle\displaystyle\frac{d^{2}q^{i}}{ds^{2}}+\Gamma^{i}_{00}\frac{dq^{0}}{ds}\frac{dq^{0}}{ds}$
$\displaystyle\displaystyle=$ $\displaystyle\displaystyle 0\ ,$ (35)

which, using the affine parametrization of the arc length with time, i.e.,
$\displaystyle ds^{2}=2(E+\eta)dt^{2}$, with
$\displaystyle(dq^{0}/dt)=2[V(q)+\eta]$ from (30), give

$\displaystyle\displaystyle\frac{d^{2}q^{0}}{dt^{2}}$
$\displaystyle\displaystyle=$ $\displaystyle\displaystyle 2\frac{d{V}}{dt}\ ,$
$\displaystyle\displaystyle\frac{d^{2}q^{i}}{dt^{2}}$
$\displaystyle\displaystyle=$ $\displaystyle\displaystyle-\frac{\partial
V}{\partial q_{i}},~{}~{}~{}~{}~{}~{}i=1,\dots,N~{},$ (36)

respectively. The first equation is the differential version of (30), and
equations (36) are Newton’s equations of motion.

####  II.2.1 Curvature of $\displaystyle(M\times\mathbb{R},G_{e})$

The basic curvature properties of the Eisenhart metric $\displaystyle G_{e}$
can be derived by means of the Riemann curvature tensor, which is found to
have the non-vanishing components

$R_{0i0j}=\frac{\partial_{i}\partial_{j}V}{(2V+2\eta)^{2}}-\frac{3(\partial_{i}V)(\partial_{j}V)}{(2V+2\eta)^{3}}\
,$ (37)

whence, after contraction, using $\displaystyle G^{00}=2V+2\eta$ the
components of the Ricci tensor are found to be

$\displaystyle\displaystyle R_{kj}$ $\displaystyle\displaystyle=$
$\displaystyle\displaystyle\frac{\partial_{k}\partial_{j}V}{(2V+2\eta)}-\frac{3(\partial_{k}V)(\partial_{j}V)}{(2V+2\eta)^{2}}\
,$ $\displaystyle\displaystyle R_{00}$ $\displaystyle\displaystyle=$
$\displaystyle\displaystyle\frac{\triangle V}{(2V+2\eta)^{2}}-\frac{3\|\nabla
V\|^{2}}{(2V+2\eta)^{3}}\ ,$ (38)

where $\displaystyle\triangle V=\sum_{i=1}^{N}{\partial^{2}V}/{\partial
q^{i\,2}}$, and thus we find that the Ricci curvature at the point
$\displaystyle q\in M\times\mathbb{R}$ and in the direction of the velocity
vector $\displaystyle\dot{q}$ is

$K_{R}(q,\dot{q})=\triangle V+R_{ij}\dot{q}^{i}\dot{q}^{j}$ (39)

and the scalar curvature at $\displaystyle q\in M\times\mathbb{R}$ is

${\mathscr{R}}(q)=\frac{\triangle V}{(2V+2\eta)}-\frac{3\|\nabla
V\|^{2}}{(2V+2\eta)^{2}}\ .$ (40)

####  II.2.2 Geodesic Spread Equation for the Eisenhart Metric $\displaystyle
G_{e}$

Let us now give the explicit form of Eq.(3) in the case of
$\displaystyle(M\times\mathbb{R},G_{e})$, the enlarged configuration space-
time equipped with one of the Eisenhart metrics. Since the nonvanishing
Christoffel coefficients are $\displaystyle\Gamma^{i}_{00}$ and
$\displaystyle\Gamma^{0}_{0i}$, then using the affine parametrization of the
arc length with physical time, we obtain

$\displaystyle\displaystyle\frac{d^{2}J^{k}}{dt^{2}}+\frac{2(\partial_{k}V)}{2V+2\eta}\frac{dJ^{0}}{dt}+\left[\partial_{kj}^{2}V-\frac{4(\partial_{k}V)(\partial_{j}V)}{2V+2\eta}\right]J^{j}$
$\displaystyle\displaystyle=$ $\displaystyle\displaystyle 0\ ,$
$\displaystyle\displaystyle\frac{d^{2}J^{0}}{dt^{2}}-\frac{2(\partial_{i}V)\dot{q}^{i}}{2V+2\eta}\frac{dJ^{0}}{dt}-2(\partial_{i}V)\frac{dJ^{i}}{dt}-\left[\partial_{ij}^{2}V-\frac{2(\partial_{i}V)(\partial_{j}V)}{2V+2\eta}\right]\dot{q}^{i}J^{j}$
$\displaystyle\displaystyle=$ $\displaystyle\displaystyle 0\ ,$

where the indexes $\displaystyle i,j,k$ run from $\displaystyle 1$ to
$\displaystyle N$. These equations have not yet been used to tackle
Hamiltonian chaos, but are certainly worth to be investigated.

As reported in Ref.cerruti1997lyapunov , the JLC equation in Eq.(7) is rather
complicated for the kinetic energy (Jacobi) metric in (6), it considerably
simplifies to (28) for $\displaystyle(M\times{\mathbb{R}}^{2},g_{e})$, and
displays an intermediate level of complexity for
$\displaystyle(M\times\mathbb{R},G_{e})$ as shown by Eqs.(II.2.2). This is
related with a different degree of ”richness” of the geometrical properties of
the respective manifolds. It is therefore important to check whether all these
geometrical frameworks provide the same information about regular and chaotic
motions rick ; cerruti1996geometric ; cerruti1997lyapunov , a necessary
condition which a-priori could be questioned as it was done in
Ref.cuervo2015non  even though the claims of this work have been proved wrong
in loris .

##  III Order and chaos in a paradigmatic two-degrees of freedom model with
$\displaystyle(M\times\mathbb{R},G_{e})$

The first benchmarking is performed for a two-degrees of freedom system. In
this case a paradigmatic candidate is the Hénon-Heiles model described by the
Hamiltonian

${H}=\frac{1}{2}\left(p_{x}^{2}+p_{y}^{2}\right)+\frac{1}{2}\left(q_{1}^{2}+q_{2}^{2}\right)+q_{1}^{2}q_{2}-\frac{1}{3}q_{2}^{3}\
.$ (42)

In this case, the JLC equation for the Jacobi metric is exactly written in the
form

$\displaystyle\displaystyle\frac{d^{2}J^{\perp}}{ds^{2}}$
$\displaystyle\displaystyle+$
$\displaystyle\displaystyle\frac{1}{2}\left[\frac{\triangle
V}{(E-V)^{2}}+\frac{\|\nabla V\|^{2}}{(E-V)^{3}}\right]\,J=0~{},$ (43)
$\displaystyle\displaystyle\frac{d^{2}J^{\parallel}}{ds^{2}}$
$\displaystyle\displaystyle=$ $\displaystyle\displaystyle 0$ (44)

where the expression in square brackets is the scalar curvature of the
manifold $\displaystyle(M_{E},g_{J})$, $\displaystyle g_{J}$ is the metric
tensor whose components are in Eq.(6), $\displaystyle J^{\perp}$ and
$\displaystyle J^{\parallel}$ are the components of the geodesic separation
vector transversal and parallel to the velocity vector along the reference
geodesic, respectively. It is well evident that this scalar curvature is
always positive and that chaotic motions can only be the consequence of
parametric instability due to the variability of the scalar curvature along
the geodesics. At first sight, the scalar curvature of
$\displaystyle(M\times{\mathbb{R}},G_{e})$ given in Eq.(40) can take also
negative values as is shown in Figure 1. On the one side this could add
another source of dynamical instability to parametric instability, but, on the
other side, the extension of regions of negative curvature depends on the
value of the arbitrary parameter $\displaystyle\eta$ that enters the metric
$\displaystyle G_{e}$, extension that can be arbitrarily reduced making its
contribution to degree of chaoticity not intrinsic. In Figure 2 the plane
$\displaystyle(q_{1},q_{2})$ is taken as surface of section of phase space
trajectories when $\displaystyle p_{2}=0$ and $\displaystyle p_{1}>0$.

Figure 1: Configuration space of the Hénon-Heiles model. The dashed lines
represent the equipotential boundaries: $\displaystyle V(q_{1},q_{2})=0.0833$
(cyan); $\displaystyle V(q_{1},q_{2})=0.125$ (green); $\displaystyle
V(q_{1},q_{2})=0.1667$ (yellow). Left panel: $\displaystyle\eta=0.045$. Right
panel: $\displaystyle\eta=0.1667$. The scale of colours represents different
intervals of values of the scalar curvature given in Eq.(40).

At the lowest energy, $\displaystyle E=0.0833$, when all the motions are
regular, the trajectories are found to visit also regions of negative
curvature, whereas at higher energies, $\displaystyle E=0.125$ and
$\displaystyle E=0.1667$, the chaotic trajectories considered display a large
number of intersections in regions of positive curvature. In other words, the
role of negatively curved regions does not appear to play a relevant role in
determining the chaotic instability of the dynamics.

Figure 2: Superposition of the configuration space of the Hénon-Heiles model
with the surfaces of section of phase space trajectories. Red dots correspond
to the crossing of the $\displaystyle(q_{1},q_{2})$ plane when $\displaystyle
p_{2}=0$ and $\displaystyle p_{1}>0$. Upper left panel corresponds to
$\displaystyle E=0.0833$; upper right panel corresponds to $\displaystyle
E=0.125$; lower panel corresponds to $\displaystyle E=0.1667$. For all these
cases $\displaystyle\eta=0.0833$.

As a matter of fact, the comparison of the results obtained by numerically
integrating the stability equations (28), (II.2.2), and (43) along with the
equations of motion of the Hénon-Heiles model, at different energies and
initial conditions, show an excellent qualitative and quantitative agreement.
The integration of the Hamilton equations of motion is performed with a
symplectic integrator. The stability equations have been integrated with a
fourth-order Runge-Kutta scheme. The choice of the energy values follows the
historical paper by Hénon-Heiles, and the initial conditions for regular and
chaotic motions are chosen according to the selections in
Ref.cerruti1996geometric . The quantity reported in Figures 3 and 4 is

$\lambda(t)=\frac{1}{t}\log\left[\frac{\|{\dot{J}}(t)\|^{2}+\|J(t)\|^{2}}{\|{\dot{J}}(0\|^{2}+\|J(0)\|^{2}}\right]$
(45)

where the separation vector $\displaystyle J$ is in turn the solution of the
three different stability equations.

Figure 3: Numerical solutions of the tangent dynamics equation (28) (black
line) compared to the solution of equation (II.2.2) (blue line), and to the
solution of equation (43) (red line). Left panel: $\displaystyle E=0.0833$,
$\displaystyle\eta=0.0833$ and the initial condition is point ($\displaystyle
a$) of Figure 1 of cerruti1996geometric . The dashed green line is the
reference $\displaystyle t^{-1}$ slope for regular motions. Right panel:
$\displaystyle E=0.125$, $\displaystyle\eta=0.0833$ and the initial condition
is point ($\displaystyle d$) of Figure 3 of cerruti1996geometric .

The robustness of the results obtained by means of Eq.(II.2.2) for the
manifold $\displaystyle(M\times\mathbb{R},G_{e})$ with respect to different
choices of the free parameter $\displaystyle\eta$ has been checked and
confirmed. It is in particular the close agreement between the results
obtained with the Eqs.(II.2.2) and (43) which confirms that chaos stems from
parametric instability, because in the latter equation the scalar curvature is
always positive. The right panel of Figure 3 shows a clear qualitative
agreement among the three patterns $\displaystyle\lambda(t)$ but some
quantitative deviations that do not change neither with longer integrations
not by changing the value of $\displaystyle\eta$ in the case of
$\displaystyle\lambda(t)$ computed with (II.2.2). Perhaps such a discrepancy
could stem from the inhomogeneity of the chaotic layer in phase space due to
the presence of very small regular islands, inhomogeneity detected differently
by the different JLC equations. Actually, this discrepancy is no longer
observed at higher energy (right panel of Figure 4) when the chaotic layer
seems more homogeneous. The reason why the geometrization of Hamiltonian
dynamics by means of $\displaystyle(M\times\mathbb{R},G_{e})$ can be of
prospective interest relies on its intermediate geometrical ”richness”.

Figure 4: Numerical solutions of the tangent dynamics equation (28) (black
line) compared to the solution of equation (II.2.2) (blue line), and to the
solution of equation (43) (red line). Here $\displaystyle E=0.1667$,
$\displaystyle\eta=0.0833$ and the initial condition for the left panel is
point ($\displaystyle a$) of Figure 5 of cerruti1996geometric , and for the
right panel point ($\displaystyle c_{2}$) of the same Figure.

On $\displaystyle(M\times{\mathbb{R}}^{2},g_{e})$ the scalar curvature is
always vanishing, the Riemann curvature tensor is just the Hessian of the
potential and the Ricci tensor has only one non-vanishing component, to the
opposite, on $\displaystyle(M_{E},g_{J})$ the Riemann curvature tensor has
$\displaystyle{\cal O}(N^{4})$ non-vanishing components and at large
$\displaystyle N$ the scalar curvature can happen to be overwhelmingly
negative without affecting the degree of chaoticity of the dynamics. The
geometry of $\displaystyle(M\times\mathbb{R},G_{e})$ is definitely richer than
that of $\displaystyle(M\times{\mathbb{R}}^{2},g_{e})$ and less complicated
than that of $\displaystyle(M_{E},g_{J})$, therefore, and mainly at large
$\displaystyle N$, this framework can offer some computational advantage for
more refined investigations about the geometric origin of parametric
instability of the geodesics. Loosely speaking, to give an idea of what a more
refined geometrical investigation might mean, it has been shown book ; cecmar
that integrability is related with the existence of Killing tensor fields on
the mechanical manifolds, therefore the degree of breaking of the hidden
symmetries associated with Killing tensor fields could be defined,
investigated, and related with the existence of weak and strong chaos in
Hamiltonian flows.

##  IV One-dimensional $\displaystyle XY$-model in the Eisenhart metric
$\displaystyle(M\times\mathbb{R},G_{e})$

Let us now proceed to investigate how Hamiltonian chaos is described in this
geometric framework at a large number of degrees of freedom. This is shown for
a specific model, the one-dimensional classical XY model. The reason for
choosing this model is that it has a rich variety of dynamical behaviors: at
low energy it is equivalent to a collection of weakly coupled harmonic
oscillators, at asymptotically high energy it represents a set of freely
rotating spins, at intermediate energies it displays a strongly chaotic
dynamics, as witnessed by the whole spectrum of Lyapounov exponents JSP .
Moreover, for this model it was necessary to introduce an ad hoc adjustment of
an otherwise successful geometric-statistical model for the analytic
computation of the largest Lyapounov exponent CasClePet  carried on in the
framework $\displaystyle(M\times{\mathbb{R}}^{2},g_{e})$. It is thus
interesting to check whether or not another geometric framework can allow to
fix the problem more naturally.

The 1D $\displaystyle XY$ model, describes a linear chain of $\displaystyle N$
spins/rotators constrained to rotate in a plane and coupled by a nearest-
neighbour interaction. This model is formally obtained by restricting the
classical Heisenberg model with $\displaystyle O(2)$ symmetry to one spatial
dimension. The potential energy of the $\displaystyle O(2)$ Heisenberg model
is $\displaystyle V=-{\cal I}\sum_{\langle i,j\rangle}{\bf s}_{i}\cdot{\bf
s}_{j}$, where the sum is extended only over nearest-neighbour pairs,
$\displaystyle{\cal I}$ is the coupling constant, and each $\displaystyle{\bf
s}_{i}$ has unit modulus and rotates in the plane. To each “spin”
$\displaystyle{\bf s}_{i}=(\cos q_{i},\sin q_{i})$, the velocity
$\displaystyle{\bf\dot{s}}_{i}=(-\dot{q}_{i}\sin q_{i},\dot{q}_{i}\cos q_{i})$
is associated, so that $\displaystyle{H}=\sum_{i=1}^{N}\frac{1}{2}\dot{\bf
s}_{i}^{2}-{\cal I}\sum_{\langle i,j\rangle}{\bf s}_{i}\cdot{\bf s}_{j}$. The
Hamiltonian of this model is then

$H(p,q)=\sum_{i=1}^{N}\frac{p_{i}^{2}}{2}+{\cal
I}\sum_{i=1}^{N}[1-\cos(q_{i}-q_{i-1})]~{},$ (46)

The canonical coordinates $\displaystyle q_{i}$ and $\displaystyle p_{i}$ are
thus given the meaning of angular coordinates and momenta. As already
mentioned above, this Hamiltonian system has two integrable limits. In the
low-energy limit it represents a chain of harmonic oscillators, as can be seen
by expanding the potential energy in power series

${H}(p,q)\approx\sum_{i=1}^{N}\left[\frac{p_{i}^{2}}{2}+\frac{{\cal
I}}{2}(q_{i+1}-q_{i})^{2}\right]~{},$ (47)

where $\displaystyle p_{i}=\dot{q}_{i}$, whereas in the high-energy limit it
represents a system of freely rotating objects, since the kinetic energy
increases with total energy without bounds, at variance with potential energy
which is bounded from above.

###  IV.1 Numerical solution of the JLC equation for
$\displaystyle(M\times\mathbb{R},G_{e})$

Let us proceed by comparing the outcomes of the integration of the equations
(28) and (II.2.2) computed along the flow of the Hamiltonian (46). The
standard tangent dynamics equations (28) can be split as

$\displaystyle\displaystyle\dot{J}^{i}_{q}$ $\displaystyle\displaystyle=$
$\displaystyle\displaystyle J^{i}_{p}$
$\displaystyle\displaystyle\dot{J}^{i}_{p}$ $\displaystyle\displaystyle=$
$\displaystyle\displaystyle-Hess(V)_{ij}\ J^{j}_{q}$ (48)

which explicitly read as

$\displaystyle\displaystyle\dot{J}^{i}_{q}$ $\displaystyle\displaystyle=$
$\displaystyle\displaystyle J^{i}_{p}$ (49)
$\displaystyle\displaystyle\dot{J}^{i}_{p}$ $\displaystyle\displaystyle=$
$\displaystyle\displaystyle-{\cal I}\cos(q_{i-1}-q_{i})J^{i-1}_{q}+{\cal
I}[\cos(q_{i-1}-q_{i})+\cos(q_{i}-q_{i+1})]J^{i}_{q}-{\cal
I}\cos(q_{i-1}-q_{i})J^{i+1}_{q}\ ,$

whence the Largest Lyapunov Exponent is worked out by computing

$\lambda_{1}=\lim_{t\rightarrow\infty}\frac{1}{t}\log\left[\frac{\|J_{q}(t)\|^{2}+\|J_{p}(t)\|^{2}}{\|J_{q}(0)\|^{2}+\|J_{p}(0)\|^{2}}\right]\
.$ (50)

At the same time, the integration of the JLC equations (II.2.2), by setting
$\displaystyle{J}=(J^{0},J^{i})$, and choosing $\displaystyle\eta=E$, yields
another estimate of the instability exponent through the analogous definition

$\lambda_{G}=\lim_{t\rightarrow\infty}\frac{1}{t}\log\left[\frac{\|{J}(t)\|^{2}_{G_{e}}+\|\dot{{J}}(t)\|^{2}_{G_{e}}}{\|{J}(0)\|^{2}_{G_{e}}+\|\dot{{J}}(0)\|^{2}_{G_{e}}}\right]\
.$ (51)

We have solved the equations of motion of the 1D XY model (setting
$\displaystyle{\cal I}=1$) and the tangent dynamics equations (49) by using a
bi-lateral symplectic algorithm lapo . The JLC equations (II.2.2) have been
solved by using a third-order predictor-corrector algorithm. Periodic boundary
conditions have been considered. Random initial conditions have been adopted
by taking the $\displaystyle q_{i}$ randomly distributed in the interval
$\displaystyle[0,2\pi]$, and by taking the $\displaystyle p_{i}$ gaussian-
distributed and suitably scaled so as to complement with the kinetic energy
the difference between the total energy initially set and the initial value of
the potential energy resulting from the random assignment of the
$\displaystyle q_{i}$. Figure 5 shows the comparison between the results
obtained at different values of the energy density $\displaystyle\epsilon=E/N$
for $\displaystyle\lambda_{1}(\epsilon)$ and
$\displaystyle\lambda_{G}(\epsilon)$ defined above. It is well evident that
the results so obtained are globally in very good agreement. At energy
densities in the interval between $\displaystyle\epsilon\simeq 0.2$ and
$\displaystyle\epsilon\simeq 100$ the agreement is perfect, whereas at lower
energy densities, below $\displaystyle\epsilon\simeq 0.2$, small discrepancies
are found which seem due to a slower time-relaxation of
$\displaystyle\lambda_{G}(t)$ with respect to $\displaystyle\lambda_{1}(t)$.

Of course, an unavoidable check of consistency has to be performed on an
integrable dynamics. This check has been performed on the flow of the
Hamiltonian (47). The results obtained with the equations (28) and (II.2.2)
are reported in Figure (6). As expected for non-chaotic dynamics, it is found
that $\displaystyle\lambda_{1}(t)$ decays as a straight line of slope
$\displaystyle-1$ in double logarithmic scale, and
$\displaystyle\lambda_{G}(t)$ decays with an oscillating pattern with a
$\displaystyle t^{-1}$ envelope. This has been checked at different
$\displaystyle N$ and energy values. Some cases are reported in Figure 6.

Figure 5: Lyapunov Exponents $\displaystyle\lambda_{1}$ (cyan circles) and
$\displaystyle\lambda_{G}$ (green triangles) versus the energy density
$\displaystyle\epsilon$ for a system of $\displaystyle N=150$ spins. The
parameter $\displaystyle\eta$ has been set as $\displaystyle\eta=E$. Figure 6:
Lyapunov Exponents $\displaystyle\lambda_{1}(t)$ (red, green and black lines)
versus $\displaystyle\lambda_{G}(t)$ (blue, magenta and cyan lines) for a
system of $\displaystyle N=2,100,1000$ harmonic oscillators, respectively. The
black dashed line is the $\displaystyle t^{-1}$ reference slope for a regular
dynamics. Here $\displaystyle\epsilon=1$ and $\displaystyle\eta=E$.

##  V The effective scalar model for the JLC equation

In CasClePet  an effective scalar approximation of the JLC equation (7) has
been worked out under some suitable hypothesis. In a nutshell, at large
$\displaystyle N$ under an hypothesis of quasi-isotropy - meaning that a
coarse-grained mechanical manifold appears as a constant curvature isotropic
manifold - with broad spatial spectrum of curvature variations at a finer
scale, the evolution of the norm of the geodesic separation vector is
described by a stochastic oscillator equation

$\frac{d^{2}\psi(s)}{ds^{2}}+\left[\langle
k_{R}\rangle+\langle\delta^{2}k_{R}\rangle^{1/2}\eta(s)\right]\psi(s)=0$

where $\displaystyle\eta(s)$ a $\displaystyle\delta$-correlated gaussian
stochastic process of zero mean and unit variance, and

$\begin{split}\langle k_{R}\rangle&=\frac{1}{N-1}\langle K_{R}\rangle\\\
\langle\delta^{2}k_{R}\rangle^{1/2}&=\frac{1}{N-1}(\langle
K_{R}^{2}\rangle-\langle K_{R}\rangle^{2})\end{split}$

where $\displaystyle K_{R}$ is the Ricci curvature of the mechanical manifold
under consideration, and the averages are meant along a reference geodesic or
as microcanonical averages on suitable energy surface
$\displaystyle\Sigma_{E}$. By putting $\displaystyle k_{0}=\langle
k_{R}\rangle$, $\displaystyle\sigma=\langle\delta^{2}k_{R}\rangle^{1/2}$,

$\begin{split}\tau_{1}&=\Big{\langle}\frac{dt}{ds}\Big{\rangle}\frac{\pi}{2\sqrt{k_{0}+\sigma}}\\\
\tau_{2}&=\Big{\langle}\frac{dt}{ds}\Big{\rangle}\frac{k_{0}^{1/2}}{\sigma}\end{split}$
(52)

and hence defining $\displaystyle\tau^{-1}=2(\tau_{1}^{-1}+\tau_{2}^{-1})$, an
analytic expression for a geometric Largest Lyapunov Exponent is given by
CasClePet

$\displaystyle\displaystyle\lambda(k_{0},\sigma,\tau)$
$\displaystyle\displaystyle=$
$\displaystyle\displaystyle\frac{1}{2}\left(\Lambda-\frac{4k_{0}}{3\Lambda}\right)\
,$ $\displaystyle\displaystyle\Lambda$ $\displaystyle\displaystyle=$
$\displaystyle\displaystyle\left(\sigma^{2}\tau+\sqrt{\left(\frac{4k_{0}}{3}\right)^{3}+\sigma^{4}\tau^{2}}\,\right)^{1/3}\
.$ (53)

This can be applied to the geometrization on the manifold
$\displaystyle(M\times{\mathbb{R}},G_{e})$ of Hamiltonian dynamics. In this
case the Ricci curvature reads as

$\begin{split}K_{R}(s)=\frac{1}{2(E+\eta)}\left(\Delta V-\frac{3\|\nabla
V\|^{2}}{2V+2\eta}+\frac{\partial^{2}_{kj}V\,\dot{q}^{j}\dot{q}^{k}}{2V+2\eta}-\frac{3\partial_{j}V\,\dot{q}^{j}\partial_{k}V\dot{q}^{k}}{(2V+2\eta)^{2}}\right)\equiv\frac{K_{R}(t)}{2(E+\eta)}\end{split}$
(54)

and using the arc-length parametrization $\displaystyle
ds^{2}=2(E+\eta)dt^{2}$ with physical time, we can compute by means of
Eqs.(53) an analytic prediction of $\displaystyle\lambda_{G}(\epsilon)$ for
$\displaystyle(M\times{\mathbb{R}},G_{e})$ and compare it to the outcome
obtained for $\displaystyle(M\times{\mathbb{R}}^{2},g_{e})$.

The first step consists in computing the average Ricci curvature and its
variance of the two manifolds at different values of the energy density. We
can limit these computations to one single choice of $\displaystyle N$ for
which the asymptotic values of $\displaystyle\langle k_{R}\rangle$ and
$\displaystyle\langle\delta^{2}k_{R}\rangle$ are already attained (see
CasClePet ). Moreover, for non-integrable systems, after the Poincaré-Fermi
theorem, all the constant energy surface is accessible to the dynamics, and
since chaos entails phase space mixing, with sufficiently long integration
times we obtain good estimate of microcanonical averages of the observables of
interest. Figures 7 and 8 provide the comparison between $\displaystyle\langle
k_{R}\rangle$ and $\displaystyle\langle\delta^{2}k_{R}\rangle$ for the two
manifolds.

Figure 7: Average of Ricci curvature $\displaystyle\langle K_{R}\rangle$ of
$\displaystyle{M\times{\mathbb{R}}^{2}}$ (red squares) and of $\displaystyle
M\times{\mathbb{R}}$ (green triangles), respectively, vs energy density
$\displaystyle\epsilon$ for a system of $\displaystyle N=150$. Here
$\displaystyle\eta=E$.

Somewhat unexpectedly these average quantities are found to be practically
coincident, thus it is not surprising that the application of the effective
scalar model for the JLC equation - recalled above - yields outcomes in close
agreement, as shown by Figure 9.

Figure 8: Average variance of the Ricci curvature $\displaystyle\sigma_{K}$ of
$\displaystyle{M\times{\mathbb{R}}^{2}}$ (red squares) and of $\displaystyle
M\times{\mathbb{R}}$ (green triangles) vs energy density
$\displaystyle\epsilon$ for a system of $\displaystyle N=150$ particles. Here
$\displaystyle\eta=E$. Figure 9: Geometric Lyapunov Exponents
$\displaystyle\lambda$ $\displaystyle\lambda$ worked out for
$\displaystyle{M\times{\mathbb{R}}^{2}}$ (red squares) and for $\displaystyle
M\times{\mathbb{R}}$ (green triangles) vs energy density
$\displaystyle\epsilon$, for a system of $\displaystyle N=150$ particles. Here
$\displaystyle\eta=E$. Figure 10: Comparison between the two Geometric
Lyapunov Exponents $\displaystyle\lambda_{g_{e}}$ (red squares),
$\displaystyle\lambda_{G_{e}}$ (green triangles) and the standard numerical
computation of $\displaystyle\lambda_{1}$ (cyan circles) vs energy density
$\displaystyle\epsilon$ for a system of $\displaystyle N=150$. Here
$\displaystyle\eta=E$.

The comparison among the outcomes $\displaystyle\lambda_{g_{e}}(\epsilon)$,
$\displaystyle\lambda_{G_{e}}(\epsilon)$ of the ”statistical” formula (53),
and the standard computation of $\displaystyle\lambda_{1}(\epsilon)$ are
displayed in Figure 10. The discrepancy, observed approximately for
$\displaystyle\epsilon$ in the interval between $\displaystyle 0.2$ and
$\displaystyle 2$, has been given an explanation in Ref.CasClePet  where it
has been shown that the numerical distribution of the Ricci curvature of
$\displaystyle{M\times{\mathbb{R}}^{2}}$ actually displays a non-vanishing
skewness with an excess of negative values with respect to a Gaussian
distribution. This information is lost in the effective scalar model for the
JLC equation above recalled. An ad hoc displacement of $\displaystyle\langle
k_{R}\rangle$ to empirically account for the excess of negative values of
$\displaystyle K_{R}$ allowed to exactly retrieve the pattern of
$\displaystyle\lambda_{1}(\epsilon)$ by means of the scalar effective model.
A-priori the use of $\displaystyle(M\times{\mathbb{R}},G_{e})$ could have
fixed the problem more naturally but, disappointedly, this has not been the
case thus calling for an improvement of the effective scalar model, possibly
taking into account higher order moments of the Ricci curvature distribution.
Finally, it is worth to mention that the potential function of the Hamiltonian
(46) has a large number of critical points $\displaystyle q_{c}$, that is such
that $\displaystyle\nabla V(q)|_{q=q_{c}}=0$ book ; near each critical point,
in Morse chart one has $\displaystyle
V(q)=V(q_{c})-\sum_{i=1}^{k}q_{i}^{2}+\sum_{i=k+1}^{N}q_{i}^{2}$ where
$\displaystyle k$ is the Morse index of a given critical point. Now, the
neighborhoods of critical points are enhancers of chaos because using the
expression for $\displaystyle V(q)$ in Morse chart together with
$\displaystyle\nabla V(q_{c})=0$, both equations (28) and (II.2.2) diagonalize
with $\displaystyle k$ unstable components in proximity of a critical point of
index $\displaystyle k$. Morse theory relates critical points of a suitable
real valued function (here the potential function) with topological properties
of its levels sets, here of equipotential manifolds in configuration space. In
other words, the 1D XY model highlights the necessity of taking into account
also some topological property of the mechanical manifolds in order to improve
the effective scalar model for the JLC equation.

##  VI Discussion

Summarizing, the geometrization of Hamiltonian dynamics within the framework
of the configuration space-time equipped with an Eisenhart metric,
$\displaystyle(M\times{\mathbb{R}},G_{e})$, provides a correct distinction of
regular and chaotic motions and it is in qualitative and quantitative
agreement with the two other geometrization frameworks reported above. As
already remarked, the advantage of this framework could be that of an
intermediate level of complexity/richness of its geometry with respect to
$\displaystyle(M_{E},g_{J})$ and
$\displaystyle(M\times{\mathbb{R}}^{2},g_{e})$ which could be useful in
performing more elaborated investigations about the relation between geometry
and chaos.

Let us conclude with an outlook at a prospective extension to generic
dynamical systems of the geometric description of chaos in systems of
differential equations

$\dot{x}^{i}=f^{i}(x^{1},\dots,x^{N})=f^{i}(\boldsymbol{x})$ (55)

that is, also in the case of dissipative systems. By differentiation with
respect to time of Eq.(55) we get a new system of equations

$\ddot{x}^{i}=\sum_{j=1}^{N}\frac{\partial f^{i}(\boldsymbol{x})}{\partial
x^{j}}\dot{x}^{j}=\sum_{j=1}^{N}\frac{\partial f^{i}(\boldsymbol{x})}{\partial
x^{j}}f^{j}(\boldsymbol{x})$ (56)

that can be derived from the Lagrangian function

$L(\boldsymbol{x},\boldsymbol{\dot{x}})=\sum_{i=1}^{N}[{\dot{x}}^{i}-f^{i}(\boldsymbol{x})]^{2}$
(57)

and the usual Lagrange equations. To this Lagrangian $\displaystyle
L(\boldsymbol{x},\boldsymbol{\dot{x}})$ one associates a metric function
homogeneous of degree one in the velocities

$\Lambda(x^{a},{\dot{x}}^{a})=L(x^{i},{\dot{x}}^{i}/{\dot{x}}^{0}){\dot{x}}^{0}\
,\hskip 14.22636pta=0,1,\dots,N;\ i=1,\dots,N$ (58)

involving an extra velocity $\displaystyle\dot{x}^{0}$; through this metric
function a metric tensor expressed as

$g_{ab}(\boldsymbol{x},\boldsymbol{\dot{x}})=\frac{1}{2}\frac{\partial^{2}\Lambda^{2}}{\partial\dot{x}^{a}\partial\dot{x}^{b}}$
(59)

provides the tangent bundle of the configuration space of the system (55) with
a Finslerian structure. The geodesics of this space, minimizing the functional
$\displaystyle\int_{\tau_{0}}^{\tau_{1}}\Lambda(x^{a},{\dot{x}}^{a})d\tau$,
are given by marco ; rund

$\frac{d^{2}x^{a}}{ds^{2}}+\gamma^{a}_{bc}(\boldsymbol{x},\boldsymbol{\dot{x}})\frac{dx^{b}}{ds}\frac{dx^{c}}{ds}=0$
(60)

where $\displaystyle\gamma^{a}_{bc}(\boldsymbol{x},\boldsymbol{\dot{x}})$ are
the connection coefficients derived from the velocity dependent metric
$\displaystyle g_{ab}(\boldsymbol{x},\boldsymbol{\dot{x}})$, and coincide with
the solutions of Eqs.(56). Then a geodesic deviation equation is defined also
on Finsler manifolds and relates stability/instability of the geodesics with
the curvature properties of the space marco . This approach certainly deserves
to be investigated to tackle chaotic dynamics of dissipative systems with the
same methodological approach successfully applied to Hamiltonian systems.

Acknowledgments

M.P. participated in this work within the framework of the project MOLINT
which has received funding from the Excellence Initiative of Aix-Marseille
University - A*Midex, a French “Investissements d’Avenir” programme.

## References

  * [1] D. V. Anosov.  Geodesic flows on closed Riemannian manifolds with negative curvature.  Proc. Steklov Math. Inst., 90:1–235, 1967.
  * [2] L. Casetti.  Efficient symplectic algorithms for numerical simulations of Hamiltonian flows.  Physica scripta, 51(1):29, 1995.
  * [3] L. Casetti, C. Clementi, and M. Pettini.  Riemannian theory of Hamiltonian chaos and Lyapunov exponents.  Physical Review E, 54(6):5969, 1996.
  * [4] M. Cerruti-Sola, R. Franzosi, and M. Pettini.  Lyapunov exponents from geodesic spread in configuration space.  Physical Review E, 56(4):4872, 1997.
  * [5] M. Cerruti-Sola and M. Pettini.  Geometric description of chaos in two-degrees-of-freedom Hamiltonian systems.  Physical Review E, 53(1):179, 1996.
  * [6] C. Clementi and M. Pettini.  A geometric interpretation of integrable motions.  Celestial Mechanics and Dynamical Astronomy, 84(3):263–281, 2002\.
  * [7] E. Cuervo-Reyes and R. Movassagh.  Non-affine geometrization can lead to non-physical instabilities.  Journal of Physics A: Mathematical and Theoretical, 48(7):075101, 2015.
  * [8] L. Di Cairano, M. Gori, and M. Pettini.  Coherent Riemannian-geometric description of Hamiltonian order and chaos with Jacobi metric.  Chaos: An Interdisciplinary Journal of Nonlinear Science, 29(12):123134, 2019.
  * [9] L. P. Eisenhart.  Dynamical trajectories and geodesics.  Annals of Mathematics, pages 591–606, 1929.
  * [10] J. Guckenheimer and P. Holmes.  Nonlinear oscillations, dynamical systems and bifurcations of vector fields.  Appl. Math. Sci. Series, 42, 1983.
  * [11] M. Hénon and C. Heiles.  The applicability of the third integral of motion: some numerical experiments.  The Astronomical Journal, 69:73, 1964.
  * [12] N. S. Krylov.  Works on the foundations of statistical physics.  Princeton Univ. Press, 1979.
  * [13] A. Lichnerowicz and T Teichmann.  Théories relativistes de la gravitation et de l’électromagnétisme.  PhT, 8(10):24, 1955.
  * [14] A. J. Lichtenberg and M. A. Lieberman.  Regular and chaotic dynamics.  Springer-Verlag, Berlin, 1992.
  * [15] R. Livi, M. Pettini, S. Ruffo, and A. Vulpiani.  Chaotic behavior in nonlinear Hamiltonian systems and equilibrium statistical mechanics.  Journal of statistical physics, 48(3-4):539–559, 1987.
  * [16] The natural and elegant geometric setting of Hamiltonian dynamics is provided by symplectic geometry. This geometrical framework is very powerful to study, for example, symmetries. However, symplectic manifolds are not endowed with a metric, and without a metric we do not know how to measure the distance between two nearby phase space trajectories and thus to study their stability/instability through the time evolution of such a distance.
  * [17] M. Pettini.  Geometrical hints for a nonperturbative approach to Hamiltonian dynamics.  Physical Review E, 47(2):828, 1993.
  * [18] M. Pettini.  Geometry and topology in Hamiltonian dynamics and statistical mechanics, volume 33.  Springer Science & Business Media, 2007.
  * [19] M. Pettini and R. Valdettaro.  On the Riemannian description of chaotic instability in Hamiltonian dynamics.  Chaos: An Interdisciplinary Journal of Nonlinear Science, 5(4):646–652, 1995.
  * [20] H. Poincaré.  Les méthodes nouvelles de la mécanique céleste, volume 3.  Blanchard, Paris, 1987.
  * [21] H. Rund.  The differential geometry of Finsler spaces, volume 101.  Springer Science & Business Media, 2012.
  * [22] S Wiggins.  Global bifurcations and Chaos.  Applied Mathematial Sciences, 73, 1988.


# Scattering Transparency of Clouds in Exoplanet Transit Spectra

Bhavesh Jaiswal1,2 and Tyler D Robinson3,4,5,6
1Space Astronomy Group, U R Rao Satellite Centre, Bengaluru 560037, India
2Dept. of Physics, Indian Institute of Science, Bengaluru 560012, India
3Department of Planetary Sciences; Lunar and Planetary Laboratory, University
of Arizona, Tucson, AZ 85721, USA
4Department of Astronomy and Planetary Science, Northern Arizona University,
Flagstaff, AZ 86011, USA
5Habitability, Atmospheres, and Biosignatures Laboratory, University of
Arizona, Tucson, AZ 85721, USA
6NASA Nexus for Exoplanet System Science Virtual Planetary Laboratory,
University of Washington, Box 351580, Seattle, WA 98195, USA E-mail:
<EMAIL_ADDRESS>
(Accepted XXX. Received YYY; in original form ZZZ)

###### Abstract

The presence of aerosols in an exoplanet atmosphere can veil the underlying
material and can lead to a flat transmission spectrum during primary transit
observations. In this work, we explore forward scattering effects from super-
micron sized aerosol particles present in the atmosphere of a transiting
exoplanet. We find that the impacts of forward scattering from larger aerosols
can significantly impact exoplanet transits and the strength of these effects
can be dependent on wavelength. In certain cloud configurations, the forward-
scattered light can effectively pass through the clouds unhindered, thus
rendering the clouds transparent. The dependence of the aerosol scattering
properties on wavelength can then lead to a positive slope in the transit
spectrum. These slopes are characteristically different from both Rayleigh and
aerosol absorption slopes. As examples, we demonstrate scattering effects for
both a rocky world and a hot Jupiter. In these models, the predicted spectral
slopes due to forward scattering effects can manifest in the transit spectrum
at the level of $\sim$10s to $\sim$100s of parts per million and, hence, could
be observable with NASA’s James Webb Space Telescope .

###### keywords:

Exoplanets — Radiative transfer

††pubyear: 2023††pagerange: Scattering Transparency of Clouds in Exoplanet
Transit Spectra–Scattering Transparency of Clouds in Exoplanet Transit Spectra

##  1 Introduction

Transiting exoplanets provide an opportunity to study the atmospheric
structure of these distant worlds through transit spectroscopy. The small
amount of stellar light which passes through the narrow annulus of the
exoplanet atmosphere in the transit geometry interacts with the atmosphere and
thus can carry information about the atmospheric physical and chemical state
(e.g. Seager & Sasselov, 2000; Brown, 2001; Hubbard et al., 2001).
Spectroscopic signatures of the atmosphere for a transiting exoplanet are
subtle, in part due to the small size of planets as compared to stars and due
to the relatively small covering area of the atmospheric annulus as compared
to the planetary disk. In the last two decades, exoplanet transit spectroscopy
has moved from noisy detections of simple atmospheric constituents (e.g.
Charbonneau et al., 2002; Tinetti et al., 2007; Swain et al., 2008) to, more
recently, searches for trends (e.g. Crossfield, 2015; Sing et al., 2016) in
the transit spectra of distinct types of worlds. Promisingly, early results
from James Webb Space Telescope (JWST) already show detections of atmospheric
species from transit spectra (JWST Transiting Exoplanet Community Early
Release Science Team et al., 2023; Alderson et al., 2023; Rustamkulov et al.,
2023).

$\begin{array}[]{rl}\lx@intercol\hfil\includegraphics[width=390.25534pt]{forward_scattering_illustration.png}\hfil\lx@intercol\\\
\\\ \\\ \\\
\includegraphics[width=216.81pt]{scattering_phase_function_0p6micron_wavelength.png}&\includegraphics[width=216.81pt]{scattering_phase_function_5micron_wavelength.png}\\\
\end{array}$

Figure 1: [Top] Cloud scattering in an exoplanet atmosphere. Strong forward
scattering from the cloud particles is shown for the slant incident direction.
Shorter wavelengths are more-strongly scattered in the forward direction than
longer wavelengths.[Bottom] The Mie scattering phase function for various
sizes ($r_{\mathrm{eff}}$, in legends) of water cloud particles shown at two
wavelengths: $\lambda$=0.6 $\upmu$m (left) and $\lambda$=5 $\upmu$m (right)

$\begin{array}[]{rl}\includegraphics[width=216.81pt]{Mie_vs_HG.png}&\includegraphics[width=216.81pt]{single_scattering_albedo.png}\\\
\end{array}$

Figure 2: [Left] Scattering phase function obtained from Mie theory and its
approximation by the Henyey-Greenstein function. $P_{\mathrm{Mie}}$ is
calculated at a wavelength of 0.6 $\upmu$m for water ice particle with
$r_{\mathrm{eff}}$ of 20 $\upmu$m. $P_{\mathrm{HG}}$ is calculated using Eq. 5
for $\theta_{s}$ of $None$. [Right] Single scattering albedo ($\omega$) for
H2O ice and Mg2SiO4 cloud particles.

Unfortunately, only a minority of the thousands of currently-known exoplanets
are suitable for atmospheric reconnaissance. In many of these cases, existing
observations indicate the presence of clouds or hazes (Pont et al., 2008; Sing
et al., 2009; Bean et al., 2010; Knutson et al., 2014; Benneke et al., 2019;
Wakeford et al., 2019; JWST Transiting Exoplanet Community Early Release
Science Team et al., 2023). Generally, these aerosols have presented
themselves either directly as slopes of decreasing transit depth with
increasing wavelength (potentially characteristic of sub-micron particles
(e.g. Lecavelier Des Etangs et al., 2008)) or somewhat indirectly through a
continuum absorption that suppresses gas absorption features (e.g. Kreidberg
et al., 2014). There have been limited observations of wavelength dependent
extinction (Benneke et al., 2019) which can also arise due to clouds. The
composition of exoplanet clouds and hazes largely remains a mystery as any
characteristic spectral signatures are either weak or outside the spectral
range of observations with the pre-JWST facilities, albeit the early broadband
JWST results are encouraging (Miles et al., 2022). Nevertheless, studies have
predicted the kinds of spectral slopes or absorption features that could be
expected for hazes in hot Jupiter atmospheres (Wakeford & Sing, 2015; Pinhas &
Madhusudhan, 2017; Ohno & Kawashima, 2020). We find that the studies which
have investigated the forward scattering from clouds so far (see for e.g.
Sengupta et al., 2020) have almost always considered sub-micron sized
scatterers which have a characteristic negative spectral slope resembling
Rayleigh-like scattering. The modeling investigations, on the other hand,
spanning hot Jupiters to sub-Neptunes have shown the potential for higher-
altitude clouds and hazes at sizes extending beyond 1 $\upmu$m and maybe even
beyond 10$\upmu$m (Powell et al., 2018; Charnay et al., 2021). Additionally,
the clouds on temperate rocky worlds can also be expected to grow much beyond
1 $\upmu$m in resemblance to clouds on Earth (see for e.g. Donovan et al.,
2001). In such a scenarios, it is imminent to explore the spectral effect of
forward scattering from super-micron sized aerosols.

Unlike smaller (i.e., sub-micron haze) particles, the super-micron sized
aerosol particles can have very strong forward scattering at wavelengths where
the size parameter of particles ($2\pi r_{\mathrm{eff}}/\lambda$) is much
larger than unity. This scattering behaviour is typical of Mie scattering.
Depending upon the clouds and the transit configuration, this forward
scattering can be so strong that deeper atmospheric layers within the cloud
can be probed, as compared to a more-isotropically scattering aerosol scenario
(de Kok & Stam, 2012; Robinson, 2017; García Muñoz & Cabrera, 2018). Thus, the
upper reaches of some clouds may appear transparent at shorter wavelengths but
transition to opaque conditions at longer wavelengths as the size parameter
changes from much larger than unity to of order (or smaller than) unity, as
predicted in Robinson (2017). This variation in transparency can occur over a
large spectral range, making the effect more-straightforward to detect.

In what follows, we investigate the spectral scattering effects of aerosols
for transiting exoplanets. For the majority of the simulations considered here
we omit gas absorption to focus on the impacts on the continuum set by the
cloud/haze. This work doesn’t seek to delineate between clouds and hazes and
explores a range of aerosol profiles. For convenience, though, we use the term
"clouds" throughout. In Section 2 of this paper we present the theory of
forward scattering for large cloud particles in the transit geometry. Section
3 discusses a Monte Carlo radiative transfer model which is used here for
more-detailed simulations. In Section 4 we present the results of simulations
for prescribed cloud structures in model atmospheres for a temperate planet
TRAPPIST-1e and a hot Jupiter HD 189733b. Finally, we briefly discuss the
feasibility and limitations of the simulations in section 5 and conclude in
section 6.

##  2 Forward Scattering Theory

Mie scattering from particles whose size parameter is larger than the
wavelength has a characteristic scattering phase function where a significant
amount of light is scattered in the forward direction (Hansen & Travis, 1974).
The scattering phase function is extremely large at $None$ scattering angle
(i.e., the forward-scattering direction) and decreases rapidly for larger
angles. Figure 1 shows the Mie scattering-derived phase function of water
droplets of various sizes at wavelengths of 0.6 $\upmu m$ and 5 $\upmu m$
where the sharp rise of the scattering phase function at forward scattering
angles is clearly seen for shorter wavelengths. Importantly, the forward-
scattering angle is probed in the exoplanet transit geometry where stellar
photons can be scattered in the same direction of travel from the limb of the
exoplanet towards the direction of the observer (de Kok & Stam, 2012). As is
apparent in Figure 1, light at shorter wavelengths is forward-scattered more
effectively than at the longer wavelengths. The angular extent of the star
apparent from the planet is an important parameter here as the incident light
is sampled from all the parts of the stellar disc. The forward scattering
phase function essentially gets convolved over the angular size of the star
($\sim$$$ for HD 189733b and $\sim$$$ for TRAPPIST-1e system).

$\begin{array}[]{rl}\includegraphics[width=216.81pt]{fxomega_waterice2.png}&\includegraphics[width=216.81pt]{g_lambda_waterice.png}\\\
\includegraphics[width=216.81pt]{fxomega_forsterite2.png}&\includegraphics[width=216.81pt]{g_lambda_forsterite.png}\\\
\end{array}$

Figure 3: Spectral variations of $f\times\omega$ and $g$ for various
$\theta_{\mathrm{s}}$. [Top] For TRAPPIST-1e with water ice particles with
$r_{\mathrm{eff}}$=20 $\upmu$m. [Bottom] For HD 189733b with silicate
particles with $r_{\mathrm{eff}}$=5 $\upmu$m. Also marked are the location of
$\theta_{s}$ corresponding to TRAPPIST-1e and HD 189733b.

For applications in planetary atmospheres, the shape of the scattering phase
function ($P(\theta)$; where $\theta$ is the scattering angle) is often
approximated with a Henyey-Greenstein (H-G) scattering phase function, whose
analytic nature makes it useful for radiative transfer calculations. The H-G
phase function (see e.g. Vasquez, M. et al., 2013; Millar-Blanchaer et al.,
2015; Robinson, 2017, for its application to exoplanets) is characterised by a
single parameter — $g$, which controls the shape of the phase function and the
strength of the forward peak — and is given by,

$P_{\rm
HG}(\theta,g)=\frac{1}{4\pi}\frac{1-g^{2}}{\left(1+g^{2}-2g\cos\theta\right)^{3/2}}\
.$ (1)

It is worth noting here that both, the HG parameter $g$ and the phase function
$P_{HG}$ are wavelength specific. Following Robinson et al. (2017)), aerosol
forward scattering can be incorporated into the transmission of stellar
intensity along the limb of a transiting exoplanet by writing the equation of
transfer as,

$\frac{dI}{d\tau_{\mathrm{eff}}}=-I\ ,$ (2)

where I is the radiant intensity along the path and $d\tau_{\mathrm{eff}}$ is
the differential “effective” aerosol optical depth along the slant path. Now,
we define a correction factor that is calculated by integrating the scattering
phase function in the forward direction,

$f=2\pi\int_{0}^{\theta_{\mathrm{S}}}P(\theta)\sin(\theta)\,d(\theta).$ (3)

This correction factor is dictated by the angular extent of the star
($\theta_{\mathrm{S}}$) as seen from the planet (García Muñoz & Cabrera,
2018). Considering this correction factor, effective optical depth is designed
to account for a lossless process where slant-path light that interacts with
the medium is simply scattered into the forward direction, with,

$d\tau_{\mathrm{eff}}=(1-f\omega)d\tau_{\mathrm{ext}}\ ,$ (4)

where $\omega$ is the single scattering albedo. Note that Equation 2 ignores
the multiple-scattering source term which comes from the isotropic component
of the phase function and is negligible for transiting geometry (explained in
Robinson et al., 2017). As seen in Figure 1, the scattering phase function can
vary drastically with wavelength, thereby leading to variations in both $f$
and $d\tau_{\mathrm{eff}}$ across the spectral range. Shorter wavelengths, due
to a steeper rise of the phase function, can have a significantly smaller
$d\tau_{\mathrm{eff}}$ causing a transparency in the clouds not encountered at
longer wavelengths. The wavelength-dependent factor $f\times\omega$ is a
measure of this transparency effect.

Alongside the more-intuitive analytic theory described above, we explore more-
detailed scattering solutions using a Monte Carlo-based approach. We adopt the
publicly-available scaTran package (discussed in Robinson, 2017), which
includes scattering, absorption, and refraction effects for exoplanet
atmospheres in the transiting geometry. Mie scattering by aerosol particles is
considered using the Monte Carlo (MC) algorithm by sampling the H-G phase
function, where the wavelength-dependent H-G asymmetry parameter is designed
to reproduce the forward scattering peak from the Mie-derived phase function.
Thus, at each wavelength we determine the value of $g$ that satisfies,

$\int_{0}^{\theta_{S}}P_{\rm
Mie}(\theta)\sin\theta\,d\theta=\int_{0}^{\theta_{S}}P_{\rm
HG}(\theta,g)\sin\theta\,d\theta\ ,$ (5)

where $P_{\mathrm{Mie}}$ is calculated using Mie theory. Figure 2 shows
$P_{\mathrm{Mie}}$ and $P_{\mathrm{HG}}$ for a water cloud particle (with
$r_{\mathrm{eff}}$ = 20 $\upmu$m) which satisfies Equation 5. It is noteworthy
that our choice of ’g’ parameter depends upon integrating the phase function
for the stellar angular size ($\theta_{S}$) with an intention to preserve the
total scattering in the direction of the star. As is seen in the figure, this
HG phase function underestimates the forward scattering by an order or
magnitude at $\sim$$$ but it overestimates it proportionally at $\sim$$$ and
hence maintaining the total scattering in the direction of the star. Here, a
choice of a multi-parameter HG phase function might be more useful for a
better fit to the forward part of the Mie scattering phase function. However
it will add more parameters (and complexity) in our model and hence we
restrict the present analysis to a single parameter HG function. The single
scattering albedo, $\omega$, is also shown in Figure 2 for the water and
silicate clouds. The optical constants for the scattering calculations of the
clouds are considered from Budaj et al. (2015).

##  3 Radiative Transfer Simulations

We prescribe a simple cloud model for full-scale MC simulations, as is shown
in Figure 1. This model contains a layer of cloud particles spread uniformly
around the terminator. In the transit geometry, we explore the effect of
clouds mainly with three parameters: cloud base height ($Z_{\mathrm{c}}$),
cloud vertical extent ($\Delta{Z}_{\mathrm{c}}$) and cloud slant scattering
optical depth ($\tau$, at $\lambda$=0.6 $\upmu$m), which are all depicted in
Figure 1. We study two cases of planets: a rocky world - TRAPPIST-1e (Gillon
et al., 2017), with water ice clouds and a hot Jupiter - HD 189733b (Bouchy et
al., 2005), with silicate clouds. The clouds can exist in a variety of
conditions on temperate as well as hot jupiter planets. For simplicity, we
consider uniformly distributed clouds having a log-normal size distribution
which is given as

$n(r)=\frac{1}{r}\exp\left(-0.5\frac{(\ln{r}-\ln{r_{g}})^{2}}{(\ln{\sigma})^{2}}\right)\
.$ (6)

Here, $n(r)$ is the number density of cloud particles with radius $r$ where
the modal radius ($r_{g}$) and width of the distribution ($\sigma$) are
related to effective radius $r_{\mathrm{eff}}$ and effective variance
$v_{\mathrm{eff}}$ as:

$r_{\mathrm{eff}}=\frac{r_{g}}{\exp(-2.5(\ln{\sigma})^{2})}\ ,$ (7)

$v_{\mathrm{eff}}=\exp((\ln{\sigma})^{2})-1\ .$ (8)

The cloud vertical optical depth ($\tau$) is considered as decreasing
exponentially above the cloud base ($Z_{c}$) i.e. for $z>Z_{c}$

$\tau(z)=\tau(Z_{c})\exp{\left(-\frac{z}{Z_{H}}\right)}\ ,$ (9)

where $Z_{\mathrm{H}}$ is the scale height of the cloud. Below the cloud base
the cloud optical depth is zero.

In order to demonstrate the effect, a fiducial value of $r_{\mathrm{eff}}$ of
5 $\upmu$m (similar to silicate clouds of Powell et al., 2018) and 20 $\upmu$m
(similar to water ice clouds of Donovan & Lammeren, 2001) is used for HD
189733b silicate clouds and TRAPPIST-1e water ice clouds respectively. The key
model parameters adopted in the simulations are given in Table 1.

Next, the evolution of spectral slope is discussed using two approaches – the
analytical approach (eq. 2-3) and a full scale 3D Monte Carlo (MC) scattering
approach of scaTran. It is noteworthy that the analytical approach ignores the
multiple scattering of photons whereas the MC approach does not. However, the
effect of multiple scattering term is negligible in the analytical approach
and, as we show later, the simulation results of both the approaches lead to a
positive slope in the spectrum.

Figure 4: Manifestation of the scattering transparency effect in the transit
spectrum of TRAPPIST-1e system having a water cloud with $\tau$=1 and
$Z_{\mathrm{c}}$=100 km. The transit spectrum is simulated with MC
calculations for cloud particles having an $r_{\mathrm{eff}}$=20 $\upmu$m. The
dashed lines are the spectra obtained for constant values of $g$. The actual
spectra, after considering the spectral variations of $g$, are shown in bold
lines. They are obtained by interpolating among the dashed spectra. The blue
curves do not consider the spectral variations of $\omega$ whereas the red
curves consider the realistic variation $\omega$ for water clouds which leads
to the spectral features of water droplets. Gas absorption is not considered.
The spectral rise towards shorter wavelengths in the dashed curves is due to
Rayleigh scattering. A color gradient for $f\times\omega$, which is a measure
of scattering transparency, is shown in the background. This gradient of
$f\times\omega$ is same as that shown in figure 3 for TRAPPIST-1e system.

###  3.1 The Analytical approach

The spectral variation of the theoretically calculated parameters
$f\times\omega$ and $g$ is shown in Figure 3. Results are shown for the water
ice cloud and silicate cloud particle for a range of $\theta_{s}$ spanning
that of TRAPPIST-1e to HD 189733b. The variation of $f\times\omega$ and $g$
with wavelength as well as angular size ($\theta_{s}$) gives an intuitive
understanding of the results via eq. 1-3. Significantly larger values of
$f\times\omega$ and $g$ at shorter wavelengths indicate the importance of a
wavelength-dependent scattering transparency effect. It arises due to the
sharp increase of the size parameter ($2\pi r_{\mathrm{eff}}/\lambda$) towards
short wavelengths which leads to a strong forward scattering phase function.
Faintly visible vertical stripes in the $f\times\omega$ plot are a result of
the spectral variation in $\omega$, as seen in Figure 2. The vertical stripes
in the $g$ plots are a result of the variation in the real and imaginary part
of the refractive index of the cloud particles. The variation of
$f\times\omega$ with wavelength is essentially the variation in the effective
optical depth $d\tau_{\mathrm{eff}}$ (via eq. 4) which manifests as a spectral
slope in the transit spectrum.

###  3.2 The Monte Carlo approach

After creating the model, the values of $g$ and $\omega$ are calculated for
the given $r_{\mathrm{eff}}$ used in the model for all the wavelengths (as
shown in Figure 3) . The parameters $\tau$, $Z_{\mathrm{c}}$ and
$Z_{\mathrm{H}}$ are left free for exploration. The cloud vertical extent
($\Delta{Z}_{\mathrm{c}}$) is always fixed at 4$\times Z_{\mathrm{H}}$.
Multiple instances of the scaTran model are run, each time for a fixed value
of $g$. Each instance takes about $\sim 10$ minutes to run on a quad-core
machine. Then, to produce a spectrum having a variation in $g$ across
wavelengths, we interpolate on these spectra. For TRAPPIST-1e, the values of
$g$ varies from 0.987 (at $\lambda$=0.5 $\upmu$m) to 0.81 (at $\lambda$=20
$\upmu$m) for $r_{\mathrm{eff}}$ of 20 $\upmu$m and for HD 189733b this
variation is from 0.98 (at $\lambda$=0.5 $\upmu$m) to 0.55 (at $\lambda$=20
$\upmu$m) for $r_{\mathrm{eff}}$ of 5 $\upmu$m.

Figure 4 shows an example of how the wavelength-dependent scattering
transparency is evolved from the interpolating the spectra, for the case of
TRAPPIST-1e. The spectra obtained with various constant values of $g$ (as
would be done without considering the scattering transparency effect) are
shown with dashed lines. Considering the variation of $g$ with wavelength, the
interpolated spectrum is shown with bold lines. As the transit depth is much
lower at shorter wavelengths than at longer wavelengths the spectrum shows a
positive slope in wavelength. This spectral shape closely resembles the
variation of $f\times\omega$ with wavelength (shown as a color gradient in the
background). To demonstrate the case of a scattering cloud without any
spectral variation in scattering albedo ($\omega$), we show another spectrum
(in blue) which has $\omega$ set to a fixed value of unity across the spectral
range. Here too, we see the similar variation of transit depth in the spectrum
(i.e. smaller transit depths at shorter wavelengths). In general, it can be
said that the clouds can be more transparent at shorter wavelengths than at
longer wavelengths if the scattering particles are sufficiently large with
respect to the wavelength.

Table 1: Key Model Parameters | TRAPPIST-1e | HD 189733b
---|---|---
Model base pressure (Pa) | $10^{5}$ | $10^{6}$
Model top pressure (Pa) | $10^{-5}$ | $10^{-4}$
Isothermal Temperature (K) | 275 | 1500
Atmospheric Composition | N2 | H2 & He (0.85:0.15)
Cloud | H2O | Mg2SiO4 (Forsterite)
Cloud $r_{\mathrm{eff}}(v_{\mathrm{eff}})$ | $20\upmu m(2)$ | $5\upmu m(0.5)$
Layers | 180 | 126
Gravity (m s-2) | 8.9 | 21.8
Star radius (R$\odot$) | 0.12 | 0.78
Orbital radius (AU) | 0.029 | 0.03
Planet radius (R$\bigoplus$) | 0.92 | 12.75

###  3.3 The planet model

As discussed in sections above, the impact of this effect on the overall
spectrum has a dependence on $r_{\mathrm{eff}}$ of cloud and also
$\theta_{\mathrm{star}}$. The fixed parameters of the model for both the
planets are presented in table 1 for the case of TRAPPIST-1e and HD 189733b.
The stellar and planetary radius for HD 189733b system is taken from Gaia
Collaboration et al. (2018) and Bonomo et al. (2017) respectively. Whereas,
those for TRAPPIST-1e system are taken from Agol et al. (2021). The planets
are considered to have an vertically isothermal atmosphere with the gravity
and radius defined at the altitude of model base pressure. The atmosphere is
divided in 180 and 126 vertical layers respectively for TRAPPIST-1e and HD
189733b. The choice of the number of layers is mainly driven by the need to
vertically resolve the extent of cloud having the smallest scale height.

##  4 Results

The scattering transparency effect can be produced by a model which can sample
the forward scattering phase function (such as the backward Monte Carlo
approach of scaTran) and not by models which work by simple geometric ray
tracing unless we factor in the non-grey variation of optical depth as
mentioned in eq. 4 (see section 4.1). As discussed in Robinson (2017) the $MC$
approach considers a Monte Carlo based multiple scattering in the clouds
whereas the $geometric$ approach considers the cloud scattering optical depth
as extinction along straight rays. First we show the spectrum produced by both
the methods i.e. $MC$ approach and $geometric$ approach (without considering
the non-grey variation in optical depth). Later, we show that considering a
variation in optical depth can also lead to a scattering transparency
spectrum.

The results for various configurations of clouds are presented in figure 5.
The top panel shows the cloud spectrum for various cloud scattering optical
depths ($\tau$) for the TRAPPIST-1e and HD 189733b configurations. Here, we do
not consider any gas absorption in the medium and show purely the effect of
clouds on the spectrum and the spectral features are due to spectral
variations of $\omega$ (as in figure 2). However, we do include the Rayleigh
scattering from the background gas for both the planets and also the H2 and He
collision induced absorptions for HD 189733b. The cloud scale height is
considered to be 1 pressure scale height for both the cases but the cloud base
altitude ($Z_{\mathrm{c}}$) is at 100 km ($\sim 10$ microbar) and 3000 km
($\sim 0.1$ millibar) for TRAPPIST-1e and HD 189733b respectively. The
difference in the $geometric$ and $MC$ cases is clearly observable. The
$geometric$ cases mostly produce a flat continuum (except the negative slopes
at shorter wavelengths which is due to the Rayleigh scattering from gas)
whereas in the $MC$ case the continuum has a positive slope due to the
scattering transparency effect. This effect is more pronounced in the hot
Jupiter configuration for two reasons: firstly, a larger planet overlap area
with the star offers large proportion of scattering particles and secondly, a
larger $\theta_{\mathrm{s}}$ for hot Jupiters allows sampling from extended
regions of the star (see the variation of $f\times\omega$ with
$\theta_{\mathrm{s}}$ in figure 3). Though the scattering transparency effect
acts to reduce the transit depth towards shorter wavelengths, at very short
wavelengths the Rayleigh scattering by molecules starts to increase the
transit depth. For smaller values of $\tau$ (thin clouds), there is a
competition between the two processes which can lead to an interesting ’U’
shape in the spectrum at around $\sim$1 $\upmu$m. The color gradient shown in
the background of the figures is the theoretically estimated values of
$f\times\omega$ which signifies the strength of this effect. One can notice
small variations in $f\times\omega$ over the spectrum accompanying the peaks
in the spectrum which is due to the variations of $\omega$.

$\begin{array}[]{rl}\includegraphics[width=216.81pt]{Trappist_spectra.png}&\includegraphics[width=216.81pt]{HD189733b_spectra.png}\\\
\lx@intercol\hfil\includegraphics[width=346.89731pt]{HD189733b_spectrum_gas_0p03AU_CF.png}\hfil\lx@intercol\end{array}$

Figure 5: Scattering transparency effect for various conditions. [Top]
Spectral behaviour of this effect for TRAPPIST-1e and HD 189733b geometry.
Results are shown for different slant optical depths $\tau$ for the
$geometric$ (dashed lines) and $MC$ (solid lines) cases. [Bottom] The spectrum
of HD 189733b considering absorption by gas as well as clouds, for $\tau=3$.

$\begin{array}[]{rl}\includegraphics[width=216.81pt]{TD_diff_TRAPPIST1e.png}&\includegraphics[width=216.81pt]{TD_diff_HD189.png}\\\
\end{array}$

Figure 6: Variation in the transit depth (TD) difference due to various cloud
parameter are shown for TRAPPIST-1e (TD($\lambda$=5 $\upmu$m) -
TD($\lambda$=0.6 $\upmu$m)) and HD 189733b (TD($\lambda$=5 $\upmu$m) -
TD($\lambda$=1.6 $\upmu$m)). Figure 7: A comparison of transit spectra of HD
189733b simulated with analytical approach (dashed lines) and Monte Carlo
approach (solid lines) shown for different values of slant optical depth.

Finally, we show this effect after including the absorption and scattering by
both clouds as well as gas. For this, we include H2O gas absorption opacities
(Rothman et al., 2010) for HD 189733b. Considering the same conditions as in
panel (b) for the cloud $\tau$=3 case, we add H2O with a volume mixing ratio
of 0.0004 in all the layers. The resultant spectra are shown in the bottom
panel of figure 5 for $geometric$ and $MC$ cases. The $geometric$ case
produces an almost flat spectrum with H2O absorption features whereas the $MC$
case produces a scattering slope continuum with H2O absorption features twice
as deep. Towards shorter wavelengths, where scattering transparency is high
(high values of $f\times\omega$), the photons can penetrate through the cloud
layers when the forward scattering is strong. This allows photons to penetrate
to deeper layers where the gas absorption is high. The absorption features of
gas (H2O) absorption are seen in the spectrum. For comparison, we also show
the spectrum of a $cloudfree$ atmosphere. Removing the high altitude cloud
leads to deeper gas absorption bands in the transit spectrum as the photons
are able to travel to deeper regions of the atmosphere. It is noteworthy that
the positive spectral slope in the transit spectrum is caused by the H2O gas
opacity. A comparison of the absorption bands in the $MC$ case and the
$cloudfree$ case reveals how the scattering transparency effect, in $MC$ case,
can cause a gradual decrease in the strength of gas absorption from
visible/NIR to MIR regions. This variation in the gas absorption band across
the spectrum would also be an indicator of scattering transparency of clouds.

We explore parameter space of the cloud profiles for the scattering
transparency spectral slopes in figure 6. Considering the NIRSpec (onboard
JWST) spectral range of 0.6 $\upmu$m to 5 $\upmu$m, we consider the transit
depth difference at 5 $\upmu$m and 0.6 $\upmu$m as a potential metric for
predicting the slopes in the JWST observations. Other JWST instruments (i.e.,
MIRI) can observe at longer wavelengths and could provide additional
constraints on slopes, although noise levels can often be higher in such mid-
infrared observations. We explore the effect of clouds using 3 parameters:
Cloud base altitude, Cloud scale height and optical depth (slant $\tau$ at
cloud base, as in figure 1). Since all of these parameters are poorly-
constrained for the exoplanets and given the conditions in which the
atmospheres exist, we consider a large range in these parameters. The cloud
base is simulated for a range of altitudes. Some of the high altitude cloud
cases shown may seem unimaginable (such as clouds at high altitudes having
large scale heights) but inspired by a few observations of high altitude haze
(Montmessin et al., 2006; Estrela et al., 2021) we consider them for the
completeness of this study. The cloud scale heights are considered for a large
range of $\sim$ 0.1 to $\sim$ 2 $\times$ the atmospheric scale heights ($H$).
It is noteworthy that the cloud vertical profile is considered to be
exponential for all the cases and the maximum cloud extent is cutoff at
$4\times$ cloud scale height. In our approach to quantify the scattering
spectral slope with the difference of two spectral points in the transit
spectrum, we do not wish to contaminate it with any kind of absorption slopes.
So, to study purely the effect of scattering on cloud parameters, we do not
consider any absorption from either the cloud particles or the surrounding gas
(similar to $\omega=1$ case of figure 4). The results are shown for
TRAPPIST-1e and HD 189733b transit geometry. Due to strong Rayleigh scattering
at lower wavelengths of our HD 189733b model the lower wavelength point in the
Transit Depth difference plots is chosen to be 1.6 $\upmu$m instead of 0.6
$\upmu$m. We find that the clouds having, either high cloud base altitudes
($Z_{\mathrm{c}}$) or large scale heights ($Z_{\mathrm{H}}$) are favourable to
show stronger effects in the scattering transparency slope. That is because
increasing either $Z_{\mathrm{c}}$ or $Z_{\mathrm{H}}$ increases the overlap
area of the planetary clouds with the star in a transit geometry leading to an
increase in the proportion of the scattered photons which will eventually
contribute to the scattering transparency slope in the spectrum. Further, the
clouds having low optical depth and low scale heights are also detectable if
the cloud is at a higher altitude and has a sufficiently large spatial extent
around the terminator. This scenario increases the overlap area of the clouds,
because of the regions where stellar rays penetrate below the cloud base
altitude (i.e. for $z<Z_{\mathrm{c}}$). For $\tau$ larger than $\sim$10 and
$\sim$20 the spectrum remains flat for all the cases for TRAPPIST-1e and HD
189733b respectively. Within the considered parameter space, we observe that
the maximum transit depth difference can reach to no more than 50 ppm and 300
ppm for the two planet cases respectively. Given that the JWST error bars on
these observations are expected to be $\sim$10-20 ppm (Greene et al. (2016)),
these slopes could be observable with JWST on rocky worlds and hot Jupiters if
such clouds are present.

###  4.1 Scattering transparency spectrum with Analytical approach

As discussed earlier, we have shown the results of scattering transparency
slopes using a full scale 3D Monte Carlo scattering approach. This approach is
accurate but can be time consuming and can take several hours to produce one
spectrum. For this reason, it becomes useful to rely on analytical approach in
situations requiring quick results. Though analytical approach is not as
accurate as an MC calculation it can be a good replacement to MC calculations
in various situations. In order to demonstrate the analytical approach we
calculate transit spectrum for a non-grey cloud extinction optical depth. The
optical depth is varied as per eq. 4, across the spectrum. The scattering
transparency spectra obtained with this approach is shown in figure 7 along
with the spectra obtained with $MC$ approach. As is seen, the analytical
approach produces spectra which closely resemble to those of $MC$ approach and
reproduce the scattering transparency slopes. Though the trends in the spectra
resemble well, there is a minor difference in the spectra produce by two
approaches. This difference is minimum for very small ($\tau<1$) and for very
large ($\tau>1$) slant optical depths. For moderate optical depths
($\tau\sim$1), the difference can be as large as $\sim$100 ppm. It happens
because the proportion of cloud-scattered photons, reaching the stellar disc,
is maximised at moderate optical depths. This leads to an increased difference
in the analytical vs. $MC$ approach at moderate optical depths. Though a $MC$
approach to such calculations is an accurate method, we encourage the
community to make use of analytical approach for quick estimations.

##  5 Discussion

The larger and hotter planets such as hot Jupiters may be the best targets to
detect the effect of scattering transparency. On smaller terrestrial planets,
like TRAPPIST-1 planets, it may be difficult to detect the scattering
transparency effect unless the clouds are formed at unusually high altitudes
(with respect to terrestrial standards). This is also the reasoning behind our
choice of 100 km cloud base altitude for TRAPPIST-1e in figure 5. Here the
pressure at the cloud base altitude is $\sim 10$ microbars, which is
approximately same as for the high altitude carbon dioxide clouds of Mars.

The transit method can only probe the clouds in the terminator region of the
planet. So far we have considered uniform clouds on terminators however this
condition may not be met in all the targets. The presence of patchy clouds
(e.g. MacDonald & Madhusudhan, 2017), for example, can lead to a slope in the
spectrum which may not be as high as for the uniform clouds. Similarly, for
the planets which may have a difference in cloud properties in east vs. west
limb of the planet (see e.g. Powell et al., 2019) the scattering slopes may be
an aggregate of the two cases of clouds (i.e. on east and west limb).

As we demonstrate here the effect of non-grey Mie scattering, we recommend the
future studies of forward scattering from cloud particles to consider this
effect for simulations especially for super-micron sized cloud particles. One
can start by first estimating the correction term (1-$f\times\omega$) and
check its spectral dependence. If the spectral dependence of this term is
relatively large, one can either calculate the spectrum using analytical
approach or by a full Monte Carlo approach. It may be possible for future
studies to improve upon the present work by including scattering phase
function of non-spherical cloud grain and also by considering a two term
Henyey-Greenstein phase function for a better approximation of forward
scattering phase function.

##  6 Conclusion

Clouds can appear in a wider variety of conditions than considered here. But
with the specific conditions considered here, one can say that moderate slant
optical depth (1-10), large scale heights and larger $r_{\mathrm{eff}}$
increase our chances of detecting the spectral slopes for a given planet-star
system. Many of the assumptions considered here, like a log-normal particle
size distribution, width of the distribution, uniformity of clouds on East and
West limbs are chosen for brevity but we emphasize that, if all else being
equal, larger particles tend to amplify the scattering transparency at shorter
wavelengths. This leads to a positive slope in the broadband transit spectrum.
This slope exists for strongly forward scattering particles (usually with
$r_{\mathrm{eff}}\gtrsim 1$ $\upmu$m). For particles with $r_{\mathrm{eff}}<1$
$\upmu$m the transparency slope starts to reduce significantly as the
scattering phase function starts to shrink in the forward direction.

The detection of scattering transparency slopes in the transit spectrum will
be an independent indicator of the presence of the clouds having super-micron
sized aerosols. In this work, we uncover a viable framework which can be used
to predict the effect of scattering transparency over a broad spectral range.
Using this framework, it might be possible to infer the important physical
properties like $r_{\mathrm{eff}}$ and $\tau$ for the clouds. With JWST
operational now, we not only can probe the hot Jupiters with better
sensitivities but also, for a first time, a large number of cooler worlds. The
broadband transit observations of JWST for a variety of exoplanets will allow
us a chance to validate the predictions of scattering transparency slopes.
Clouds having large optical depths or smaller scale heights will eventually
lead to a flat spectrum irrespective of cloud particle size. However, given
the variety of exoplanet atmospheres which will be probed by JWST with its
large wavelength coverage, we envisage the detection of the scattering
transparency effect on several targets.

## Acknowledgement

We acknowledge many comments and suggestions from the reviewer which have
helped us to improve the quality of the manuscript. BJ acknowledges the
initial guidance and discussions with Nikku Madhusudhan which have led to this
investigation.

## Data Availability

The data underlying this article will be shared on reasonable request to the
corresponding author.

## References

  * Agol et al. (2021) Agol E., et al., 2021, The Planetary Science Journal, 2, 1
  * Alderson et al. (2023) Alderson L., et al., 2023, Nature, 614, 664
  * Bean et al. (2010) Bean J. L., Miller-Ricci Kempton E., Homeier D., 2010, Nature, 468, 669
  * Benneke et al. (2019) Benneke B., et al., 2019, Nature Astronomy, 3, 813
  * Bonomo et al. (2017) Bonomo A. S., et al., 2017, A&A, 602, A107
  * Bouchy et al. (2005) Bouchy F., et al., 2005, A&A, 444, L15
  * Brown (2001) Brown T. M., 2001, ApJ, 553, 1006
  * Budaj et al. (2015) Budaj J., Kocifaj M., Salmeron R., Hubeny I., 2015, MNRAS, 454, 2
  * Charbonneau et al. (2002) Charbonneau D., Brown T. M., Noyes R. W., Gilliland R. L., 2002, ApJ, 568, 377
  * Charnay et al. (2021) Charnay B., Blain D., Bézard B., Leconte J., Turbet M., Falco A., 2021, A&A, 646, A171
  * Crossfield (2015) Crossfield I. J. M., 2015, PASP, 127, 941
  * Donovan & Lammeren (2001) Donovan D., Lammeren A., 2001, Journal of Geophysical Research, 106, 27425
  * Donovan et al. (2001) Donovan D. P., et al., 2001, J. Geophys. Res., 106, 27,449
  * Estrela et al. (2021) Estrela R., Swain M. R., Roudier G. M., West R., Sedaghati E., Valio A., 2021, AJ, 162, 91
  * Gaia Collaboration et al. (2018) Gaia Collaboration et al., 2018, A&A, 616, A1
  * García Muñoz & Cabrera (2018) García Muñoz A., Cabrera J., 2018, MNRAS, 473, 1801
  * Gillon et al. (2017) Gillon M., et al., 2017, Nature, 542, 456
  * Greene et al. (2016) Greene T. P., Line M. R., Montero C., Fortney J. J., Lustig-Yaeger J., Luther K., 2016, ApJ, 817, 17
  * Hansen & Travis (1974) Hansen J. E., Travis L. D., 1974, Space Sci. Rev., 16, 527
  * Hubbard et al. (2001) Hubbard W. B., Fortney J. J., Lunine J. I., Burrows A., Sudarsky D., Pinto P., 2001, ApJ, 560, 413
  * JWST Transiting Exoplanet Community Early Release Science Team et al. (2023) JWST Transiting Exoplanet Community Early Release Science Team et al., 2023, Nature, 614, 649
  * Knutson et al. (2014) Knutson H. A., et al., 2014, ApJ, 794, 155
  * Kreidberg et al. (2014) Kreidberg L., et al., 2014, Nature, 505, 69
  * Lecavelier Des Etangs et al. (2008) Lecavelier Des Etangs A., Pont F., Vidal-Madjar A., Sing D., 2008, A&A, 481, L83
  * MacDonald & Madhusudhan (2017) MacDonald R. J., Madhusudhan N., 2017, MNRAS, 469, 1979
  * Miles et al. (2022) Miles B. E., et al., 2022, arXiv e-prints, p. arXiv:2209.00620
  * Millar-Blanchaer et al. (2015) Millar-Blanchaer M. A., et al., 2015, ApJ, 811, 18
  * Montmessin et al. (2006) Montmessin F., et al., 2006, Icarus, 183, 403
  * Ohno & Kawashima (2020) Ohno K., Kawashima Y., 2020, ApJ, 895, L47
  * Pinhas & Madhusudhan (2017) Pinhas A., Madhusudhan N., 2017, Monthly Notices of the Royal Astronomical Society, 471, 4355
  * Pont et al. (2008) Pont F., Knutson H., Gilliland R. L., Moutou C., Charbonneau D., 2008, MNRAS, 385, 109
  * Powell et al. (2018) Powell D., Zhang X., Gao P., Parmentier V., 2018, ApJ, 860, 18
  * Powell et al. (2019) Powell D., Louden T., Kreidberg L., Zhang X., Gao P., Parmentier V., 2019, ApJ, 887, 170
  * Robinson (2017) Robinson T. D., 2017, The Astrophysical Journal, 836, 236
  * Robinson et al. (2017) Robinson T. D., Fortney J. J., Hubbard W. B., 2017, ApJ, 850, 128
  * Rothman et al. (2010) Rothman L. S., et al., 2010, J. Quant. Spectrosc. Radiative Transfer, 111, 2139
  * Rustamkulov et al. (2023) Rustamkulov Z., et al., 2023, Nature, 614, 659
  * Seager & Sasselov (2000) Seager S., Sasselov D. D., 2000, ApJ, 537, 916
  * Sengupta et al. (2020) Sengupta S., Chakrabarty A., Tinetti G., 2020, ApJ, 889, 181
  * Sing et al. (2009) Sing D. K., Désert J. M., Lecavelier Des Etangs A., Ballester G. E., Vidal-Madjar A., Parmentier V., Hebrard G., Henry G. W., 2009, A&A, 505, 891
  * Sing et al. (2016) Sing D. K., et al., 2016, Nature, 529, 59
  * Swain et al. (2008) Swain M. R., Vasisht G., Tinetti G., 2008, Nature, pp 329–331
  * Tinetti et al. (2007) Tinetti G., et al., 2007, Nature, 448, 169
  * Vasquez, M. et al. (2013) Vasquez, M. Schreier, F. Gimeno García, S. Kitzmann, D. Patzer, B. Rauer, H. Trautmann, T. 2013, A&A, 557, A46
  * Wakeford & Sing (2015) Wakeford H. R., Sing D. K., 2015, A&A, 573, A122
  * Wakeford et al. (2019) Wakeford H. R., Wilson T. J., Stevenson K. B., Lewis N. K., 2019, Research Notes of the American Astronomical Society, 3, 7
  * de Kok & Stam (2012) de Kok R. J., Stam D. M., 2012, Icarus, 221, 517