# Moduli stabilization with bulk scalar in nested doubly warped braneworld
model

Arko Bhaumik<EMAIL_ADDRESS>Physics and Applied Mathematics Unit
Indian Statistical Institute, Kolkata-700108, India  Soumitra SenGupta
<EMAIL_ADDRESS>School of Physical Sciences
Indian Association for the Cultivation of Science, Kolkata-700032, India

###### Abstract

We examine the modulus stabilization mechanism of a warped geometry model with
nested warping. Such a model with multiple moduli is known to offer a possible
resolution of the fermion mass hierarchy problem in the standard model. A six
dimensional doubly warped braneworld model under consideration admits two
distinct moduli, with the associated warp factors dynamically generating
different physical mass scales on four 3-branes. In order to address the
hierarchy problem related to the Higgs mass, both moduli need to be stabilized
around their desired values without any extreme fine tuning of parameters. We
show that it is possible to stabilize them simultaneously due to the
appearence of an effective 4D moduli potential, which is generated by a single
bulk scalar field having non-zero VEVs frozen on the 3-branes. We also discuss
how the entire mechanism can possibly be understood from a purely
gravitational point of view, with higher curvature $f(R)$ contributions in the
bulk automatically providing a scalar degree of freedom that can serve as the
stabilizing field in the Einstein frame.

###### Contents

  1. 1 Introduction
  2. 2 Review of doubly warped braneworld model
  3. 3 Modulus stabilization in doubly warped model
    1. 3.1 Dynamics of the bulk field
    2. 3.2 The effective potential
    3. 3.3 Stabilizing $k$
    4. 3.4 Stabilizing $c$
    5. 3.5 Phenomenological implications
  4. 4 Insight from higher curvature gravity
  5. 5 Discussions

##  1 Introduction

The genesis of extra dimensional theories in physics dates back to the early
1920s, marked by Kaluza and Klein’s original attempt to unify classical
electrodynamics and Einstein gravity. In recent decades, higher-dimensional
braneworld models have attracted fresh interest [1, 2, 3, 4, 5, 6, 7, 8, 9,
10, 11, 12, 13, 14], arising as potential candidates for addressing the gauge
hierarchy problem which is inherently linked to the problem of the Higgs mass.
As the only fundamental scalar in the Standard Model (SM), the Higgs should
experience large radiative corrections to its mass due to its self-coupling
and couplings to other matter and gauge fields. These corrections are expected
to flow up to the ultraviolet cut-off scale of the underlying quantum field
theory. Protecting the Higgs mass from such Planck-scale corrections and
keeping it safely within the TeV scale requires an extremely precise fine
tuning of parameters, which leads to a naturalness problem within the SM. In
this sense, the discovery of a Higgs boson as light as
$125\>\textrm{GeV}/c^{2}$ [15, 16, 17] has simultaneously validated the last
major prediction of the SM and exposed its central conundrum.

The warped braneworld model proposed by Randall and Sundrum [9] offers an
elegant explanation of the large mass hierarchy, based on the premise of a
non-factorizable warped five-dimensional spacetime geometry. The extra
dimension is taken to be an orbifold with topology $S^{1}/\mathbb{Z}_{2}$. Two
3-branes are located at the opposite boundaries of the bulk, which is a slice
of $\textrm{AdS}_{5}$. One of these branes (the “TeV brane”) is identified
with our visible universe, on which all the SM fields are assumed to be
confined to. Gravity alone is assumed to permeate through the bulk. This
places the apparent scale of gravity $(M_{Pl})$ on the visible brane very
close to the fundamental scale $(M_{5})$. In contrast, the energy density in
the bulk induces an exponential warping along the extra dimension, which
causes the physical mass $m_{H}^{(ph)}$ of the Higgs on the TeV brane to be
suppressed from its fundamental mass $m_{H}\sim O(M_{5})$ according to
$m_{H}^{(ph)}=m_{H}e^{-\pi kr_{c}}$, where $k\sim M_{5}$ and $r_{c}$ is the
brane separation (which plays the role of the extra dimensional modulus). For
$kr_{c}\sim 11.54$, the electroweak scale is dynamically generated from the
fundamental scale. But the original model contains no mechanism to stabilize
the magnitude of $kr_{c}$ around its desired value. This problem was addressed
by Goldberger and Wise [18], who showed that $kr_{c}$ can be stabilized
appropriately by introducing a massive scalar field $\phi$ in the bulk with
quartic interactions localized on the branes. Integrating out the extra
dimensions results in an effective potential $V_{eff}(r_{c})$ which stabilizes
the modulus without any excessive fine tuning. This mechanism further allows
one to interpret the modulus $r_{c}$ as the vacuum expectation value (VEV) of
a dynamical radion field $\rho(x^{\mu})$, i.e.,
$r_{c}=\left\langle\rho\right\rangle$, with the effective potential
$V_{eff}(\rho)$ driving $\rho$ to settle at its minimum [19]. The mass of the
radion and its coupling to matter fields on the visible brane are suggested to
be $\mathcal{O}(\textrm{TeV})$, rendering it detectable, in principle, at
present generation colliders. Notably, the radion is expected to be lighter
than the lowest-lying KK excitations of generic bulk fields, which are
typically above $\mathcal{O}(\textrm{TeV})$. This should make the radion the
lightest detectable signature of the warped higher dimensional world.

There exists a wealth of work on the dynamics of various bulk fields and KK
gravitons in warped spacetime, alongside a variety of generalizations with
their own phenomenological and cosmological features (for a small body of
examples, see [20, 21, 22, 23, 24, 25, 26, 27, 28, 29, 30, 31, 32, 33, 34, 35,
36]). The original Goldberger-Wise analysis has also been extended along
several avenues, eg. by assuming finite values of the brane coupling constants
[37], incorporating back reaction of the bulk scalar on the metric [38, 39,
40], and dynamical stabilization schemes in cosmological backdrops [41, 42,
43]. However, in spite of its theoretical appeal, the original 5D RS model has
been increasingly challenged by experimental data in recent years. The mass of
the first excitation of the graviton KK tower in the 5D setting [23] is
suggested to be $\mathcal{O}(\textrm{TeV})$. Moreover, the graviton KK modes
are expected to couple to visible brane matter with amplified strength, as the
lowest-lying excitations cluster close to the TeV brane. These considerations
subsequently inspired several unsuccessful attempts to detect the signatures
of the graviton modes through various channels at the LHC [44, 45, 46, 47, 48,
49]. The absence of any such evidence till date has imposed serious
constraints on the parameter space of the model. In particular, a small
hierarchy $m_{H}/M_{5}\sim 10^{-2}$ is required to explain the null results,
indicating the potential appearance of new physics at least two orders of
magnitude below the fundamental Planck scale. This is particularly problematic
since this hierarchy requires $r_{c}^{-1}$ (which plays the role of the cut-
off above which new physics is expected to appear) to be lowered by two orders
of magnitude as well. But the sensitive dependence of the exponential warp
factor on $r_{c}$ severely restricts this possibility, thereby limiting the
efficacy of the five-dimensional setup in addressing the gauge hierarchy
problem.

##  2 Review of doubly warped braneworld model

Several higher dimensional extensions of the original RS model have been
proposed [50, 51, 52, 53], with most of them introducing additional orbifolds
with topology $S^{1}/\mathbb{Z}_{2}$ besides the first one. Some of these
models, generalized to six dimensions, have interesting cosmological features
stemming from dynamical stabilization of the extra spacelike dimensions [54].
A particularly interesting extension [55], capable of addressing the
aforementioned issues, emerges in the form of a doubly warped braneworld with
topology $\left[\mathcal{M}(1,3)\times S^{1}/\mathbb{Z}_{2}\right]\times
S^{1}/\mathbb{Z}_{2}$. It results in a “brane-box” configuration, with four
4-branes forming the “walls” of the box and four 3-branes located at the
“vertices” where adjacent 4-branes intersect. The bulk is a slice of
$\textrm{AdS}_{6}$ with cosmological constant $\Lambda_{6}\sim-\tilde{M}^{6}$,
where $\tilde{M}$ is the fundamental Planck scale. The complete bulk-brane
action $(S)$, comprised of the bulk Einstein-Hilbert action $(S_{6})$ and the
brane tension terms $(S_{5})$, is

$S=S_{6}+S_{5}$ $S_{6}=\int
d^{4}xdydz\sqrt{-g_{6}}\left(\dfrac{1}{2}\tilde{M}^{4}R_{6}-\Lambda_{6}\right)$
$\begin{split}S_{5}=&\int
d^{4}xdydz\sqrt{-g_{5}}\left[V_{1}(z)\delta(y)+V_{2}(z)\delta(y-\pi)\right]\\\
&+\int
d^{4}xdydz\sqrt{-\bar{g}_{5}}\left[V_{3}(y)\delta(z)+V_{4}(y)\delta(z-\pi)\right]\end{split}$
(1)

where $R_{6}$ is the six dimensional Ricci scalar, and $y$ and $z$ are angular
coordinates charting the extra dimensions. For the current purpose, the
visible 3-brane has been assumed to be devoid of any matter field, which would
have otherwise introduced a further ($S_{4}$) term as the contribution of the
matter Lagrangian. In a multiply warped setting, the brane tensions can, in
general, be functions of the bulk coordinates. This is a departure from the 5D
case, where the presence of only one extra dimension rules out any such
dependence at the very outset. Using an RS-like ansatz to solve the resulting
Einstein equations, one obtains the doubly warped metric

$ds^{2}=\dfrac{\textrm{cosh}^{2}(kz)}{\textrm{cosh}^{2}(k\pi)}\left(e^{-2c|y|}\eta_{\mu\nu}dx^{\mu}dx^{\nu}+R_{y}^{2}dy^{2}\right)+r_{z}^{2}dz^{2}$
(2)

where $R_{y}$ and $r_{z}$ are the radii of the orbifolds, and $c$ and $k$ are
constants interconnected through

$c=\dfrac{R_{y}k}{r_{z}\textrm{cosh}(k\pi)}\qquad;\qquad
k=r_{z}\sqrt{-\dfrac{\Lambda_{6}}{10\tilde{M}^{4}}}$ (3)

This peculiar relation makes it clear that in absence of a considerably large
hierarchy between $R_{y}$ and $r_{z}$, both $c$ and $k$ cannot be
simultaneously large. Instead, any one of the two distinct regimes $c>k$ or
$c<<k$ can emerge. It is not entirely possible to get rid of a little
hierarchy between $R_{y}$ and $r_{z}$, but with suitable choice of $c$ and
$k$, it can be ensured that $R_{y}/r_{z}$ does not exceed $\mathcal{O}(10)$.
The physical mass scale $m_{H}^{(ph)}$ on the 3-brane located at the vertex
$(y_{i},z_{j})$ is related to the fundamental scale $m_{H}$ by

$m_{H}^{(ph)}=m_{H}e^{-cy_{i}}\dfrac{\textrm{cosh}(kz_{j})}{\textrm{cosh}(k\pi)}$
(4)

In the $c>k$ regime, the choice of $c\sim 10$ and $k\sim 0.1$ yields warp
factors which together generate the TeV scale on the two 3-branes at
$\\{y=\pi,z=0\\}$ and $\\{y=\pi,z=\pi\\}$, one of which can be identified with
the visible brane. The other couple of 3-branes experiences negligible warping
and remains at the Planck scale. This clustering of multiple 3-branes around
each of the scales is a salient property of generalized RS models with nested
warping, one having important phenomenological implications.

Once the metric is fully solved, the brane tensions can be derived using
appropriate junction conditions at the orbifold fixed points; $V_{1}$ and
$V_{2}$ retain their $z$-dependence, whereas $V_{3}$ and $V_{4}$ turn out to
be independent of $y$.

$V_{1}(z)=-V_{2}(z)=8\tilde{M}^{2}\sqrt{-\dfrac{\Lambda_{6}}{10}}\textrm{sech}(kz)$
(5) $V_{3}(y)=0\qquad;\qquad
V_{4}(y)=-8\tilde{M}^{2}\sqrt{-\dfrac{\Lambda_{6}}{10}}\textrm{sech}(k\pi)$
(6)

The coordinate dependence of the former pair can be shown to emerge from
suitable scalar fields confined to the corresponding 4-branes. In the large
$k$ regime, one needs to account for a phantom scalar field on the $y=\pi$
brane. This field may serve as a potential candidate for dark energy on the
TeV branes. But the presence of such a phantom scalar renders the setup an
effective field theory as opposed to a fundamental one. No such phantom field
is needed in the large $c$ regime.

The doubly warped model has numerous phenomenological advantages over its
progenitor model. Firstly, owing to the presence of two extra dimensions, the
first excited graviton KK mode turns out to be considerably heavier than that
of the 5D case. Moreover, coupling between graviton KK modes and SM fields on
the visible brane is largely suppressed compared to the 5D model. Taken
together, these features can satisfactorily explain the non-detection of KK
gravitons at the LHC so far, without recourse to any small hierarchy between
$m_{H}$ and $\tilde{M}$ [56]. At the same time, significant portions of the
parameter space of the extended model remain accessible to the LHC, allowing
them to be explored in future runs [57]. Secondly, the doubly warped model can
offer an explanation of the mass hierarchy among the SM fermions [55]. In the
$c>k$ regime, the $\mathcal{O}$(TeV) 4-brane at $y=\pi$ intersects two other
4-branes at $z=0$ and $z=\pi$. Assuming SM fermions to be described by five
dimensional fields confined to the $y=\pi$ brane implies natural
$\mathcal{O}$(TeV) masses for the fermions. In addition to the $z$-dependent
bulk wavefunction on this 4-brane, the fermionic fields can have kinetic terms
on the 3-branes at the two intersection points. These boundary terms can
modify the fermion-scalar Yukawa coupling on the $\\{y=\pi,z=0\\}$ and
$\\{y=\pi,z=\pi\\}$ 3-branes, thereby causing a splitting among the effective
fermion masses. This splitting is arguably small, as the natural mass scales
of these 3-branes are already clustered close to each other around the TeV
scale.

##  3 Modulus stabilization in doubly warped model

Like the 5D RS model, the action of the 6D scenario contains no dynamics which
can stabilize the extra dimensional moduli around their desired values. In
absence of such an underlying stabilizing mechanism, the braneworld model
alone cannot be considered adequate. Motivated by the success of the 5D
Goldberger-Wise mechanism, it becomes natural to seek a similar approach for
the 6D model that might stabilize both $c$ and $k$ (or equivalently, $R_{y}$
and $r_{z}$) simultaneously. While other phenomenological aspects of multiply
warped spacetimes are well-studied [58, 59, 60, 61], there has been little
work in this direction so far, with the notable exception of [62]. The latter
study proposes disjoint stabilization mechanisms for $c$ and $k$ with the help
of two separate bulk and brane-localized fields in the $c>k$ regime. In the
other regime, taking $c<<k$ makes the metric almost conformally flat,
wherefore only $r_{z}$ needs to be stabilized satisfactorily, with $R_{y}$
either left unstabilized (which is justified because of the negligibly small
value of $c$) or stabilized with the help of another brane-localized field.
Whether stabilization of both moduli can be achieved with the help of a single
bulk scalar field, has, however, been an open question so far. In this paper,
we attempt to address this very question.

###  3.1 Dynamics of the bulk field

For our current purpose, we choose to work in the $c>k$ regime, which, as
noted already, obviates the need to introduce a phantom scalar (which raises
phenomenological difficulties) on the $y=\pi$ brane to explain the coordinate
dependence of its tension. Additionally, the relative smallness of $k$ allows
us to view the entire setup as a not-too-large departure from the 5D model,
helping identify the key points of deviation from the latter clearly. A study
of the complementary $c<<k$ domain using tools from supersymmetric quantum
mechanics appears in [62].

Analogous to the original Goldberger-Wise scenario, we consider a bulk field
propagating freely through the extra dimensional $y-z$ bulk, and interacting
only at the locations of the four 3-branes through quartic self-interaction
terms. The action $(S_{GW})$ of the bulk field is given by an immediate
generalization of the original 5D Goldberger-Wise action.

$S_{GW}=-\dfrac{1}{2}\int
d^{4}x\int\limits_{-\pi}^{+\pi}dy\int\limits_{-\pi}^{+\pi}dz\sqrt{-g_{6}}\left(g_{6}^{AB}\partial_{A}\phi\partial_{B}\phi+m^{2}\phi^{2}\right)-\sum_{i=1}^{4}\int
d^{4}x\sqrt{\tilde{g}_{4}^{(i)}}\lambda_{i}\left(\phi^{2}-u_{i}^{2}\right)^{2}$
(7)

Here, $\tilde{g}_{4}^{(i)}$ is the induced metric on the $i^{th}$ 3-brane,
with $u_{i}$ being the VEV of the bulk field (with mass dimension
$[u_{i}]=+2$) and $\lambda_{i}$ the corresponding coupling constant
($[\lambda_{i}]=-4$). With its dynamics confined to the higher dimensional
space, the bulk field is essentially “frozen” on the corner branes. Defining
$a(y)=\textrm{exp}(-c|y|)$ and $b(z)=\textrm{sech}(k\pi)\textrm{cosh}(kz)$, we
extremize this action with respect to $\phi$ to obtain the following equation
of motion.

$\begin{split}&-\dfrac{1}{R_{y}^{2}}\partial_{y}\left(a^{4}b^{3}\partial_{y}\phi\right)-\dfrac{1}{r_{z}^{2}}\partial_{z}\left(a^{4}b^{5}\partial_{z}\phi\right)+m^{2}a^{4}b^{5}\phi\\\
&+\dfrac{4a^{4}b^{4}}{R_{y}r_{z}}\left[\lambda_{1}\phi\left(\phi^{2}-u_{1}^{2}\right)\delta(y)\delta(z)+\lambda_{2}\phi\left(\phi^{2}-u_{2}^{2}\right)\delta(y-\pi)\delta(z)\right]\\\
&+\dfrac{4a^{4}b^{4}}{R_{y}r_{z}}\left[\lambda_{3}\phi\left(\phi^{2}-u_{3}^{2}\right)\delta(y)\delta(z-\pi)+\lambda_{4}\phi\left(\phi^{2}-u_{4}^{2}\right)\delta(y-\pi)\delta(z-\pi)\right]=0\end{split}$
(8)

Away from the boundaries, the contributions of the interaction terms vanish.
Assuming a separable solution of the form $\phi(y,z)=\phi_{1}(y)\phi_{2}(z)$,
the bulk equation of motion can be reduced to the following pair of uncoupled
ODEs, where $-4\alpha^{2}$ is the separation constant.

$-\dfrac{1}{R_{y}^{2}a^{4}\phi_{1}}\dfrac{d}{dy}\left(a^{4}\dfrac{d\phi_{1}}{dy}\right)=\dfrac{1}{r_{z}^{2}b^{3}\phi_{2}}\dfrac{d}{dz}\left(b^{5}\dfrac{d\phi_{2}}{dz}\right)-m^{2}b^{2}=-4\alpha^{2}$
(9)

An interesting physical interpretation of $\alpha$ emerges immediately. From
the equation for $\phi_{1}$, it is apparent that $2\alpha$ plays the role of
the mass of the $\phi_{1}$ field. Also, as the $k\to 0$ (and simultaneously
$r_{z}\to 0$) limit implies $b(z)\to 1$ and mathematically reduces the 6D
metric to the familiar singly warped form, it must also enforce $2\alpha\to m$
in order to reduce this equation to the 5D Goldberger-Wise bulk equation of
motion. As for the equation for $\phi_{2}$, it is satisfied trivially in this
limit, as $\phi_{2}(z)\to 1$ leads merely to $m^{2}=4\alpha^{2}$. However,
from a physical standpoint, one cannot of course allow $k$ to be arbitrarily
small, as any effective field theory in the semiclassical approach strictly
remains valid only for $r_{z}\geq M^{-1}$, with quantum gravity effects
dominating for $r_{z}<M^{-1}$. In this sense, the $k\to 0$ limit is not a
physically tenable one, but comes with a lower cut-off regulated by
$\Lambda_{6}$. The bottom line is simply that for sufficiently small $k$, one
can expect $\alpha$ to be reasonably close to $m$, which is a fact that proves
useful in due course.

eq. 9 can be solved exactly for the component fields $\phi_{1}(y)$ and
$\phi_{2}(z)$, leading to the following general solution for the bulk field
under the condition of $\mathbb{Z}_{2}$ orbifold symmetry.

$\phi(y,z)=e^{2c|y|}\left(Ae^{\nu c|y|}+Be^{-\nu
c|y|}\right)\left[DP_{n}^{l}(\textrm{tanh}(k|z|))+EQ_{n}^{l}(\textrm{tanh}(k|z|))\right]\textrm{sech}^{\frac{5}{2}}(kz)$
(10)

where $A$, $B$, $D$ and $E$ are four arbitrary constants (arising on account
of each equation being a second order ODE), and $P_{n}^{l}$ and $Q_{n}^{l}$
are associated Legendre functions of the first and second kind respectively.
The quantities $\nu$, $n$ and $l$ are defined as

$\nu=2\sqrt{1+\dfrac{R_{y}^{2}\alpha^{2}}{c^{2}}}$ (11)
$n=2\sqrt{1+\dfrac{\alpha^{2}r_{z}^{2}\textrm{cosh}^{2}(k\pi)}{k^{2}}}-\dfrac{1}{2}$
(12) $l=\dfrac{1}{2}\sqrt{25+\dfrac{4m^{2}r_{z}^{2}}{k^{2}}}$ (13)

Assuming the magnitudes of the brane coupling constants to be large (i.e.
$\lambda_{i}\to\infty$) allows the identification of simple energetically
favourable configurations to serve as boundary conditions. The very structure
of the interaction terms, given by
$\lambda_{i}(\phi^{2}-u_{i}^{2})^{2}\delta(y-y_{0})\delta(z-z_{0})$,
necessitates $\phi(y_{0},z_{0})\to u_{i}$ as such a configuration on the
corner brane at $(y_{0},z_{0})$. This approach leads to the following four
equations, where the shorthand $\tau_{z}=\textrm{tanh}(k|z|)$ has been
introduced for convenience.

$\phi_{1}(0)\phi_{2}(0)=(A+B)\left[DP_{n}^{l}(0)+EQ_{n}^{l}(0)\right]\approx
u_{1}$ (14) $\phi_{1}(\pi)\phi_{2}(0)=e^{2c\pi}(Ae^{\nu c\pi}+Be^{-\nu
c\pi})\left[DP_{n}^{l}(0)+EQ_{n}^{l}(0)\right]\approx u_{2}$ (15)
$\phi_{1}(0)\phi_{2}(\pi)=(A+B)\left[DP_{n}^{l}(\tau_{\pi})+EQ_{n}^{l}(\tau_{\pi})\right]\textrm{sech}^{\frac{5}{2}}(k\pi)\approx
u_{3}$ (16) $\phi_{1}(\pi)\phi_{2}(\pi)=e^{2c\pi}(Ae^{\nu c\pi}+Be^{-\nu
c\pi})\left[DP_{n}^{l}(\tau_{\pi})+EQ_{n}^{l}(\tau_{\pi})\right]\textrm{sech}^{\frac{5}{2}}(k\pi)\approx
u_{4}$ (17)

eq. 14$-$eq. 17 are not all linearly independent, and as such, only allow
three of the constants to be solved in terms of the remaining fourth. This
poses no problem though, as the forms of the solutions still allow $\phi(y,z)$
to be determined uniquely. In the large $c$ regime, the normalized solutions
of the individual component fields are

$\phi_{1}(y)=\sqrt{u_{1}}e^{2c|y|}(c_{1}e^{\nu c|y|}+c_{2}e^{-\nu c|y|})$ (18)
$\phi_{2}(z)=\sqrt{u_{1}}\left[c_{3}P_{n}^{l}(\tau_{z})+c_{4}Q_{n}^{l}(\tau_{z})\right]\textrm{sech}^{\frac{5}{2}}(kz)$
(19)

where $c_{2}=1$, and the other three constants are given by

$c_{1}=e^{-(\nu+2)c\pi}\left[\dfrac{u_{2}}{u_{1}}-e^{-(\nu-2)c\pi}\right]$
(20)
$c_{3}=\dfrac{Q_{n}^{l}(\tau_{\pi})-\left(\frac{u_{4}}{u_{2}}\right)Q_{n}^{l}(0)\textrm{cosh}^{\frac{5}{2}}(k\pi)}{P_{n}^{l}(0)Q_{n}^{l}(\tau_{\pi})-Q_{n}^{l}(0)P_{n}^{l}(\tau_{\pi})}$
(21)
$c_{4}=\dfrac{\left(\frac{u_{4}}{u_{2}}\right)P_{n}^{l}(0)\textrm{cosh}^{\frac{5}{2}}(k\pi)-P_{n}^{l}(\tau_{\pi})}{P_{n}^{l}(0)Q_{n}^{l}(\tau_{\pi})-Q_{n}^{l}(0)P_{n}^{l}(\tau_{\pi})}$
(22)

It is interesting to note that, apart from the overall normalization factor of
$u_{1}$, the dynamics of $\phi(y,z)$ is controlled not by the absolute values
of the VEVs but only by their ratios. This feature is reminiscent of the 5D
mechanism. Further, eq. 18$-$eq. 22 dictate $u_{4}/u_{2}=u_{3}/u_{1}$. So the
two ratios appearing explicitly in the solution are sufficient to completely
specify all the six possible ratios among the four VEVs, thus ruling out any
ambiguity.

###  3.2 The effective potential

Having obtained the solution of the bulk field, one needs to substitute it
back in equation eq. 7 and integrate over $y$ and $z$ in order to obtain the
stabilizing potential $V_{eff}(c,k)$, or equivalently, $V_{eff}(R_{y},r_{z})$.
From the resultant effective 4D action $S_{eff}$, the definition
$S_{eff}=-\int d^{4}xV_{eff}(c,k)$ allows the potential to be read off
directly. In the large $\lambda_{i}$ limit, this amounts to evaluating only
the bulk contribution. First, let us define the following dimensionless
parameters.

$\mu_{1}=\dfrac{mr_{z}}{k}\quad;\quad\mu_{2}=\dfrac{\alpha r_{z}}{k}$ (23)

These quantities roughly estimate the ratios between the mass parameters of
the bulk field and the fundamental Planck scale. From the semiclassical
standpoint, the admissible range of each ratio should be $0<\mu_{i}<1$. In
terms of these ratios, the parameters from eq. 11$-$eq. 13 can be re-expressed
as $\nu=2\sqrt{1+\mu_{2}^{2}\textrm{cosh}^{2}(k\pi)}$, $n=\nu-0.5$, and
$l=5/2\sqrt{1+4\mu_{1}^{2}/25}$. Upon substituting the solution from eq. 10 in
the action, the integral over $y$ can be readily evaluated, but the presence
of the special functions makes the $z$-integral analytically intractable.
Defining the dimensionless potential $\tilde{V}_{eff}=V_{eff}/u_{1}^{2}$, and
making use of eq. 3 to eliminate $R_{y}$ and $r_{z}$ in favour of $c$ and $k$,
the following form emerges for the potential

$\begin{split}\tilde{V}_{eff}(c,k)=&\left[\dfrac{c_{1}^{2}(\nu+2)^{2}}{2\nu}e^{2\nu
c\pi}+\dfrac{c_{2}^{2}(\nu-2)^{2}}{2\nu}-2\pi
c_{1}c_{2}(\nu^{2}-4)c\right]\times F_{1}(k)\\\
&+\left(\dfrac{c_{1}^{2}}{2\nu}e^{2\nu c\pi}+\dfrac{c_{2}^{2}}{2\nu}+2\pi
c_{1}c_{2}c\right)\times F_{2}(k)\end{split}$ (24)

where the functions $F_{1}(k)$ and $F_{2}(k)$, arising out of the
$z$-integral, are given explicitly by

$F_{1}(k)=k\textrm{sech}^{4}(k\pi)\int\limits_{-\pi}^{+\pi}dz\left(c_{3}P_{n}^{l}(\tau_{z})+c_{4}Q_{n}^{l}(\tau_{z})\right)^{2}\textrm{sech}^{2}(kz)$
(25)
$\begin{split}F_{2}(k)=\dfrac{1}{4}k\textrm{sech}^{4}(k\pi)\int\limits_{-\pi}^{+\pi}dz&\left[2(l-n-1)(c_{3}P_{n+1}^{l}(\tau_{z})+c_{4}Q_{n+1}^{l}(\tau_{z}))+(2n-3)\tau_{z}(c_{3}P_{n}^{l}(\tau_{z})+c_{4}Q_{n}^{l}(\tau_{z}))\right]^{2}\\\
&+k\textrm{sech}^{4}(k\pi)\int\limits_{-\pi}^{+\pi}dz\>\mu_{1}^{2}(c_{3}P_{n}^{l}(\tau_{z})+c_{4}Q_{n}^{l}(\tau_{z}))^{2}\end{split}$
(26)

There is little choice but to evaluate $F_{1}(k)$ and $F_{2}(k)$ numerically
for different values of $\mu_{1}$, $\mu_{2}$, and $u_{4}/u_{2}$, which, by
this point, serve as three of the fundamental parameters of the model. The
fourth parameter $u_{2}/u_{1}$, contained in the coefficient $c_{1}$,
contributes chiefly to the $c$-dependence of the potential. Equipped with the
potential given by eq. 24$-$eq. 26, we are in a position to demonstrate the
existence of a simultaneous minimum of $\tilde{V}_{eff}$ in $c$ and $k$ over
some region of the parameter space that doesn’t require excessive fine tuning.

###  3.3 Stabilizing $k$

Owing to the explicit form of $c_{1}$ from eq. 20, both the linear and
quadratic terms in $c_{1}$ appearing in the coefficients of $F_{1}(k)$ and
$F_{2}(k)$ in eq. 24 are suppressed at least by $\mathcal{O}(e^{-4c\pi})$. For
$c\sim 10$, this suppression factor is nearly of order $10^{-68}$. So the
leading order terms within both coefficients are the $c_{2}^{2}$ terms which
suffer no such suppression. In order to study the $k$-dependence of the
potential, it suffices to approximate $\tilde{V}_{eff}$ by taking only the
dominant contributions offered by these two terms. Recalling $c_{2}=1$, the
potential effectively reduces to a single-variable function of $k$ given by

$\tilde{V}_{eff}(c,k)\to\dfrac{1}{2}F(k)\>\approx\>\dfrac{1}{4\sqrt{1+\mu_{2}^{2}\textrm{cosh}^{2}(k\pi)}}\left[4\left(\sqrt{1+\mu_{2}^{2}\textrm{cosh}^{2}(k\pi)}-1\right)^{2}F_{1}(k)+F_{2}(k)\right]$
(27)

Establishing the stabilization of $k$ is thus tantamount to locating a
suitable minimum of $F(k)$. The plots in figure 1 show the behaviour of $F(k)$
for various combinations of the parameters.

Figure 1: Profiles of $F(k)$ vs $k$ for VEV ratios $u_{4}/u_{2}=1.3$ (solid
blue), $u_{4}/u_{2}=1.5$ (dotted green), and $u_{4}/u_{2}=1.7$ (dashed red),
and a few representative values of $\mu_{1}$ and $\mu_{2}$ for each case.

For all the choices of parameters, $F(k)$ admits a broad minimum at some
$k<0.5$. The location of this minimum is physically important, as the range
$0.1\leq k\leq 0.6$ is of particular interest to scenarios attempting to
explain the non-detectability of graviton KK modes at the LHC [56][57], or
exploring the phenomenology of off-brane SM fields that extend into the bulk
[60][61]. As special cases, it can be checked that $\mu_{2}\to 0$ produces an
asymptotically decaying $F(k)$ with no finite minimum, whereas $u_{4}=u_{2}$
admits only $k=0$ as the global minimum. As expected, the most pronounced
dependence of the minimum is on the VEV ratio $u_{4}/u_{2}$, which governs the
potential at the leading order.

###  3.4 Stabilizing $c$

The previous plots reveal the existence of a fairly large parameter space
which can stabilize $k$ around its desired value. But $\tilde{V}_{eff}(c,k)$
also needs to stabilize $c$, which plays the leading role in determining the
degree of warping. To this end, we need to reinstate the previously neglected
$c$-dependence in $\tilde{V}_{eff}(c,k)$ in order to locate a suitable minimum
along $c$. Owing to the form of $c_{1}$ in eq. 20, one might expect such a
minimum to occur at $c_{1}=0$, which would be close in spirit to the 5D
Goldberger-Wise result. Upon closer inspection though, it becomes clear that
the existence of such a minimum depends crucially on the value of $\mu_{2}$.
For considerably small values of $\mu_{2}$ which render
$\nu\>\approx\>2+\epsilon_{k}$, where
$\epsilon_{k}=\mu_{2}^{2}\textrm{cosh}^{2}(k\pi)$ is small enough so that
$\mathcal{O}(\epsilon_{k}^{2})$ and higher are negligible compared to unity,
the deviation $\delta\tilde{V}_{eff}$ for any arbitrarily small deviation
$c=c_{0}-\delta$ from the alleged minimum $c_{0}$ can be estimated as

$\delta\tilde{V}_{eff}|_{c_{0}}\sim
e^{-4\pi(c_{0}-\delta)-2\epsilon_{k}c_{0}\pi}\left[25\left(e^{\epsilon_{k}\delta\pi}-1\right)^{2}-8\pi(c_{0}-\delta)e^{\epsilon_{k}\delta\pi}\left(e^{\epsilon_{k}\delta\pi}-1\right)\left(1-6\epsilon_{k}\right)\right]$
(28)

For excessively small $\epsilon_{k}$ and vanishingly small $\delta$, this
quantity is negative. This clearly rules out a minimum at $c_{0}$. The trouble
can be traced to the existence of the $\mathcal{O}(c_{1})$ terms in
$\tilde{V}_{eff}(c,k)$. While the 5D stabilization mechanism ensured automatic
cancellation of the pair of $\mathcal{O}(c_{1})$ terms in the potential, the
current model offers no such way out by default. However, if some particular
combination of parameters can result in such cancellation (without excessive
fine tuning of course), the issue would be resolved. Once $k$ has been
stabilized at some specific $k_{min}$, the condition for the
$\mathcal{O}(c_{1})$ terms to cancel each other translates to

$\left[\nu(k_{min})^{2}-4\right]F_{1}(k_{min})-F_{2}(k_{min})=0$ (29)

As $F_{1}(k)$ and $F_{2}(k)$ themselves depend on all three parameters in a
rather complicated manner, this equation is transcendental and cannot be
solved analytically. However, it can be checked numerically that for every
choice of $\mu_{1}$, the allowed parameter space contains approximate
solutions to eq. 29, some of which are presented in Table 1.

$\mu_{1}$ | $\mu_{2}$ | $u_{4}/u_{2}$ | $k_{min}$ | $\nu(k_{min})$ | $F_{1}(k_{min})$ | $F_{2}(k_{min})$ | $u_{1}/u_{2}$
---|---|---|---|---|---|---|---
0.40 | 0.24 | 1.122 | 0.119 | 2.065 | 0.697 | 0.185 | 20
0.50 | 0.30 | 1.160 | 0.123 | 2.102 | 0.741 | 0.307 | 40
0.60 | 0.36 | 1.200 | 0.124 | 2.145 | 0.777 | 0.469 | $\sim 10^{2}$
0.70 | 0.42 | 1.248 | 0.128 | 2.197 | 0.833 | 0.689 | $\sim 10^{3}$
0.80 | 0.48 | 1.303 | 0.132 | 2.256 | 0.899 | 0.980 | $\sim 10^{4}$
Table 1: A few examples of critical parameter values approximately satisfying
eq. 29 and allowing simultaneous stabilization of both moduli around
$k_{min}\sim 0.1$ and $c_{min}\sim 10$.

Each of the combinations ensures that the non-negative
$\mathcal{O}(c_{1}^{2})$ terms dominate over the $\mathcal{O}(c_{1})$ terms by
nearly two orders of magnitude, leading to a proper minimum $c_{min}$
(satisfying $c_{1}=0$) with the following approximate form

$c_{min}\>\approx\>\dfrac{1}{\pi\left[\nu(k_{min})-2\right]}\textrm{ln}\left(\dfrac{u_{1}}{u_{2}}\right)$
(30)

One can immediately estimate the corresponding $u_{1}/u_{2}$ ratio which
produces the desired value of $c_{min}\sim 11.54$. As given in Table 1, its
magnitude depends strongly on $\mu_{1}$. In order to stabilize $k$ between
$0.1$ and $0.6$, the minimum admissible magnitude of $u_{1}/u_{2}$ is
$\mathcal{O}(10)$, which is slightly larger than the requirement in 5D. For
larger values of $\mu_{1}$, a minor hierarchy between $u_{1}$ and $u_{2}$
becomes increasingly prominent. This hierarchy has its physical origin in the
fact that we have essentially made the component fields $\phi_{1}(y)$ and
$\phi_{2}(z)$ yield two distinct warping scales $k\sim 0.1$ and $c\sim 10$ in
spite of having chosen their mass parameters ($m$ and $\alpha$) to be of
comparable magnitudes! The somewhat large $u_{1}/u_{2}$ ratio is precisely the
price we need to pay for that. As a reality check, it is instructive to check
that choosing $u_{1}\sim u_{2}$ in eq. 30 would have resulted in $c_{min}\sim
0.1$, just as this argument suggests.

The mutual cancellation of the $\mathcal{O}(c_{1})$ terms involves a certain
degree of fine tuning among the three parameters which control $F_{1}(k)$ and
$F_{2}(k)$. While the situation is clearly more delicate than the 5D case,
this tuning need not be very extreme. The $\mathcal{O}(c_{1}^{2})$ terms can
dominate and provide $c$ with a stable minimum if the cancellation of the
$\mathcal{O}(c_{1})$ terms in eq. 29 is accurate roughly up to
$\mathcal{O}(10^{-2})$, which is sufficient to suppress the latters’
contribution by $\mathcal{O}(10^{-2})$. For $c\sim 10$, the critical values of
the parameters need to be accurate at most up to $\mathcal{O}(10^{-3})$, while
larger values of $c$ are somewhat more likely to ameliorate the situation due
to increased suppression. This is no worse than the fine tuning associated
with the choice of $m/M$ for a given VEV ratio, that is required to generate
the TeV scale accurately in the 5D case. Moreover, for some fixed $\nu(k)$,
the tuning associated with $u_{2}/u_{1}$ can be even more lenient due to the
logarithmic dependence of $c_{min}$ on $u_{2}/u_{1}$. So simultaneous
stabilization of $c$ and $k$ depends on a small degree of fine tuning among
the parameters $\mu_{1}$, $\mu_{2}$, and $u_{4}/u_{2}$, which constitutes a
crucial aspect of the extended Goldberger-Wise mechanism in the doubly warped
scenario. This feature is not surprising since one should physically expect a
two-level tuning for the stabilization of two distinct moduli in a spacetime
with nested warping. As the stabilization of $k$ alone requires negligible
tuning (as evident from Fig. 1), the stabilization of $c$ justifiably involves
both levels. It might be interesting to investigate if incorporating back
reaction or quartic self-interaction terms within the bulk can lead to
improvements for this tuning requirement.

###  3.5 Phenomenological implications

The minor hierarchy between $u_{1}$ and $u_{2}$ (or equivalently between
$u_{3}$ and $u_{4}$) can be better interpreted if we arrange the magnitudes of
the four VEVs in proper order. As the boundary conditions require
$u_{1}/u_{2}=u_{3}/u_{4}$, the VEVs obtained in each case satisfy

$u_{3}>u_{1}>u_{4}>u_{2}$ (31)

The first two values are very close to each other, as are the last two, with
the aforesaid hierarchy pushing the pairs apart. Interestingly, eq. 31
reflects the order of physical mass scales on the corresponding corner branes
(on which these classical values of $\phi$ are defined) for large $c$ and
small $k$, as can be checked using eq. 4 . This can be understood physically
as the massive bulk scalar, once frozen on the boundary branes, naturally
tends to act against the warping induced by the bulk energy density.
Consequently, the resulting mass scales on corner branes associated with
larger VEVs can be generally expected to be higher than those with smaller
VEVs. This feature is visible in the 5D Goldberger-Wise mechanism as well,
albeit in a greatly tempered form. In the present case, it is more pronounced
as it also incorporates the clustering effect observed among the physical mass
scales.

Due to the smallness of $k$, phenomenological features associated with $c$
should be sufficiently similar to those of the 5D model. At first glance, this
might make the smallness of $u_{2}$ compared to $u_{1}$ for larger values of
$\mu_{1}$ appear alarming, as the mass of the associated radion in the 5D
model is proportional to the ratio between the classical value of $\phi$ on
the TeV brane and the fundamental scale. In particular, one might worry that
such small $u_{2}$ could bring the predicted radion mass down to the MeV scale
or lower, thereby giving rise to phenomenological issues. It must be borne in
mind, however, that we have neglected the back reaction of the bulk field. On
occasions, this can cause underestimation of the radion mass by a few orders
of magnitude [38][39]. Furthermore, unlike the 5D model, one can exploit the
6D setup immensely by allowing SM fields to extend into the bulk. This
possibility, which constitutes one of the most interesting features of the
doubly warped model, might imply further corrections to the radion mass due to
interactions of $\phi$ with these fields. There can also be significant mixing
between $\phi$ and the scalar degrees of freedom which give rise to coordinate
dependent 4-brane tensions. Any realistic study of radion phenomenology needs
to address these questions in order to estimate the magnitude of the radion
mass accurately.

Besides, due to the presence of two distinct moduli in the model, one should
physically expect the appearance of two fundamental radions with different
masses and widths. Owing to the smallness of $k$, the second radion should not
suffer any significant suppression from the fundamental scale. This is echoed
by the form of the stabilizing potential $\tilde{V}_{eff}(c,k)$ from eq. 24,
where $F(k)$ is usually $\mathcal{O}(1)$ around its minimum. However, the
aforementioned subtleties need to be accounted for in this case too before its
phenomenology can be deduced conclusively. Due to the evidently complicated
nature of such a study, we defer it to a future work.

##  4 Insight from higher curvature gravity

One of the most compelling advantages of having the size of the extra
dimensions set by a single bulk field (as opposed to separate bulk and brane
localized fields) is that it allows a purely gravitational interpretation of
the stabilizing mechanism. Conventional wisdom suggests that the Einstein-
Hilbert action, which provides an effective low energy description of gravity,
needs to be amended with additional higher curvature terms respecting
diffeomorphism invariance at sufficiently high energy scales. The warped
geometry model which is being considered here has a large cosmological
constant ( $\sim M_{P}$ ) in the bulk and as a result the inclusion of higher
curvature terms is a natural choice. Two broad classes of such higher
curvature theories are the quasi-linear Lanczos-Lovelock models and $f(R)$
models. While Lanczos-Lovelock models enjoy the benefit of being naturally
ghost-free [63, 64, 65], the mathematically simpler $f(R)$ models, equipped
with specific conditions to ensure freedom from ghosts, pass some of the
cosmological tests [66, 67, 68, 69, 70, 71]. Furthermore, any given $f(R)$
action typically admits a dual scalar-tensor representation [72, 73, 74, 75].
In the so-called Einstein frame (related to the Jordan frame through a
conformal transformation), the situation is equivalent to that of a scalar
field $\tilde{\phi}$ coupled minimally to gravity, alongside a potential
$U(\tilde{\phi})$ whose form is determined by the functional form of $f(R)$.
For singularity-free metrics which are not experiencing rapid evolution, this
equivalence holds physically [76, 77, 78, 79, 80].

In recent works [81][82], it has been shown how such a scalar degree of
freedom, arising solely from gravity in the 5D RS model, can play the role of
the bulk field in the Goldberger-Wise scheme, thus obviating the need to
introduce the latter by hand. By choosing $f(R)=R+aR^{2}-|b|R^{4}$ (where $a$
and $b$ are coupling constants satisfying $a>0$ and $a>|b|$ to ensure freedom
from ghosts), the potential can be arranged to contain both quadratic and
quartic terms, with the latter encapsulating the effects of back reaction.
Although the technique can be readily extended to spacetimes of arbitrary
dimensionality, the inclusion of back reaction quickly renders multiply warped
settings intractable in their full generality. In the following analysis, we
discuss an extension of this technique to the doubly warped model, by
retaining only the lowest-order correction term in $f(R)$. To that end, we
choose $f(R)=R+aR^{2}$, which provides the action

$\mathcal{A}=\int
d^{6}x\sqrt{-g}\left[\dfrac{f(R)}{2\kappa_{6}^{2}}-\Lambda_{6}\right]=\int
d^{6}x\sqrt{-g}\left[\dfrac{1}{2\kappa_{6}^{2}}\left(R+aR^{2}\right)-\Lambda_{6}\right]$
(32)

where $\kappa_{6}$ is the six-dimensional gravitational constant. In order to
arrive at the Einstein frame, one conventionally makes a detour [75] through
the intermediate Jordan frame representation

$\mathcal{A}=\int d^{6}x\sqrt{-g}\left[\dfrac{1}{2\kappa_{6}^{2}}\\{\psi
R-V(\psi)\\}-\Lambda_{6}\right]\quad,\quad
V(\psi)=\chi(\psi)\psi-f(\chi(\psi))$ (33)

where one introduces the auxiliary field $\chi$ and defines
$\psi=f^{\prime}(\chi)$. For $f^{\prime\prime}(\chi)\neq 0$, the equation of
motion for $\chi$ from eq. 33, i.e. the on-shell condition, imposes $\chi=R$.
This makes eq. 33 equivalent to the original $f(R)$ action in eq. 32. In order
to reduce this to the minimally coupled Einstein frame representation, we
apply the conformal transformation

$g_{AB}\to\tilde{g}_{AB}=\sqrt{\psi}g_{AB}\quad,\quad\kappa_{6}\tilde{\phi}=\dfrac{\sqrt{5}}{2}\>\textrm{ln}\>\psi=\dfrac{\sqrt{5}}{2}\>\textrm{ln}\>f^{\prime}(R)$
(34)

where the last equality follows from the on-shell condition. The action then
transforms to

$\mathcal{A}=\int
d^{6}x\sqrt{-\tilde{g}}\left[\left(\dfrac{\tilde{R}}{2\kappa_{6}^{2}}-\Lambda_{6}\right)-\dfrac{1}{2}\partial_{A}\tilde{\phi}\partial^{A}\tilde{\phi}-U(\tilde{\phi})\right]\quad,\quad
U(\tilde{\phi})=\dfrac{Rf^{\prime}(R)-f(R)}{2\kappa_{6}^{2}(f^{\prime}(R))^{\frac{3}{2}}}$
(35)

Substituting $f(R)=R+aR^{2}$ in eq. 34 and inverting the relation yields
$R(\tilde{\phi})$, which can be plugged immediately in eq. 35 to obtain
$U(\tilde{\phi})$ as follows.

$R=\dfrac{1}{2a}\left(e^{\frac{2}{\sqrt{5}}\kappa_{6}\tilde{\phi}}-1\right)\quad\implies\quad
U(\tilde{\phi})=\dfrac{1}{8a\kappa_{6}^{2}}\left(e^{\frac{1}{\sqrt{5}}\kappa_{6}\tilde{\phi}}-2e^{-\frac{1}{\sqrt{5}}\kappa_{6}\tilde{\phi}}+e^{-\frac{3}{\sqrt{5}}\kappa_{6}\tilde{\phi}}\right)$
(36)

The minimum of $U(\tilde{\phi})$ occurs at $\tilde{\phi}=0$, as evident from
$U^{\prime}(0)=0$ and $U^{\prime\prime}(0)>0$. Expanding $U(\tilde{\phi})$
about this minimum, the leading order non-vanishing contribution comes from
the quadratic term, with all subsequent terms increasingly suppressed by
higher powers of $\kappa_{6}$.

$U(\tilde{\phi})\>\approx\>\dfrac{1}{8a\kappa_{6}^{2}}\left[\left(1+\dfrac{\kappa_{6}\tilde{\phi}}{\sqrt{5}}+\dfrac{\kappa_{6}^{2}\tilde{\phi}^{2}}{10}\right)-2\left(1-\dfrac{\kappa_{6}\tilde{\phi}}{\sqrt{5}}+\dfrac{\kappa_{6}^{2}\tilde{\phi}^{2}}{10}\right)+\left(1-\dfrac{3\kappa_{6}\tilde{\phi}}{\sqrt{5}}+\dfrac{9\kappa_{6}^{2}\tilde{\phi}^{2}}{10}\right)\right]=\dfrac{\tilde{\phi}^{2}}{10a}$
(37)

Having identified $1/5a$ with the mass squared of the scalar mode
$\tilde{\phi}$, the action from eq. 35 reduces precisely to the sum of the
Einstein-Hilbert action and the Goldberger-Wise action in the 6D bulk, with
$\tilde{\phi}$ in the Einstein frame remarkably playing the role of the bulk
field.

While [82] uses this approach to deal with the full back-reacted problem (by
including the $R^{4}$ term in the $f(R)$ Lagrangian) for the 5D RS model in
the Jordan and Einstein frames separately, such an exact treatment is too
complicated to be feasible in the 6D case. Instead, with the results from the
gravitational sector at hand, one can use heuristic arguments to bridge the
gap with the existing bulk field method from the preceding sections. As the
physical origin of the higher curvature correction(s) can be traced to the
bulk energy density, it stands to reason that the resulting scalar
$\tilde{\phi}$ should be an explicit function of the coordinates $y$ and $z$
only. With the metric solution in the Einstein frame given by eq. 2, the
equation of motion for $\tilde{\phi}$ is identical to eq. 8 away from the
boundaries. This prevents $\tilde{\phi}$ from having any non-trivial dynamics
on the 3-branes, which is in agreement with the bulk field prescription. As an
immediate corollary, we obtain four constant values of $\tilde{\phi}$ (having
mass dimension $+2$) serving as fixed boundary values on the four branes. The
results are completely analogous to eq. 14$-$eq. 17, with the only difference
being that the current approach renders these boundary conditions exact,
whereas in the earlier formulation they were valid only in the large coupling
regime. With the full solution eq. 18$-$eq. 22 for $\tilde{\phi}$ subsequently
at hand, the rest of the analysis can proceed exactly as before, leading to
the stabilization of both moduli under appropriate choices of parameters.

##  5 Discussions

The prospect of stabilizing both moduli of a doubly warped Randall-Sundrum
braneworld model using a single bulk scalar field has been studied. Such a
mechanism is crucial for a complete resolution of the gauge hierarchy problem
in a higher dimensional scenario, and needs to supplement the gravitational
part of the action giving rise to the warped metric. While the approach taken
here is essentially a direct generalization of the Goldberger-Wise mechanism
to six dimensions, the presence of nested warping brings out additional
subtleties and constraints on the parameter space. As demonstrated, these
constraints do not necessarily involve any extreme fine tuning of the
fundamental parameters. In the $c>k$ regime, which is phenomenologically
preferred as it requires no brane-localized phantom field to explain the
coordinate dependence of the 4-brane tensions, the effective potential admits
a true minimum in $k$ around $k_{min}\sim 0.1$ without any significant tuning.
The stabilization of $c$, on the other hand, requires tuning on two different
levels: firstly among the parameters $\mu_{1}$, $\mu_{2}$, and $u_{4}/u_{2}$
(for the appearance of a suitable $c_{min}$), and secondly for $u_{1}/u_{2}$
(to ensure $c_{min}\sim 12$). This can be interpreted physically as a
consequence of the doubly warped structure of the underlying spacetime. A
further departure from the singly warped model is that in order to achieve
$0.1\leq k_{min}\leq 0.6$ and $c_{min}\sim 12$, the minimum admissible
magnitude of $u_{1}/u_{2}$ is $\mathcal{O}(10)$, which is one order of
magnitude higher than the analogous requirement in case of the 5D mechanism.

The bulk scalar approach is especially attractive as it allows room for a
purely gravitational interpretation, with higher curvature contributions in
the bulk automatically giving rise to the required scalar mode and its
potential in the Einstein frame. Since the bulk of such warped geometry model
is endowed with large bulk cosmological constant , the contributions from
higher curvature terms become natural. This motivates us to include higher
curvature terms in the bulk such as in $f(R)$ model. This feature
distinguishes it from certain other stabilization schemes, eg. stabilization
of the two moduli with two distinct bulk and brane-localized fields. For
appropriate choices of $f(R)$, a variety of bulk scalar potentials can be
generated, the simplest of which is the quadratic potential considered here.
While the singly warped model could be solved exactly in presence of non-
negligible back reaction, the doubly warped model becomes unsolvable as the
field equations turn out to be non-linear PDEs. Nevertheless, the well-known
$R^{2}$ correction is sufficient to make contact with the conventional
approach on physical grounds. In soothe, at sufficiently high energy scales,
gravity alone appears capable of both warping spacetime and determining the
degree of warping. This possibility is intriguing, as it may produce
observable TeV scale signatures of higher curvature gravity that can be
explored in future collider experiments. In principle, foremost among them
should be the $\mathcal{O}(\textrm{TeV})$ radion mass associated with the
larger modulus, which would be of fundamental importance in
estimating/constraining the magnitude(s) of the higher curvature coupling(s).

The current study can be extended along various avenues. As pointed out
earlier, a realistic study of radion phenomenology in a multiply warped
background needs to include interactions of the bulk field with higher
dimensional fermionic and gauge fields into account, alongside mixing with the
brane-localized scalars responsible for coordinate dependent brane tensions.
These effects, which constitute salient features of geometries with nested
warping, may introduce a plethora of non-trivial modifications as far as
collider signatures are concerned. Theoretically, it is worth investigating
how other viable choices of $f(R)$, or other classes of higher curvature
theories (eg. Einstein-Gauss-Bonnet gravity), incorporate higher order
phenomena like significant back reaction and affect the stabilization scheme.
Such studies would necessarily have to rely on numerical techniques due to the
complexity of the model. As plausible alternatives, one can also attempt to
explain the origin of the stabilizing potential from quantum and/or
thermodynamic perspectives, proceeding along the lines of [83, 84, 85, 86,
87]. Finally, it would be interesting to study the role of multiple dynamical
radion fields in various cosmological contexts as well, eg. in the
inflationary and bouncing settings, alongside their signature on the Cosmic
Microwave Background.

## Acknowledgements

Research work of AB was partly supported by MS studentship of IACS Kolkata.

## References

  * [1] I. Antoniadis, “A possible new dimension at a few TeV”, Phys. Lett. B 246, 377 (1990).
  * [2] P. Horava and E. Witten, “Heterotic and type I string dynamics from eleven-dimensions”, Nucl. Phys. B, 460, 506 (1996).
  * [3] P. Horava and E. Witten, “Eleven-dimensional supergravity on a manifold with boundary”, Nucl. Phys. B 475, 94 (1996).
  * [4] N. Arkani-Hamed, S. Dimopoulos, and G. Dvali, “The hierarchy problem and new dimensions at a millimeter”, Phys. Lett. B 429, 263 (1998).
  * [5] I. Antoniadis, N. Arkani-Hamed, S. Dimopoulos, and G. Dvali, “New dimensions at a millimeter to a Fermi and superstrings at a TeV”, Phys. Lett. B 436, 257 (1998).
  * [6] N. Arkani-Hamed, S. Dimopoulos, and G. Dvali, “Phenomenology, astrophysics, and cosmology of theories with submillimeter dimensions and TeV scale quantum gravity”, Phys. Rev. D 59, 086004 (1999).
  * [7] Z. Kakushadze and S.H.H. Tye, “Brane world”, Nucl. Phys. B 548, 180 (1999).
  * [8] N. Kaloper, “Bent domain walls as braneworlds”, Phys. Rev. D 60, 123506 (1999).
  * [9] L. Randall and R. Sundrum, “Large Mass Hierarchy from a Small Extra Dimension”, Phys. Rev. Lett. 83, 3370 (1999).
  * [10] L. Randall and R. Sundrum, “An Alternative to Compactification”, Phys. Rev. Lett. 83, 4690 (1999).
  * [11] R. Sundrum, “Effective field theory for a three-brane universe”, Phys. Rev. D 59, 085009 (1999).
  * [12] T. Appelquist, H.C. Cheng, and B.A. Dobrescu, “Bounds on universal extra dimensions”, Phys. Rev. D 64, 035002 (2001).
  * [13] T.G. Rizzo, “Probes of universal extra dimensions at colliders”, Phys. Rev. D 64, 095010 (2001).
  * [14] M. Gogberashvili, “Hierarchy problem in the shell-universe model”, Int. J. Mod. Phys. D 11, 1635 (2002).
  * [15] G. Aad et al. [ATLAS Collaboration], “Observation of a new particle in the search for the Standard Model Higgs boson with the ATLAS detector at the LHC”, Phys. Lett. B 716, 1 (2012).
  * [16] G. Aad et al. [ATLAS Collaboration], “Combined search for the Standard Model Higgs boson in pp collisions at $\sqrt{s}=7\>\textrm{TeV}$ with the ATLAS detector”, Phys. Rev. D 86, 032003 (2012).
  * [17] S. Chatrchyan et al. [CMS Collaboration], “Observation of a new boson at a mass of 125 GeV with the CMS experiment at the LHC”, Phys. Lett. B 716, 30 (2012).
  * [18] W.D. Goldberger and M.B. Wise, “Modulus Stabilization with Bulk Fields”, Phys. Rev. Lett. 83, 4922 (1999).
  * [19] W.D. Goldberger and M.B. Wise, “Phenomenology of a stabilized modulus”, Phys. Lett. B 475, 275 (1999).
  * [20] H. Davoudiasl, J. L. Hewett, and T. G. Rizzo, “Bulk gauge fields in the Randall–Sundrum model”, Phys. Lett. B 473, 43 (1999).
  * [21] O. DeWolfe, D.Z. Freedman, S.S. Gubser, and A. Karch, “Modeling the fifth dimension with scalars and gravity”, Phys. Rev. D 62, 046008 (2000).
  * [22] S. Chang, J. Hisano, H. Nakano, N. Okada, and M. Yamaguchi, “Bulk standard model in the Randall-Sundrum background”, Phys. Rev. D 62, 084025 (2000).
  * [23] H. Davoudiasl, J. L. Hewett, and T. G. Rizzo, “Phenomenology of the Randall-Sundrum Gauge Hierarchy Model”, Phys. Rev. Lett. 84, 2080 (2000).
  * [24] S.J. Huber, “Fermion masses, mixings and proton decay in a Randall–Sundrum model”, Phys. Lett. B 498, 256 (2001).
  * [25] S.J. Huber and Q. Shafi, “Higgs mechanism and bulk gauge boson masses in the Randall-Sundrum model”, Phys. Rev. D 63, 045010 (2001).
  * [26] H. Davoudiasl, J. L. Hewett, and T. G. Rizzo, “Experimental probes of localized gravity: On and off the wall”, Phys. Rev. D 63, 075004 (2001).
  * [27] H. Davoudiasl, J. L. Hewett, and T. G. Rizzo, “Brane-localized kinetic terms in the Randall-Sundrum model”, Phys. Rev. D 68, 045002 (2003).
  * [28] S. Das, D. Maity, and S. SenGupta, “Cosmological constant, brane tension and large hierarchy in a generalized Randall-Sundrum braneworld scenario”, J. High Energy Phys. 05, 042 (2008).
  * [29] R. Koley, J. Mitra, and S. SenGupta, “Fermion localization in a generalized Randall-Sundrum model”, Phys. Rev. D 79, 041902(R) (2009).
  * [30] C. Csáki, M. Graesser, C. Kolda, and J. Terning, “Cosmology of one extra dimension with localized gravity”, Phys. Lett. B 462, 34 (1999).
  * [31] T. Nihei, “Inflation in the five-dimensional universe with an orbifold extra dimension”, Phys. Lett. B 465, 81 (1999).
  * [32] C. Csáki, M. Graesser, L. Randall, and J. Terning, “Cosmology of brane models with radion stabilization”, Phys. Rev. D 62, 045015 (2000).
  * [33] S. Pal, S. Bharadwaj, and S. Kar, “Can extra-dimensional effects replace dark matter?”, Phys. Lett. B 609, 194 (2005).
  * [34] F. Chen, J.M. Cline, and S. Kanno, “Modified Friedmann equation and inflation in a warped codimension-two braneworld”, Phys. Rev. D 77, 063531 (2008).
  * [35] P. Dey, B. Mukhopadhyaya, and S. SenGupta, “Neutrino masses, the cosmological constant, and a stable universe in a Randall-Sundrum scenario”, Phys. Rev. D 80, 055029 (2009).
  * [36] N. Banerjee and T. Paul, “Inflationary scenario from higher curvature warped spacetime”, Eur. Phys. J. C 77, 672 (2017).
  * [37] A. Dey, D. Maity, and S. SenGupta, “Critical analysis of Goldberger-Wise stabilization of the Randall-Sundrum braneworld scenario”, Phys. Rev. D 75, 107901 (2007).
  * [38] T. Tanaka and X. Montes, “Gravity in the brane-world for two-branes model with stabilized modulus”, Nucl. Phys. B 582, 259 (2000).
  * [39] C. Csáki, M.L. Graesser, and G.D. Kribs, “Radion dynamics and electroweak physics”, Phys. Rev. D 63, 065002 (2001).
  * [40] T. Konstandin, G. Nardini, and M. Quiros, “Gravitational backreaction effects on the holographic phase transition”, Phys. Rev. D 82, 083513 (2010).
  * [41] A. Mazumdar and A. Pérez-Lorenzana, “A dynamical stabilization of the radion potential”, Phys. Lett. B 508, 340 (2001).
  * [42] S. Chakraborty and S. SenGupta, “Radion cosmology and stabilization”, Eur. Phys. J. C 74, 3045 (2014).
  * [43] I. Banerjee, S. Chakraborty, and S. SenGupta, “Radion induced inflation on nonflat brane and modulus stabilization”, Phys. Rev. D 99, 023515 (2019).
  * [44] G. Aad et al. [ATLAS Collaboration], “Search for extra dimensions using diphoton events in 7 TeV proton–proton collisions with the ATLAS detector”, Phys. Lett. B 710, 538 (2012).
  * [45] G. Aad et al. [ATLAS Collaboration], “Search for extra dimensions in diphoton events from proton–proton collisions at $\sqrt{s}=7\>\textrm{TeV}$ in the ATLAS detector at the LHC”, New J. Phys. 15, 043007 (2013).
  * [46] G. Aad et al. [ATLAS Collaboration], “Search for high-mass dilepton resonances in pp collisions at $\sqrt{s}=8\>\textrm{TeV}$ with the ATLAS detector”, Phys. Rev. D 90, 052005 (2014).
  * [47] V. Khachatryan et al. [CMS Collaboration], “Search for massive resonances decaying into pairs of boosted bosons in semi-leptonic final states at $\sqrt{s}=8\>\textrm{TeV}$”, J. High Energy Phys. 1408, 174 (2014).
  * [48] G. Aad et al. [ATLAS Collaboration], “Search for high-mass diphoton resonances in pp collisions at $\sqrt{s}=8\>\textrm{TeV}$ with the ATLAS detector”, Phys. Rev. D 92, 032004 (2015).
  * [49] [CMS Collaboration], “Search for High-Mass Diphoton Resonances in pp Collisions at $\sqrt{s}=8\>\textrm{TeV}$ with the CMS Detector”, CMS-PAS-EXO-12-045 (2015).
  * [50] S. Randjbar-Daemi and M. Shaposhnikov, “On some new warped brane world solutions in higher dimensions”, Phys. Lett. B 491, 329 (2000).
  * [51] P. Kanti, R. Madden, and K.A. Olive, “6-dimensional brane world model”, Phys. Rev. D 64, 044021 (2001).
  * [52] T. Gherghetta and A. Kehagias, “Anomaly-Free Brane Worlds in Seven Dimensions”, Phys. Rev. Lett.90, 101601 (2003).
  * [53] N. Kaloper, “Origami world”, J. High Energy Phys. 05, 061 (2004).
  * [54] B. Cuadros-Melgar and E. Papantonopoulos, “Need of dark energy for dynamical compactification of extra dimensions on the brane”, Phys. Rev. D 72, 064008 (2005).
  * [55] D. Choudhury and S. SenGupta, “Generalization of the Randall-Sundrum warped braneworld model to higher dimensions”, Phys. Rev. D 76, 064030 (2007).
  * [56] M.T. Arun, D. Choudhury, A. Das, and S. SenGupta, “Graviton modes in multiply warped geometry”, Phys. Lett. B 746, 266 (2015).
  * [57] M.T. Arun P. Saha, “Gravitons in multiply warped scenarios: At 750 GeV and beyond”, Pramana 88, 93 (2017).
  * [58] A. Das, R.S. Hundi, and S. SenGupta, “Bulk Higgs and gauge fields in a multiply warped braneworld model”, Phys. Rev. D 83, 116003 (2011).
  * [59] S. Chakraborty and S. SenGupta, “Bulk scalar field in warped extra dimensional models”, Phys. Rev. D 89, 126001 (2014).
  * [60] M.T. Arun and D. Choudhury, “Bulk gauge and matter fields in nested warping: I. The formalism”, J. High Energy Phys. 2015, 202 (2015).
  * [61] M.T. Arun and D. Choudhury, “Bulk gauge and matter fields in nested warping: II. Symmetry breaking and phenomenological consequences”, J. High Energy Phys. 2016, 133 (2016).
  * [62] M.T. Arun and D. Choudhury, “Stabilization of moduli in spacetime with nested warping and the UED”, Nucl. Phys. B 923, 258 (2017).
  * [63] T. Padmanabhan and D. Kothawala, “Lanczos–Lovelock models of gravity”, Phys. Rep. 531, 115 (2013).
  * [64] S. Chakraborty and T. Padmanabhan, “Geometrical variables with direct thermodynamic significance in Lanczos-Lovelock gravity”, Phys. Rev. D 90, 084021 (2014).
  * [65] S. Chakraborty, “Lanczos-Lovelock gravity from a thermodynamic perspective”, J. High Energy Phys. 08, 029 (2015).
  * [66] S. Nojiri, S.D. Odintsov, and S. Ogushi, “Cosmological and black hole brane-world universes in higher derivative gravity”, Phys. Rev. D 65, 023521 (2001).
  * [67] S. Capozziello and A. Troisi, “Parametrized post-Newtonian limit of fourth order gravity inspired by scalar-tensor gravity”, Phys. Rev. D 72, 044022 (2005).
  * [68] T.P. Sotiriou, “The nearly Newtonian regime in non-linear theories of gravity”, Gen. Rel. Grav. 38, 1407 (2006).
  * [69] S. Capozziello, A. Stabile, and A. Troisi, “Newtonian limit of $f(R)$ gravity”, Phys. Rev. D 76, 104019 (2007).
  * [70] S. Nojiri and S.D. Odintsov, “Unifying inflation with $\Lambda$CDM epoch in modified $f(R)$ gravity consistent with Solar System tests”, Phys. Lett. B 657, 238 (2007).
  * [71] L. Pogosian and A. Silvestri, “Pattern of growth in viable $f(R)$ cosmologies”, Phys. Rev. D 77, 023503 (2008). [Erratum: Phys. Rev. D 81, 049901 (2010)].
  * [72] J.D. Barrow and S. Cotsakis, “Inflation and the conformal structure of higher-order gravity theories”, Phys. Lett. B 214, 515 (1988).
  * [73] S. Capozziello, R. de Ritis, and A.A. Marino, “Some aspects of the cosmological conformal equivalence between the ‘Jordan frame’ and the ‘Einstein frame’”, Class. Quant. Grav. 14, 3243 (1997).
  * [74] A.D. Felice and S. Tsujikawa, “$f(R)$ Theories”, Living Rev. Rel. 13, 3 (2010).
  * [75] T.P. Sotiriou and V. Faraoni, “$f(R)$ theories of gravity”, Rev. Mod. Phys. 82, 451 (2010).
  * [76] I. Antoniadis, J. Rizos, and K. Tamvakis, “Singularity-free cosmological solutions of the superstring effective action”, Nucl. Phys. B 415, 497 (1994).
  * [77] R. Catena, M. Pietroni, and L. Scarabello, “Einstein and Jordan frames reconciled: A frame-invariant approach to scalar-tensor cosmology”, Phys. Rev. D 76, 084039 (2007).
  * [78] S. Capozziello, P. Martin-Moruno, and C. Rubano, “Physical non-equivalence of the Jordan and Einstein frames”, Phys. Lett. B 689, 117 (2010).
  * [79] T. Chiba and M. Yamaguchi, “Conformal-frame (in)dependence of cosmological observations in scalar-tensor theory”, JCAP 1310, 040 (2013).
  * [80] S. Bahamonde, S.D. Odintsov, V.K. Oikonomou, and M. Wright, “Correspondence of $F(R)$ gravity singularities in Jordan and Einstein frames”, Ann. Phys. 373, 96 (2016).
  * [81] S. Chakraborty and S. SenGupta, “Solving higher curvature gravity theories”, Eur. J. Phys. C 76, 552 (2016).
  * [82] S. Chakraborty and S. SenGupta, “Gravity stabilizes itself”, Eur. J. Phys. C 77, 573 (2017).
  * [83] I. Brevik, K.A. Milton, S. Nojiri, and S.D. Odintsov, “Quantum (in)stability of a brane-world $\textrm{AdS}_{5}$ universe at nonzero temperature”, Nucl. Phys. B 599, 305 (2001).
  * [84] I. Brevik, “Can the Hierarchy Problem be Solved by Finite-Temperature Massive Fermions in the Randall-Sundrum Model?”, arXiv:hep-th/0106225 (2001).
  * [85] Z. Chacko, R.K. Mishra, and D. Stolarski, “Dynamics of a stabilized radion and duality”, J. High Energy Phys. 09, 121 (2013).
  * [86] S. Anand, D. Choudhury, A.A. Sen, and S. SenGupta, “Geometric approach to modulus stabilization”, Phys. Rev. D 92, 026008 (2015).
  * [87] K. Fujikura, Y. Nakai, and M. Yamada, “A more attractive scheme for radion stabilization and supercooled phase transition”, J. High Energy Phys. 02, 111 (2020).


we deduce that $h_{l}$ is actually $L^{2}(\Omega)$-differentiable, with an
$L^{2}(\Omega)$-derivative denoted as $\dot{h}_{l}$, and that the derivative
process $\dot{h}_{l}$ is Gaussian, as it is an $L^{2}(\Omega)$-limit of
Gaussian random variables, satisfying
$\operatorname{\mathbb{E}}_{f\sim\Pi_{n}}\dot{h}_{l}(x)\dot{h}_{l}(y)=\left\langle\frac{\partial\widetilde{K}_{l}}{\partial
x}(x,\cdot),\frac{\partial\widetilde{K}_{l}}{\partial
x}(y,\cdot)\right\rangle_{\widetilde{\mathbb{H}}_{l}}$

Having established the existence of an $L^{2}(\Omega)$-derivative
$\dot{h}_{l}$ of the process $h_{l}$, we would like now to show that
$\dot{h}_{l}$ possesses a $(\gamma-1)$-regular version for every $\gamma<\nu$.
For this, we would like to apply Lemma 25 to $\dot{h}_{l}$. Notice that, by
isometry, for all $h>0$ we have

$\displaystyle\operatorname{\mathbb{E}}_{f\sim\Pi_{n}}\left|\dot{h}_{l}(x)-\dot{h}_{l}(y)\right|^{2}$
$\displaystyle=\left\lVert\frac{\partial\widetilde{K}_{l}}{\partial
x}(y,\cdot)-\frac{\partial\widetilde{K}_{l}}{\partial
x}(x,\cdot)\right\rVert_{\widetilde{\mathbb{H}}_{l}}^{2}$ (168)
$\displaystyle\leq
3\left\lVert\frac{\widetilde{K}_{l}(y+\delta,\cdot)-\widetilde{K}_{l}(y,\cdot)}{h}-\frac{\partial\widetilde{K}_{l}}{\partial
x}(y,\cdot)\right\rVert_{\widetilde{\mathbb{H}}_{l}}^{2}$ (169)
$\displaystyle~{}+3\left\lVert\frac{\widetilde{K}_{l}(x+\delta,\cdot)-\widetilde{K}_{l}(x,\cdot)}{h}-\frac{\partial\widetilde{K}_{l}}{\partial
x}(x,\cdot)\right\rVert_{\widetilde{\mathbb{H}}_{l}}^{2}$ (170)
$\displaystyle~{}+3\left\lVert\frac{\widetilde{K}_{l}(x+\delta,\cdot)-\widetilde{K}_{l}(x,\cdot)}{h}-\frac{\widetilde{K}_{l}(y+\delta,\cdot)-\widetilde{K}_{l}(y,\cdot)}{h}\right\rVert_{\widetilde{\mathbb{H}}_{l}}^{2}.$
(171)

Therefore, by the same arguments as above, we have

$\displaystyle\left(\\!\operatorname{\mathbb{E}}_{f\sim\Pi_{n}}\left|\dot{h}_{l}(x)\\!-\\!\dot{h}_{l}(y)\right|^{2}\right)^{1/2}\\!\\!\\!\\!\\!\\!\\!$
$\displaystyle\lesssim\liminf_{h\to
0}\left\lVert\frac{\widetilde{K}_{l}(x+\delta,\cdot)\\!-\\!\widetilde{K}_{l}(x,\cdot)}{h}\\!-\\!\frac{\widetilde{K}_{l}(y+\delta,\cdot)\\!-\\!\widetilde{K}_{l}(y,\cdot)}{h}\right\rVert_{\widetilde{\mathbb{H}}_{l}}$
(172) $\displaystyle\leq\liminf_{h\to
0}\sup_{\left\lVert\varphi\right\rVert_{\widetilde{\mathbb{H}}_{l}}=1}\int_{0}^{1}\left|\varphi^{\prime}(x+th)-\varphi^{\prime}(y+th)\right|\mathop{}\\!\mathrm{d}t$
(173) $\displaystyle\leq\liminf_{h\to
0}C_{\nu}\left|x-y\right|^{\nu-1}=C_{\nu}\left|x-y\right|^{\nu-1}$ (174)

where the transition from the second-to-last line to the last line is similar
to (163)–(165).

Therefore, we can apply Lemma 25 to $\dot{h}_{l}$, and find a version
$\tilde{h}_{l}^{\prime}$ of $\dot{h}_{l}$ with sample paths in
$\mathcal{CH}^{\alpha-1}\left(\mathcal{V}_{l}\right)$ almost surely for all
$\alpha<\nu$, such that

$\displaystyle\operatorname{\mathbb{E}}_{f\sim\Pi_{n}}\big{\lVert}\tilde{h}_{l}^{\prime}\big{\rVert}_{\mathcal{CH}^{\alpha-1}\left(\mathcal{V}_{l}\right)}^{2}\leq
C_{\nu,\alpha}<+\infty$ $\displaystyle\alpha<\nu.$ (175)

This gives Hölder regularity of the respective derivatives: we now integrate
these to obtain a Hölder-regular version of the process itself. Take any
$c_{l}\in(0,1)$ and consider
$\tilde{h}_{l}=h_{l}(c_{l})+\int_{c_{l}}^{\cdot}\tilde{h}_{l}^{\prime}(t)\mathop{}\\!\mathrm{d}t$.
Then since $\tilde{h}_{l}^{\prime}$ is almost surely in
$\mathcal{CH}^{\alpha-1}\left(\mathcal{V}_{l}\right)$, $\tilde{h}_{l}$ is has
sample paths almost surely in
$\mathcal{CH}^{\alpha}\left(\mathcal{V}_{l}\right)$. Moreover, it is easy to
check using our previous results that $\tilde{h}_{l}$ has an
$L^{2}(\Omega)$-derivative given by $\tilde{h}_{l}^{\prime}$. This implies
that $\tilde{h}_{l}$ is a version of $h_{l}$.

To conclude the argument, we construct $\tilde{f}$ from the obtained parts, by
pulling $\tilde{h}_{l}$ back from the charts to the manifold. Consider
$\tilde{f}=\sum_{l=1}^{L}\tilde{h}_{l}\circ\phi_{l}$. Arguing as in the case
$\nu\leq 1$, we see that $\tilde{f}$ is a version of $f$ with
$\mathcal{CH}^{\alpha}\left(\mathcal{M}\right)$ sample paths for every
$\alpha<\nu$, and for every $\alpha<\nu$ we have
$\operatorname{\mathbb{E}}_{f\sim\Pi_{n}}\lVert\tilde{f}\rVert_{\mathcal{CH}^{\alpha}\left(\mathcal{M}\right)}^{2}\leq
C_{\alpha,\nu}<+\infty$. This gives the claim. ∎

With this, it is easy to prove that all Matérn Gaussian processes considered
in this paper can be seen as Gaussian random elements in the Banach space
$\big{(}\mathcal{C}\left(\mathcal{M}\right),\left\lVert\cdot\right\rVert_{\mathcal{C}(\mathcal{M})}\big{)}$
of continuous functions on $\mathcal{M}$. This allows us to use the same proof
scheme as in [56] through the control of the _concentration functions_ defined
in Appendix B. It is also important that we work with Gaussian random elements
in $\mathcal{C}\left(\mathcal{M}\right)$—and not only with the classical
notion of Gaussian process, as the concentration functions are defined using
the _Gaussian random element RKHS_ defined in [59], which can potentially be
different from the classical RKHS. Fortunately, when the process is a Gaussian
random element in $\mathcal{C}\left(\mathcal{M}\right)$, [59, Theorem 2.1]
implies that the two notions of RKHS coincide.

######  Corollary 27.

The intrinsic Matérn Gaussian processes of Definition 4, their truncated
versions as in Section 3.2 as well as the extrinsic Matérn Gaussian processes
of Definition 5 are Gaussian random elements in
$\big{(}\mathcal{C}\left(\mathcal{M}\right),\left\lVert\cdot\right\rVert_{\mathcal{C}(\mathcal{M})}\big{)}$.

###### Proof.

By Lemma 15, it suffices to show that the processes have almost surely
continuous sample paths. Since Euclidean Matérn Gaussian processes have
continuous sample paths, this implies the same for their restrictions, the
extrinsic Matérn Gaussian processes on $\mathcal{M}$. For the intrinsic Matérn
process, we apply Lemma 26. ∎

Using Lemma 26 and known properties of Euclidean Matérn processes, we now
show, in a sense, that all of the Matérn processes presented in this paper are
sub-Gaussian, in a manner which holds uniformly with respect to the truncation
parameter in the case of the truncated intrinsic Matérn process, and live in
Hölder spaces with appropriate exponents. This result is used to control
Hölder norms when going from the error at input locations to the
$L^{2}(p_{0})$-error. We use the notation $\Pi_{n}$ to emphasize that the
prior depends on $n$ when we consider a truncated intrinsic Matérn process.

######  Lemma 28.

For $f\sim\Pi_{n}$ with $\Pi_{n}$ the prior in either Definition 4, Section
3.2 or Definition 5, for every $\nu>0$ and
$\gamma<\nu,\gamma\notin\mathbb{Z}_{>0}$, there exists a constant
$\sigma\left(f\right)=\sigma_{\gamma}\left(f\right)$ independent of $n$ we
have for $x>0$ that

$\operatorname{\mathbb{P}}\left(\left\lVert
f\right\rVert_{\mathcal{CH}^{\gamma}\left(\mathcal{M}\right)}>\left(x+1\right)\sigma\left(f\right)\right)\leq
2e^{-x^{2}/2}.$ (176)

###### Proof.

We start with the restriction $f$ of an extrinsic Matérn process $\tilde{f}$
to $\mathcal{M}$, as in Definition 5. By [56, Section 3.1], for every
$\gamma<\nu$ we have $\tilde{f}\in\mathcal{CH}^{\gamma}\left([0,1]^{D}\right)$
almost surely. By [20, Lemma I.7], for every $\gamma<\nu$, $\tilde{f}$ is a
Gaussian random element in the Banach space
$\mathcal{CH}^{\gamma}\left([0,1]^{D}\right)$. In particular, by the
Borell–TIS inequality [20, Proposition I.8] we have for $x>0$ that

$\operatorname{\mathbb{P}}\left(\big{\lVert}\tilde{f}\big{\rVert}_{\mathcal{CH}^{\gamma}\left([0,1]^{D}\right)}>\left(x+1\right)\sigma\big{(}\tilde{f}\big{)}\right)\leq
2e^{-x^{2}/2}$ (177)

where
$\sigma\big{(}\tilde{f}\big{)}=\left(\operatorname{\mathbb{E}}\,\big{\lVert}\tilde{f}\big{\rVert}_{\mathcal{CH}^{\gamma}\left([0,1]^{D}\right)}^{2}\right)^{1/2}<\infty$.
Since $\mathcal{M}$ is smooth, the restriction $f$ also satisfies for $x>0$
the expression

$\operatorname{\mathbb{P}}\left(\left\lVert
f\right\rVert_{\mathcal{CH}^{\gamma}\left(\mathcal{M}\right)}>\left(x+1\right)\sigma\left(f\right)\right)\leq
2e^{-x^{2}/2}$ (178)

perhaps for a possibly larger constant $\sigma\left(f\right)$. Finally, the
case of the intrinsic Matérn process $f\sim\Pi_{n}$ truncated at
$J_{n}\in\mathbb{Z}_{>0}\cup\left\\{\infty\right\\}$ follows in the same
manner, as we have shown in Lemma 26 that $\sup_{n\geq
1}\operatorname{\mathbb{E}}_{f\sim\Pi_{n}}\left\lVert
f\right\rVert_{\mathcal{CH}^{\alpha}\left(\mathcal{M}\right)}^{2}\leq
C_{\alpha,\nu}$. ∎

##  Appendix E Small Ball Asymptotics

Here, we bound the probability that a Matérn Gaussian process lies in the
$\varepsilon$-ball with respect to the $L^{\infty}$-norm for small
$\varepsilon$. For a Banach space $\mathbb{B}$, an element $x\in\mathbb{B}$
and a number $r\in\mathbb{R}_{>0}$, let us denote the closed $r$-ball around
$x$ by $B_{r}^{\mathbb{B}}(x)$. We start by an upper bound on the metric
entropy of Sobolev balls on $\mathcal{M}$ with respect to the uniform norm.

######  Lemma 29 (Entropy of Sobolev balls).

For any $s>0$ define the $\varepsilon$-covering number of $A\subseteq
H^{s}(\mathcal{M})$ with respect to the norm
$\lVert\cdot\rVert_{L^{\infty}\left(\mathcal{M}\right)}$ by

$N\left(\varepsilon,A,\left\lVert\cdot\right\rVert_{L^{\infty}\left(\mathcal{M}\right)}\right)=\operatorname*{\arg\min}_{J\in\mathbb{Z}_{>0}}\bigg{\\{}\exists
h_{1},..,h_{J}\in
A:A\subset\bigcup\limits_{j=1}^{J}B^{L^{\infty}(\mathcal{M})}_{\varepsilon}(h_{j})\bigg{\\}}.$
(179)

Then for any $s>d/2$, there exist $C,\varepsilon_{0}>0$ such that for every
$\varepsilon\leq\varepsilon_{0}$

$\ln
N\left(\varepsilon,B^{H^{s}(\mathcal{M})}_{1}(0),\lVert\cdot\rVert_{L^{\infty}\left(\mathcal{M}\right)}\right)\leq
C\varepsilon^{-\frac{d}{s}},$ (180)

where the left-hand side of the inequality above, as a function of
$\varepsilon$, is called the _metric entropy_ of the Sobolev ball
$B^{H^{s}(\mathcal{M})}_{1}(0)$ with respect to the uniform norm
$\lVert\cdot\rVert_{L^{\infty}\left(\mathcal{M}\right)}$.

###### Proof.

Using charts we will reduce the problem to the entropy of the unit ball of the
Sobolev space $H^{s}\left([0,1]^{d}\right)$ for which the upper bound is
known. Take $f\in B^{H^{s}(\mathcal{M})}_{1}(0)$ and consider approximations
of $f$ by $g$ of the form

$g=\sum_{l=1}^{L}\chi_{l}\left(g_{l}\circ\phi_{l}\right)$ (181)

for some functions $g_{l}:\mathcal{V}_{l}\to\mathbb{R}$ where
$\mathcal{V}_{l}\subseteq\mathbb{R}^{d}$. We have

$\displaystyle\big{\lVert}f-g\big{\rVert}_{L^{\infty}(\mathcal{M})}$
$\displaystyle=\big{\lVert}\sum_{l=1}^{L}\chi_{l}\left(g_{l}\circ\phi_{l}-f\right)\big{\rVert}_{L^{\infty}\left(\mathcal{M}\right)}\leq\sum_{l=1}^{L}\big{\lVert}\chi_{l}\left(g_{l}\circ\phi_{l}-f\right)\big{\rVert}_{L^{\infty}\left(\mathcal{U}_{l}\right)}$
(182)
$\displaystyle\leq\sum_{l=1}^{L}\big{\lVert}g_{l}\circ\phi_{l}-f\big{\rVert}_{L^{\infty}\left(\mathcal{U}_{l}\right)}\leq\sum_{l=1}^{L}\big{\lVert}g_{l}-f\circ\phi_{l}^{-1}\big{\rVert}_{L^{\infty}\left(\mathcal{V}_{l}\right)}$
(183) $\displaystyle\leq L\max_{1\leq l\leq
L}\big{\lVert}g_{l}-f\circ\phi_{l}^{-1}\big{\rVert}_{L^{\infty}\left([0,1]^{d}\right)}.$
(184)

This means that to approximate $f$ by $g$ uniformly on $\mathcal{M}$, we need
to choose the functions $g_{l}$ that approximate $f\circ\phi_{l}^{-1}$ well
with respect to the uniform norm on $[0,1]^{d}$.

Next, we show that the functions $f\circ\phi_{l}^{-1}$ are contained in an
Euclidean Sobolev ball of radius $R$, with $R$ depending only on $\nu$ and the
atlas. To do this we use [24, Lemma 2.1] to write

$\displaystyle\Big{\lVert}f\circ\phi_{l}^{-1}\Big{\rVert}_{H^{s}\left([0,1]^{d}\right)}$
$\displaystyle=\Big{\lVert}\sum_{l^{\prime}=1}^{L}\left(\chi_{l^{\prime}}f\right)\circ\phi_{l}^{-1}\Big{\rVert}_{H^{s}\left([0,1]^{d}\right)}\leq\sum_{l^{\prime}=1}^{L}\Big{\lVert}\left(\chi_{l^{\prime}}f\right)\circ\phi_{l}^{-1}\Big{\rVert}_{H^{s}\left([0,1]^{d}\right)}$
(185)
$\displaystyle=\sum_{l^{\prime}=1}^{L}\Big{\lVert}\left(\chi_{l^{\prime}}f\right)\circ\phi_{l^{\prime}}^{-1}\circ\phi_{l^{\prime}}\circ\phi_{l}^{-1}\Big{\rVert}_{H^{s}\left([0,1]^{d}\right)}$
(186)
$\displaystyle\lesssim\sum_{l^{\prime}=1}^{L}\Big{\lVert}\left(\chi_{l^{\prime}}f\right)\circ\phi_{l^{\prime}}^{-1}\Big{\rVert}_{H^{s}\left([0,1]^{d}\right)}\lesssim\left\lVert
f\right\rVert_{H^{s}\left(\mathcal{M}\right)}.$ (187)

Note, importantly, that the remark just above [24, Lemma 2.1] allows us to
consider Besov spaces $B_{2,2}^{s}$ coinciding with the Sobolev spaces $H^{s}$
instead of the Besov spaces $B_{2,\infty}^{s}$—to get from the second line to
the third. Note also that the constant hidden behind the notation $\lesssim$
in the last line where we use Theorem 13 is the radius $R$. Without loss of
generality we assume $R=1$. By the Euclidean counterpart of the result we are
proving [21, Theorem 4.3.36], we have

$\ln
N\left(\varepsilon,B^{H^{s}[0,1]^{d}}_{1}(0),\lVert\cdot\rVert_{L^{\infty}\left([0,1]^{d}\right)}\right)\lesssim\varepsilon^{-\frac{d}{s}}.$
(188)

Let $h_{1},..,h_{J}\in B^{H^{s}([0,1]^{d})}_{1}(0)$ be such that
$B^{H^{s}([0,1]^{d})}_{1}(0)\subset\bigcup_{j=1}^{J}B^{L^{\infty}([0,1]^{d})}_{\varepsilon/L}(h_{k})$.
Then for any $f\in B^{H^{s}(\mathcal{M})}_{1}(0)$ there exists a sequence
$\left\\{j_{l}\right\\}_{l=1}^{L}\subseteq\left\\{1,..,J\right\\}$ such that

$\big{\lVert}f-\sum_{l=1}^{L}\chi_{l}\left(h_{j_{l}}\circ\phi_{l}\right)\big{\rVert}_{L^{\infty}(\mathcal{M})}<L\frac{\varepsilon}{L}=\varepsilon.$
(189)

This shows that
$N\big{(}\varepsilon,B^{H^{s}(\mathcal{M})}_{1}(0),\left\lVert\cdot\right\rVert_{L^{\infty}\left(\mathcal{M}\right)}\big{)}\leq
LJ$, where $L$ is just the number of charts, thereby proving the claim. ∎

For the related _diffusion spaces_ [16], the RKHS corresponding to the heat
(diffusion) kernels, [12] uses the results of [14] to bound the entropy in
terms of a wavelet frame instead of relying on charts. We believe this
alternative proof scheme should work in our case as well. However, we could
not, to the best of our effort, get a tight-enough bound for the Sobolev
spaces by directly using the results of [14] and therefore we chose to rely on
charts instead.

Having established regularity properties for our prior processes, we now turn
to the _small ball problem_ : we want to find sharp lower bounds on
$\operatorname{\mathbb{P}}(\left\lVert
f\right\rVert_{L^{\infty}(\mathcal{M})}<\varepsilon)$ where $f\sim\Pi$ is our
prior process. This will be crucial in order to control the _concentration
functions_ used in Appendix B. In fact, it is well-known that this problem is
closely related to the estimation of the metric entropy of the unit ball of
the RKHS of $f$ with respect to the uniform norm: see [31] for details. Since
we have already characterized the RKHS of our processes in Proposition 24 and
Lemma 21, we are able to lower bound the small-ball probabilities. The
technicality here involves getting a bound uniform in the truncation parameter
for the truncated intrinsic Matérn process, as the truncated Matérn process is
a sequence of priors rather than a fixed prior.

######  Lemma 30.

If $f\sim\Pi_{n}$ where $\Pi_{n}$ is the prior in either Definitions 4 and 3.2
or Definition 5 with smoothness parameter $\nu>d/2$, then there exist two
constants $C,\varepsilon_{0}>0$ that do not depend on $n$ such that for all
$\varepsilon\leq\varepsilon_{0}$ we have
$-\ln\operatorname{\mathbb{P}}\left(\left\lVert
f\right\rVert_{L^{\infty}(\mathcal{M})}<\varepsilon\right)\leq
C\varepsilon^{-\frac{d}{\nu}}$.

###### Proof.

Because the processes are Gaussian random elements in
$\mathcal{C}\left(\mathcal{M}\right)$ by Lemma 15, their stochastic process
RKHS given by Proposition 24 coincide with their _Gaussian random element
RKHS_ defined in [59]. Hence, for the non-truncated intrinsic and the
extrinsic Matérn processes, the result follows by a direct application of
Lemma 29 and [31, Theorem 1.2]. For the intrinsic Matérn process truncated at
$J_{n}$ it is not immediately clear that the constants $C,\varepsilon_{0}$ can
be taken independent of $n$, so we go through the proof of [31, Proposition
3.1] to see that this is, in fact, the case. We first need a crude upper bound
of the form

$-\ln\operatorname{\mathbb{P}}\left(\left\lVert
f\right\rVert_{L^{\infty}(\mathcal{M})}<\varepsilon\right)\leq
c\varepsilon^{-c}$ (190)

for some possibly large constant $c>0$. To get such a bound, we use [12,
Proposition 3] which shows the existence of a universal constant $C>0$ such
that for all $\varepsilon\leq\min(1,4\sigma\left(f\right))$

$-\ln\operatorname{\mathbb{P}}\left(\left\lVert
f\right\rVert_{L^{\infty}(\mathcal{M})}<\varepsilon\right)\leq
Cn\left(\varepsilon\right)\ln\left(\frac{6n\left(\varepsilon\right)\max(1,\sigma\left(f\right))}{\varepsilon}\right)$
(191)

where
$\sigma\left(f\right)=\left(\operatorname{\mathbb{E}}_{f\sim\Pi_{n}}\left\lVert
f\right\rVert_{L^{\infty}(\mathcal{M})}^{2}\right)^{1/2}$ and
$n\left(\varepsilon\right)$ is defined in the following way in [12, page 684]
using auxiliary quantities $l_{J}$ that are defined in [31, page 1562], namely

$\displaystyle n(\varepsilon)$ $\displaystyle=\max\left\\{J\geq
0:4l_{J}(f)\geq\varepsilon\right\\},$ (192) $\displaystyle l_{J}(f)$
$\displaystyle=\inf\Big{\\{}\Big{(}\operatorname{\mathbb{E}}_{f\sim\Pi_{n}}\Big{\lVert}\sum\nolimits_{j\geq
J}\varepsilon_{j}h_{j}\Big{\rVert}_{L^{\infty}(\mathcal{M})}^{2}\Big{)}^{1/2}:f\overset{(d)}{=}\sum\nolimits_{j\geq
1}\varepsilon_{j}h_{j}\Big{\\}}$ (193)

with $\overset{(d)}{=}$ standing for the equality in distributions and the
infimum being taken over all possible decompositions $\sum_{j\geq
1}\varepsilon_{j}h_{j}$ with $h_{j}\in\mathcal{C}\left(\mathcal{M}\right)$,
$\varepsilon_{j}$ being a sequence of IID $\operatorname{N}(0,1)$ random
variables, and the series being required to converge uniformly almost
surely.999We consider $\sum_{j\geq 1}$, unlike $\sum_{j\geq 0}$ frequently
used above, in order to follow the respective references.

The function
$f=\sum_{j=1}^{J_{n}+1}\left(\frac{2\nu}{\kappa^{2}}+\lambda_{j-1}\right)^{-\frac{\nu+d/2}{2}}\varepsilon_{j}f_{j-1}$
is a valid decomposition. Therefore

$l_{J}(f)\leq\left(\operatorname{\mathbb{E}}_{f\sim\Pi_{n}}\Big{\lVert}\sum_{j=J}^{J_{n}+1}\left(\frac{2\nu}{\kappa^{2}}+\lambda_{j-1}\right)^{-\frac{\nu+d/2}{2}}\varepsilon_{j}f_{j-1}\Big{\rVert}_{L^{\infty}(\mathcal{M})}^{2}\right)^{1/2}.$
(194)

By the Sobolev Embedding Theorem and by Weyl’s Law, given in Result 7, for
every $d/2<\gamma<\nu$ there exists a constant $C=C_{\gamma,\mathcal{M}}$ such
that for all $J\in\mathbb{Z}_{>0}$, allowing $C$ to change from line to line,
we have

$\displaystyle\operatorname{\mathbb{E}}_{f\sim\Pi_{n}}\bigg{\lVert}\sum_{j=J}^{J_{n}+1}\left(1+\lambda_{j-1}\right)^{-\frac{\nu+d/2}{2}}\varepsilon_{j}f_{j-1}\bigg{\rVert}_{L^{\infty}(\mathcal{M})}^{2}\\!\\!\\!\\!\\!\\!\\!\\!\\!\\!\\!\\!\\!$
$\displaystyle\leq
C^{2}\operatorname{\mathbb{E}}_{f\sim\Pi_{n}}\bigg{\lVert}\sum_{j=J}^{J_{n}+1}\left(1+\lambda_{j-1}\right)^{-\frac{\nu+d/2}{2}}\varepsilon_{j}f_{j-1}\bigg{\rVert}_{H^{\gamma}\left(\mathcal{M}\right)}^{2}\\!\\!\\!\\!\\!\\!\\!\\!\\!\\!\\!\\!\\!\\!\\!\\!\\!$
(195)
$\displaystyle=C^{2}\sum_{j=J}^{J_{n}+1}\left(1+\lambda_{j-1}\right)^{-\left(\nu+d/2-\gamma\right)}$
(196) $\displaystyle\leq
C^{2}\sum_{j=J}^{J_{n}+1}j^{-\left(1+2\left(\nu-\gamma\right)/d\right)}$ (197)
$\displaystyle\leq C^{2}\sum_{j\geq
J}j^{-\left(1+2\left(\nu-\gamma\right)/d\right)}$ (198) $\displaystyle\leq
C^{2}J^{-2\left(\nu-\gamma\right)/d}.$ (199)

By choosing $J=1$ this gives us $\sigma\left(f\right)\leq C$ independent of
$n$. Moreover, by choosing $J\geq
C\varepsilon^{-\frac{d}{2\left(\nu-\gamma\right)}}$, again for a comparison
constant $C$ independent of $n$, this gives us $n\left(\varepsilon\right)\leq
C\varepsilon^{-\frac{d}{2\left(\nu-\gamma\right)}}$ for $C$ independent of
$n$. This implies using [12, Proposition 3] that

$-\ln\operatorname{\mathbb{P}}\left(\left\lVert
f\right\rVert_{L^{\infty}(\mathcal{M})}<\varepsilon\right)\leq
c\varepsilon^{-c}$ (200)

for $c>0$ independent of $n$.

With this crude bound, we can now continue the proof of [31, Proposition 3.1].
For this, we need a metric entropy estimate. For this notice that for all
$J\in\mathbb{Z}_{>0}\cup\left\\{\infty\right\\}$ we have
$B^{\mathbb{H}^{J}}_{1}(0)\subset
B^{\mathbb{H}^{\infty}}_{1}(0)=B^{H^{\nu+d/2}\left(\mathcal{M}\right)}_{1}(0)$,
and therefore using Lemma 29, we have the metric entropy estimate

$\ln
N\left(\varepsilon,B^{\mathbb{H}^{J}}_{1}(0),\left\lVert\cdot\right\rVert_{L^{\infty}\left(\mathcal{M}\right)}\right)\leq
C\varepsilon^{-\frac{d}{\nu+d/2}}$ (201)

for a constant $C>0$ independent of $J$. Therefore following the proof of
proposition 3.1 in [31] (with $J\equiv 1$) we find
$-\ln\operatorname{\mathbb{P}}\left(\left\lVert
f\right\rVert_{L^{\infty}(\mathcal{M})}<\varepsilon\right)\leq
C\varepsilon^{-\frac{d}{\nu}}$ for every $\varepsilon\leq\varepsilon_{0}$,
where $C,\varepsilon_{0}>0$ are constants independent of $n$. ∎

##  Appendix F Expressions for Pointwise Worst-case Errors

Let $k$ be a kernel on some abstract input domain $\mathcal{X}$, and let
$\mathbb{H}_{k}$ be the respective RKHS. Consider $n$ input values
$\mathbf{X}\subseteq\mathcal{X}$ and let $\sigma_{\varepsilon}^{2}>0$ be the
noise variance. Define

$\displaystyle m_{k,\mathbf{X},f,\varepsilon}(t)$
$\displaystyle=\mathbf{K}_{t\mathbf{X}}(\mathbf{K}_{\mathbf{X}\mathbf{X}}+\sigma_{\varepsilon}^{2}\mathbf{I})^{-1}\left(f(\mathbf{X})+\boldsymbol{\varepsilon}\right),$
(202) $\displaystyle v^{(\operatorname{i})}(t)=v_{k,\mathbf{X}}(t)$
$\displaystyle=k(t,t)-\mathbf{K}_{t\mathbf{X}}\left(\mathbf{K}_{\mathbf{X}\mathbf{X}}+\sigma_{\varepsilon}^{2}\mathbf{I}\right)^{-1}\mathbf{K}_{\mathbf{X}t}.$
(203)

######  Proposition 31.

With notation above

$v^{(\operatorname{i})}(t)=\sup_{f\in\mathcal{H}_{k},\left\lVert
f\right\rVert_{\mathcal{H}_{k}}\leq
1}\operatorname{\mathbb{E}}_{\boldsymbol{\varepsilon}\sim\operatorname{N}(\boldsymbol{0},\sigma_{\varepsilon}^{2}\mathbf{I})}\left|f(t)-m_{k,\mathbf{X},f,\varepsilon}(t)\right|^{2}.$
(204)

###### Proof.

To simplify notation, we shorten
$\operatorname{\mathbb{E}}_{\boldsymbol{\varepsilon}\sim\operatorname{N}(\boldsymbol{0},\sigma_{\varepsilon}^{2}\mathbf{I})}$
to $\operatorname{\mathbb{E}}$ and denote
$\boldsymbol{\alpha}=\mathbf{K}_{t\mathbf{X}}(\mathbf{K}_{\mathbf{X}\mathbf{X}}+\sigma_{\varepsilon}^{2}\mathbf{I})^{-1}$.
First of all, by direct computation,

$\displaystyle\operatorname{\mathbb{E}}m_{k,\mathbf{X},f,\varepsilon}(t)$
$\displaystyle=\boldsymbol{\alpha}f(\mathbf{X}),$ (205)
$\displaystyle\operatorname{\mathbb{E}}m_{k,\mathbf{X},f,\varepsilon}(t)^{2}$
$\displaystyle=\boldsymbol{\alpha}f(\mathbf{X})f(\mathbf{X})^{\top}\boldsymbol{\alpha}^{\top}+\sigma_{\varepsilon}^{2}\boldsymbol{\alpha}\boldsymbol{\alpha}^{\top}.$
(206)

Write

$\displaystyle\operatorname{\mathbb{E}}\left|f(t)-m_{k,\mathbf{X},f,\varepsilon}(t)\right|^{2}$
$\displaystyle=f(t)^{2}-2f(t)\operatorname{\mathbb{E}}m_{k,\mathbf{X},f,\varepsilon}(t)+\operatorname{\mathbb{E}}m_{k,\mathbf{X},f,\varepsilon}(t)^{2}$
(207)
$\displaystyle=f(t)^{2}-2f(t)\boldsymbol{\alpha}f(\mathbf{X})+\boldsymbol{\alpha}f(\mathbf{X})f(\mathbf{X})^{\top}\boldsymbol{\alpha}^{\top}+\sigma_{\varepsilon}^{2}\boldsymbol{\alpha}\boldsymbol{\alpha}^{\top}$
(208)
$\displaystyle=\left(f(t)-\boldsymbol{\alpha}f(\mathbf{X})\right)^{2}+\sigma_{\varepsilon}^{2}\boldsymbol{\alpha}\boldsymbol{\alpha}^{\top}$
(209)
$\displaystyle=\Big{\langle}k(t,\cdot)-\sum_{j=1}^{n}\alpha_{j}k(x_{j},\cdot),f\Big{\rangle}_{\mathbb{H}_{k}}^{2}+\sigma_{\varepsilon}^{2}\boldsymbol{\alpha}\boldsymbol{\alpha}^{\top}.$
(210)

As $\left\lVert
g\right\rVert_{\mathbb{H}_{k}}=\sup_{f\in\mathbb{H}_{k},\left\lVert
f\right\rVert_{\mathbb{H}_{k}}\leq 1}\left\langle
g,f\right\rangle_{\mathbb{H}_{k}}$, implying
$\sup_{f\in\mathbb{H}_{k},\left\lVert f\right\rVert_{\mathbb{H}_{k}}\leq
1}\left\langle g,f\right\rangle_{\mathbb{H}_{k}}^{2}=\left\lVert
g\right\rVert_{\mathbb{H}_{k}}^{2}$, we have

$\displaystyle\sup_{\begin{subarray}{c}f\in\mathcal{H}_{k}\\\ \left\lVert
f\right\rVert_{\mathcal{H}_{k}}\leq
1\end{subarray}}\operatorname{\mathbb{E}}\left|f(t)-m_{k,\mathbf{X},f,\varepsilon}(t)\right|^{2}$
$\displaystyle=\Big{\lVert}k(t,\cdot)-\sum_{j=1}^{n}\alpha_{j}k(x_{j},\cdot)\Big{\rVert}_{\mathbb{H}_{k}}^{2}+\sigma_{\varepsilon}^{2}\boldsymbol{\alpha}\boldsymbol{\alpha}^{\top}$
(211)
$\displaystyle=k(t,t)-2\boldsymbol{\alpha}\mathbf{K}_{\mathbf{X}t}+\underbracket{\boldsymbol{\alpha}\mathbf{K}_{\mathbf{X}\mathbf{X}}\boldsymbol{\alpha}^{\top}+\sigma_{\varepsilon}^{2}\boldsymbol{\alpha}\boldsymbol{\alpha}^{\top}}_{\vphantom{|}\mathclap{\boldsymbol{\alpha}\mathbf{K}_{\mathbf{X}t}}\vphantom{|}}$
(212)
$\displaystyle=k(t,t)-\boldsymbol{\alpha}\mathbf{K}_{\mathbf{X}t}=\underbracket{k(t,t)-\mathbf{K}_{t\mathbf{X}}(\mathbf{K}_{\mathbf{X}\mathbf{X}}+\sigma_{\varepsilon}^{2}\mathbf{I})^{-1}\mathbf{K}_{\mathbf{X}t}}_{\vphantom{|}\mathclap{v_{k,\mathbf{X}}(t)}\vphantom{|}}.$
(213)

∎

We now move to the misspecified case. Consider the RKHS $\mathcal{H}_{c}$ for
some other kernel $c:\mathcal{X}\times\mathcal{X}\to\mathbb{R}$ instead of
$\mathbb{H}_{k}$. Then, continuing from (210), write

$\sup_{\begin{subarray}{c}f\in\mathcal{H}_{c}\\\ \left\lVert
f\right\rVert_{\mathcal{H}_{c}}\leq
1\end{subarray}}\operatorname{\mathbb{E}}\left|f(t)-m_{k,\mathbf{X},f,\varepsilon}(t)\right|^{2}=\Big{\lVert}c(t,\cdot)-\sum_{j=1}^{n}\alpha_{j}c(x_{j},\cdot)\Big{\rVert}_{\mathbb{H}_{c}}^{2}+\sigma_{\varepsilon}^{2}\boldsymbol{\alpha}\boldsymbol{\alpha}^{\top}.$
(214)

The next question is how to compute the norm on the right-hand side. There is
not much hope of doing this exactly in the misspecified case: thus, we
consider approximations. To this end, we take some large set of locations
$\mathbf{X}^{\prime}\subseteq\mathcal{X}$. Then we use $\left\lVert
g\right\rVert_{\mathbb{H}_{c}}^{2}\approx
g(\mathbf{X}^{\prime})^{\top}\mathbf{C}_{\mathbf{X}^{\prime}\mathbf{X}^{\prime}}^{-1}g(\mathbf{X}^{\prime})$
for $g(\cdot)=c(t,\cdot)-\sum_{j=1}^{n}\alpha_{j}c(x_{j},\cdot)$. As a result,
we obtain the approximation

$\sup_{\begin{subarray}{c}f\in\mathcal{H}_{c}\\\ \left\lVert
f\right\rVert_{\mathcal{H}_{c}}\leq
1\end{subarray}}\operatorname{\mathbb{E}}\left|f(t)-m_{k,\mathbf{X},f,\varepsilon}(t)\right|^{2}\approx
g(\mathbf{X}^{\prime})^{\top}\mathbf{C}_{\mathbf{X}^{\prime}\mathbf{X}^{\prime}}^{-1}g(\mathbf{X}^{\prime})+\sigma_{\varepsilon}^{2}\boldsymbol{\alpha}\boldsymbol{\alpha}^{\top}=\tilde{v}_{k,c,\mathbf{X}}(t)=v^{(\operatorname{e})}(t)$
(215)

where $v^{(\operatorname{e})}(t)$ was first introduced in Section 4

To compute spatial averages of this quantity, let
$g_{t}(\cdot)=c(t,\cdot)-\sum_{j=1}^{n}\alpha_{j}c(x_{j},\cdot)$, the same as
$g$ before, but now with explicit dependence on $t$. Similarly, put
$\boldsymbol{\alpha}_{t}=\mathbf{K}_{t\mathbf{X}}(\mathbf{K}_{\mathbf{X}\mathbf{X}}+\sigma_{\varepsilon}^{2}\mathbf{I})^{-1}$.
Then

$\displaystyle g_{t}(\mathbf{X}^{\prime})$
$\displaystyle=\mathbf{C}_{\mathbf{X}^{\prime}\,t}-\mathbf{C}_{\mathbf{X}^{\prime}\,\mathbf{X}}\boldsymbol{\alpha}_{t}^{\top}=\mathbf{C}_{\mathbf{X}^{\prime}\,t}-\mathbf{C}_{\mathbf{X}^{\prime}\,\mathbf{X}}(\mathbf{K}_{\mathbf{X}\mathbf{X}}+\sigma_{\varepsilon}^{2}\mathbf{I})^{-1}\mathbf{K}_{\mathbf{X}t}$
(216) $\displaystyle
g_{t}(\mathbf{X}^{\prime})^{\top}\mathbf{C}_{\mathbf{X}^{\prime}\mathbf{X}^{\prime}}^{-1}g_{t}(\mathbf{X}^{\prime})$
$\displaystyle=\left(\mathbf{C}_{t\,\mathbf{X}^{\prime}}-\boldsymbol{\alpha}_{t}\mathbf{C}_{\mathbf{X}\,\mathbf{X}^{\prime}}\right)\mathbf{C}_{\mathbf{X}^{\prime}\mathbf{X}^{\prime}}^{-1}\left(\mathbf{C}_{\mathbf{X}^{\prime}\,t}-\mathbf{C}_{\mathbf{X}^{\prime}\,\mathbf{X}}\boldsymbol{\alpha}_{t}^{\top}\right).$
(217)

From here we can also deduce that

$\displaystyle\frac{1}{\left|\mathbf{X}^{\prime}\right|}\sum_{t\in\mathbf{X}^{\prime}}\tilde{v}_{k,c,\mathbf{X}}(t)$
$\displaystyle=\frac{1}{\left|\mathbf{X}^{\prime}\right|}\sum_{t\in\mathbf{X}^{\prime}}g_{t}(\mathbf{X}^{\prime})^{\top}\mathbf{C}_{\mathbf{X}^{\prime}\mathbf{X}^{\prime}}^{-1}g_{t}(\mathbf{X}^{\prime})$
(218)
$\displaystyle=\frac{1}{\left|\mathbf{X}^{\prime}\right|}\operatorname{tr}\left(g_{\mathbf{X}^{\prime}}(\mathbf{X}^{\prime})^{\top}\mathbf{C}_{\mathbf{X}^{\prime}\mathbf{X}^{\prime}}^{-1}g_{\mathbf{X}^{\prime}}(\mathbf{X}^{\prime})\right)$
(219)

where
$g_{\mathbf{X}^{\prime}}(\mathbf{X}^{\prime})=\mathbf{C}_{\mathbf{X}^{\prime}\,\mathbf{X}^{\prime}}\\!-\\!\mathbf{C}_{\mathbf{X}^{\prime}\,\mathbf{X}}(\mathbf{K}_{\mathbf{X}\mathbf{X}}\\!+\\!\sigma_{\varepsilon}^{2}\mathbf{I})^{-1}\mathbf{K}_{\mathbf{X}\mathbf{X}^{\prime}}$.

##  Appendix G Full Experimental Details

All of our kernels were computed using GPJax [38] and the geometric kernels
library.101010See https://gpjax.readthedocs.io and https://geometric-
kernels.github.io. We use three manifolds, each represented by a mesh: (i) a
dumbbell-shaped manifold represented as a mesh with $1556$ nodes, (ii) a
sphere represented by an icosahedral mesh with $2562$ nodes, and (iii) the
Stanford dragon mesh, preprocessed to keep only its largest connected
component, which has $100179$ nodes. For the sphere, we also considered a
finer icosahedral mesh with $10242$ nodes, but this was found to have
virtually no effect on the computed pointwise expected errors.

We use extrinsic Matérn and Riemannian Matérn kernels with the following
hyperparameters: $\sigma_{f}^{2}=1$ and $\sigma_{\varepsilon}^{2}=0.0005$. For
the truncated Karhunen–Loève expansion, we used $J=500$ eigenpairs obtained
from the mesh. We selected smoothness values to ensure norm-equivalence of the
intrinsic and extrinsic kernels’ reproducing kernel Hilbert spaces, which was
$\nu=5/2$ for the intrinsic model, and $\nu=5/2+d/2$ for the extrinsic model,
where $d$ is the manifold’s dimension. We used different length scales for
each manifold: $\kappa=200$ for the dumbbell, $\kappa=0.25$ for the sphere,
and $\kappa=0.05$ for the dragon, selected to ensure that the Gaussian
processes were neither approximately constant, nor white-noise-like. We
considered data sizes of $N=50$, $N=500$, and $N=1000$, respectively, for the
dumbbell, sphere, and dragon, sampled uniformly from the mesh’s nodes, which
in each case resulted in a reasonably-uniform distribution of points across
the manifold. Finally, for the extrinsic pointwise error approximation, we
used a subset $\mathbf{X}^{\prime}$ uniformly sampled from each mesh’s nodes,
of size equal to the data size. For each respective test set, we used the full
mesh. Each experiment was repeated for $10$ different seeds.

To set the length scales for the extrinsic process, we used maximum marginal
likelihood optimization on the full data, except for the dumbbell whose full
data size is small and for which we instead generated a larger set consisting
of $500$ points. We optimzied only the length scale, leaving all other
hyperparameters fixed. We used ADAM with a learning rate of $0.005$, and an
initialization equal to the length scale $\kappa$ of the intrinsic model,
except for the dumbbell where this lead to divergence and we instead used an
initial value of $\kappa/4$. We ran the optimizer for a total of $1000$ steps.
With these settings, we found empirically that maximum marginal likelihood
optimization always converged.