$\lambda^{2}=g(y;L)h(y,\Gamma_{\pm};L)\big{[}\partial_{y}X_{2}(y,\Gamma_{\pm};L)\pm\partial_{y}X_{1}(y,\Gamma_{\pm};L)].$

Since clearly $g(y;L)h(y,\Gamma;L)>0$, we need to study the sign of

$\mathcal{E}_{\pm}(y;L)=\partial_{y}X_{2}(y,\Gamma_{\pm}(y;L);L)\pm\partial_{y}X_{1}(y,\Gamma_{\pm}(y;L);L)=\pm\partial_{y}\mathcal{F}_{\pm}(y,\Gamma_{\pm};L)$

for $L\in[L_{\min},L_{\max}]$ and $\hat{y}_{\min}L\leq y\leq\frac{L}{2}$.

On the other hand, from $\mathcal{F}_{\pm}(y,\Gamma_{\pm}(y;L);L)=0$, we have
that

$\partial_{y}\mathcal{F}_{\pm}(y,\Gamma_{\pm}(y;L);L)+\partial_{\Gamma}\mathcal{F}_{\pm}(y,\Gamma_{\pm}(y;L);L)\partial_{y}\Gamma_{\pm}(y;L)=0$

and then

$\mathcal{E}_{\pm}(y;L)=\mp\partial_{\Gamma}\mathcal{F}_{\pm}(y,\Gamma_{\pm}(y;L);L)\partial_{y}\Gamma_{\pm}(y;L).$
(5.21)

######  Proposition 5.15.

We have that

  1. 1.

For $L\in[L_{\min},L_{\max}]$ and
$\Gamma\in[L\hat{\Gamma}_{\min}^{+},L\hat{\Gamma}_{\max}^{+}(L)]$, the unique
fixed point of the form $(0,y)$ satisfies that
$y=y_{+}(\Gamma;L):=L\hat{y}_{+}\left(\frac{\Gamma}{L};L\right)$ and it is a
saddle.

  2. 2.

For $L\in[L_{\min},L_{\max}]$ and
$\Gamma\in[L\hat{\Gamma}_{\min}^{-}(L),L\hat{\Gamma}_{\max}^{-}(L))$, the
unique fixed point of the form $(\pi,y)$ satisfies that
$y=y_{-}(\Gamma;L):=L\hat{y}_{-}\left(\frac{\Gamma}{L};L\right)$ and it is a
center.

  3. 3.

For $L\in(L_{*}(s_{\max}),L_{\max}]$ and
$\Gamma\in[L\hat{\Gamma}_{\max}^{-}(L),L\hat{\Gamma}_{*}(L))$ there are two
fixed points of the form $(\pi,y)$ with
$y=y_{-}^{1,2}(\Gamma;L):=L\hat{y}_{-}^{1,2}\left(\frac{\Gamma}{L};L\right)$.
Assume that $y^{1}<y_{\max}<y^{2}$. Then $(\pi,y_{-}^{1}(\Gamma;L))$ is a
center whereas $(\pi,y_{-}^{2}(\Gamma;L))$ is a saddle.

  4. 4.

When $L\in[L_{*}(s_{\max}),L_{\max}]$ and $\Gamma=L\hat{\Gamma}_{*}^{-}(L)$
the unique fixed point $(\pi,y_{-}(\Gamma;L))$ is parabolic.

###### Proof.

It is clear that if $\mathcal{E}_{\pm}>0$, the corresponding fixed point is a
saddle and if $\mathcal{E}_{\pm}<0$, it is a center.

Note that (5.21) implies
$\mathcal{E}_{\pm}=\mp\partial_{\Gamma}\mathcal{F}_{\pm}\partial_{y}\Gamma_{\pm}$
and Corollary 5.5 implies $\partial_{\Gamma}\mathcal{F}_{\pm}<0$. Then, when
$\partial_{y}\Gamma_{+}>0$ the fixed point will be a saddle and if
$\partial_{y}\Gamma_{\color[rgb]{1,0,0}\definecolor[named]{pgfstrokecolor}{rgb}{1,0,0}+}<0$
the fixed point will be a center. Conversely, when $\partial_{y}\Gamma_{-}>0$
the fixed point will be a center and if $\partial_{y}\Gamma_{-}<0$ the fixed
point will be a saddle.

From Proposition 5.7, we have that $\partial_{y}\Gamma_{+}>0$ which proves the
first item of the lemma.

For the values of $L,\Gamma$ in the second item, by Lemma 5.10 (see also
Figure 9), $\partial_{y}\Gamma_{-}=L\partial_{\hat{y}}\hat{\Gamma}_{-}>0$ that
implies that $\mathcal{E}_{-}<0$ and the result follows.

With respect to the third item, we have that $\hat{\Gamma}_{-}(\cdot;L)$ has a
maximum at $\hat{y}_{\max}$ so that
$\partial_{\hat{y}}\hat{\Gamma}_{-}(\hat{y};L)>0$ for $\hat{y}<\hat{y}_{\max}$
and negative if $\hat{y}>\hat{y}_{\max}$ and the result holds true since
$\hat{y}^{1}_{-}<\hat{y}_{\max}$ and $\hat{y}^{2}_{-}>\hat{y}_{\max}$.

The last item follows from the fact that, if $L=L_{*}(s)$ and
$\hat{\Gamma}=\hat{\Gamma}_{*}(L)$, then
$\partial_{\hat{y}}\hat{\Gamma}_{-}(\hat{y}_{*}(s);L_{*}(s))=0$ where we
recall that $\hat{y}_{*}(s)=\hat{y}_{\max}$ and
$\hat{\Gamma}_{*}^{-}(L)=\hat{y}_{*}(s)s$. ∎

Theorem 3.5 is a straightforward consequence of Proposition 5.15 keeping track
of the range of values of $\hat{\Gamma}_{\min,\max}^{\pm}$ and
$\hat{\Gamma}_{*}$. Indeed, we only need to rename
$\hat{\Gamma}_{1}=\hat{\Gamma}_{\min}^{-}$,
$\hat{\Gamma}_{2}=\Gamma_{\max}^{-}$,
$\hat{\Gamma}_{0}=\hat{\Gamma}_{\min}^{+}$,
$\hat{\Gamma}_{3}=\hat{\Gamma}_{\max}^{+}$,
$\hat{\Gamma}_{*}=\hat{\Gamma}_{*}^{-}$, $L_{*}=L_{*}(s_{\max})$ and to
rewrite Proposition 5.15 in the terms of Theorem 3.5, that is fixing the
values of $L$. Notice that by Remark 5.13 we have that

$\hat{\Gamma}_{1}(L)\leq\hat{\Gamma}_{0}\leq\hat{\Gamma}_{2}(L)\leq\hat{\Gamma}_{*}(L)\leq\Gamma_{3}(L),$

with the convection that $\hat{\Gamma}_{*}(L)=\Gamma_{2}(L)$ when
$L\in[L_{\min},L_{*}(s_{\max})]$. To finish, we note that by Proposition 5.7,
$\hat{\Gamma}_{3}$ is increasing and by Lemma 5.10 $\hat{\Gamma}_{1,2}$ are
decreasing.

##  6 Periodic orbits of the coplanar Hamiltonian

###  6.1 Existence of periodic orbits. Proof of Theorem 3.7

Consider the coplanar Hamiltonian
$\mathcal{H}_{\mathrm{CP}}=\mathcal{H}_{0}+\alpha^{3}\mathcal{H}_{\mathrm{CP},1}$.
We first emphasize that, by expressions of $\mathcal{H}_{0}$ in (2.14) and
${\mathcal{H}}_{\mathrm{CP},1}$ in (A.1) in Appendix A, $(\eta,\xi)=(0,0)$ is
invariant by the flow of $\mathcal{H}_{\mathrm{CP}}$ and every orbit of the
form $(0,\Gamma(t;L),0,h(t;L))$ has to be in an energy level
$\mathcal{H}_{\mathrm{CP}}(0,\Gamma,0,h)=\mathbf{E}$ for some energy
$\mathbf{E}$ with $\Gamma(t),h(t)$ satisfying the differential equations in
Theorem 3.7.

We now prove that for a certain range of energies $E$ (to be determined), the
solutions $(0,\Gamma(t;L),0,h(t;L))$ are periodic orbits which are a graph
over the variable $h$.

Fix $E\in\mathbb{R}$ and $L\in(0,L_{\max}]$. Writing $\Gamma=L\hat{\Gamma}$,

$\displaystyle E=\widehat{\mathcal{H}}_{\mathrm{CP}}(0,\hat{\Gamma},0,h):=$
$\displaystyle\frac{\rho_{0}}{16L^{6}}(1-12\hat{\Gamma}-12\hat{\Gamma}^{2})+\alpha^{3}\frac{\rho_{1}U_{2}^{0,0}}{32L^{2}}(1-12\hat{\Gamma}-12\hat{\Gamma}^{2})$
(6.1)
$\displaystyle-\alpha^{3}\frac{\rho_{1}U_{2}^{1,0}}{8L^{2}}\sqrt{(1-2\hat{\Gamma})(3+2\hat{\Gamma})}(2\hat{\Gamma}+1)\cos
h$
$\displaystyle-\alpha^{3}\frac{\rho_{1}U_{2}^{2,0}}{32L^{2}}(1-2\hat{\Gamma})(3+2\hat{\Gamma})\cos
2h.$

We impose that
$\partial_{\hat{\Gamma}}\widehat{\mathcal{H}}_{\mathrm{CP}}(0,\hat{\Gamma},0,h)\neq
0$ for all $h$, in other words $\dot{h}\neq 0$. This condition will give a set
of possible values for $L,\hat{\Gamma}$. From the differential equations in
Theorem 3.7, we need that

$\displaystyle\dot{h}=$
$\displaystyle-\frac{3\rho_{0}(1+2\hat{\Gamma})}{4L^{7}}-\alpha^{3}\frac{3\rho_{1}}{8L^{3}}U_{2}^{0,0}(1+2\hat{\Gamma})$
$\displaystyle-\alpha^{3}\frac{3\rho_{1}}{24L^{3}}\Big{(}4U_{2}^{1,0}\frac{1-4\hat{\Gamma}-4\hat{\Gamma}^{2}}{\sqrt{(1-2\hat{\Gamma})(3+2\hat{\Gamma})}}\cos
h-U_{2}^{2,0}(1+2\hat{\Gamma})\cos(2h)\Big{)}<0$

for all $h\in[0,2\pi]$ or equivalently

$\displaystyle(1+2\hat{\Gamma})+\frac{1}{2}\rho\alpha^{3}L^{4}U_{2}^{0,0}(1+2\hat{\Gamma})$
$\displaystyle+\frac{2}{3}\rho\alpha^{3}L^{4}U_{2}^{1,0}\frac{1-4\hat{\Gamma}-4\hat{\Gamma}^{2}}{\sqrt{(1-2\hat{\Gamma})(3+2\hat{\Gamma})}}\cos
h$ (6.2)
$\displaystyle-\frac{1}{6}\rho\alpha^{3}L^{4}U_{2}^{2,0}(1+2\hat{\Gamma})\cos
2h>0,$

where $\rho$ has been introduced in (4.5).

To avoid cumbersome notations, we introduce

$\displaystyle A=A(\hat{\Gamma};L)$
$\displaystyle=(1+2\hat{\Gamma})+\frac{1}{2}\rho\alpha^{3}L^{4}U_{2}^{0,0}(1+2\hat{\Gamma})$
$\displaystyle B=B(\hat{\Gamma};L)$
$\displaystyle=\frac{2}{3}\rho\alpha^{3}L^{4}U_{2}^{1,0}\frac{1-4\hat{\Gamma}-4\hat{\Gamma}^{2}}{\sqrt{(1-2\hat{\Gamma})(3+2\hat{\Gamma})}}$
$\displaystyle C=C(\hat{\Gamma};L)$
$\displaystyle=-\frac{1}{6}\rho\alpha^{3}L^{4}U_{2}^{2,0}(1+2\hat{\Gamma}),$

so that condition (6.2) reads as

$f(h)=f(h;\hat{\Gamma},L):=A+B\cos h+C\cos 2h>0.$

With respect to $h$, $f(h)$ has its global minimum either at $h=0,\pi$ or, if
$\left|\frac{B}{4C}\right|<1$ at $h_{1,2}$ satisfying

$\cos h_{1}=\cos h_{2}=-\frac{B}{4C}.$

Using that $\cos 2h=2\cos^{2}h-1$, we have that,

$f(h)=A+B\cos h+C(2\cos^{2}h-1)$

and then when $\left|\frac{B}{4C}\right|<1$

$f(h_{1})=f(h_{2})=A-\frac{B^{2}}{8C}-C.$

Since $A>0$ and $C<0$, $f(h_{1}),f(h_{2})$ are positive. Therefore we only
need to impose $f(0),f(\pi)>0$ Notice that $f(0)=A+B+C$ and $f(\pi)=A-B+C$ so
that both conditions can be written as

$A+C>|B|.$

When $|B|<4|C|$, this last condition is satisfied provided
$a_{\max}=30000\,\mathrm{km}$. Indeed, in this case

${5|C|}\leq\frac{{5}}{6}\rho\frac{a_{\max}^{3}}{a_{{\mathrm{M}}}^{3}}L_{\max}^{4}U_{2}^{2,0}(1+2\hat{\Gamma})=0.023691360650697(1+2\hat{\Gamma})<A$

and then $|B|<4|C|\leq 5|C|+C<A+C$. Summarizing we only need to impose

$A+C-|B|>0\quad\text{if }\quad\left|\frac{B}{4C}\right|\geq 1$ (6.3)

######  Lemma 6.1.

For $\hat{\Gamma}\in[0,0.49]$ and $L\in[L_{\min},L_{\max}]$, we have that
$\frac{5}{4}|B(\hat{\Gamma},L)|\leq 1<A$.

###### Proof.

It is clear that

$\frac{5}{4}|B(\hat{\Gamma},L)|\leq\frac{5}{6}\rho\frac{a_{\max}^{3}}{a_{{\mathrm{M}}}^{3}}L_{\max}^{4}U_{2}^{1,0}\frac{|1-4\hat{\Gamma}-4\hat{\Gamma}^{2}|}{\sqrt{(1-2\hat{\Gamma})(3+2\hat{\Gamma})}}\leq\frac{5}{3}\rho\frac{a_{\max}^{3}}{a_{{\mathrm{M}}}^{3}}L_{\max}^{4}U_{2}^{1,0}\frac{1}{\sqrt{3(1-2\cdot
0.49)}}.$

Computing this value, we have that $|B(\hat{\Gamma},L)|\leq 0.446157052927936$
and we are done. ∎

Lemma 6.1 implies that condition $A+C-|B|>0$ is always satisfied if
$\hat{\Gamma}\in[0,0.49]$ and $|B|\geq 4|C|$ because
$|B|-C\leq\frac{5}{4}|B|<A.$

As a consequence of the previous analysis,

$\partial_{\hat{\Gamma}}\widehat{\mathcal{H}}_{\mathrm{CP}}\neq
0\qquad\text{for}\quad
L\in[L_{\min},L_{\max}],\,h\in[0,2\pi],\,\hat{\Gamma}\in[0,0.49].$

Therefore, any energy level $\widehat{\mathcal{H}}_{\mathrm{CP}}=E$ which
belongs to the cylinder $h\in[0,2\pi],\hat{\Gamma}\in[0,0.49]$ is a closed
curve, which is a graph over $h$ and moreover its dynamics is periodic. Then,
it only remains to characterize such energy levels.

We notice that, since $L\hat{\Gamma}(t;L),h(t;L)$ satisfy the differential
equations in Theorem 3.7, one deduces that $\dot{\hat{\Gamma}}=0$ if and only
if $h=0,\pi,2\pi$. Indeed, $\dot{\hat{\Gamma}}=0$ is equivalent to

$\sin
h\left(U_{2}^{1,0}(2\hat{\Gamma}+1)+U_{2}^{2,0}\sqrt{(1-2\hat{\Gamma})(3+2\hat{\Gamma})}\cos
h\right)=0$

and since

$-\frac{U_{2}^{1,0}}{U_{2}^{2,0}}\frac{2\hat{\Gamma}+1}{\sqrt{(1-2\hat{\Gamma})(3+2\hat{\Gamma})}}\leq-\frac{U_{2}^{1,0}}{U_{2}^{2,0}\sqrt{3}}=-1.331627517097770<-1$

the conclusion is obvious. Clearly $h=\pi$ corresponds to a maximum and $h=0$
corresponds to a minimum of $\hat{\Gamma}$ as function of $h$.

We denote by $\hat{\Gamma}_{\min}(E),\hat{\Gamma}_{\max}(E)$ the values of
$\hat{\Gamma}$ such that

$\widehat{\mathcal{H}}_{\mathrm{CP}}(0,\hat{\Gamma}_{\min}(E),0,0)=E,\qquad\widehat{\mathcal{H}}_{\mathrm{CP}}(0,\hat{\Gamma}_{\max}(E),0,\pi)=E.$

The level curves $E=\widehat{\mathcal{H}}_{\mathrm{CP}}(0,\hat{\Gamma},0,h)$
cannot intersect and, since
$\partial_{\hat{\Gamma}}\widehat{\mathcal{H}}_{\mathrm{CP}}\neq 0$, we have
that
$\partial_{E}\hat{\Gamma}_{\min}(E),\partial_{E}\hat{\Gamma}_{\max}(E)\neq 0$.
Therefore, in order to check the range of $E$ allowed in our analysis, we
compute

$\mathbf{E}_{\max}=\mathbf{E}_{\max}(L)=\widehat{\mathcal{H}}_{\mathrm{CP}}(0,0,0,0),\qquad\mathbf{E}_{\min}=\mathbf{E}_{\min}=(L)\widehat{\mathcal{H}}_{\mathrm{CP}}(0,0.49,0,\pi)$

and Theorem 3.7 is proven, taking into account that, from (6.1), we easily
deduce that

$\mathbf{E}_{\max}(L)=\frac{\rho_{0}}{16L^{6}}\widehat{\mathbf{E}}_{\max}(\delta),\qquad\mathbf{E}_{\min}(L)=\frac{\rho_{0}}{16L^{6}}\widehat{\mathbf{E}}_{\min}(\delta)$

for some function $\widehat{\mathbf{E}}_{\min,\max}$ defined for
$\delta=\rho\alpha^{3}L^{4}\in[0,\delta_{\max}]$ (recall that $\delta_{\max}$
is defined in (3.3)).

To obtain the values in Remark 3.8, it is enough to recall that, for Galileo,
$a=29600\,\mathrm{km}$ (equivalently $L=1$), and therefore
$\mathbf{E}_{\min,\max}(1)$ with

$\mathbf{E}_{\max}(1)=2.477266122798186\cdot
10^{-6},\qquad\mathbf{E}_{\min}(1)=-2.558100888960067\cdot 10^{-5}.$

####  6.1.1 Proof of Corollary 3.10

We fix $\mathbf{E}$, a given energy level, and let
$\Gamma_{\mathrm{CP}}^{\mathbf{E}}$ be such that
$\mathcal{H}_{\mathrm{CP}}(0,\Gamma_{\mathrm{CP}}^{\mathbf{E}},0,0)=\mathbf{E}$.
Then, using formula (6.1) for
$\widehat{\mathcal{H}}_{\mathrm{CP}}(0,\hat{\Gamma},0,h)=\mathcal{H}_{\mathrm{CP}}(0,L\hat{\Gamma},0,h)$
one has that

$\frac{16L^{6}}{\rho_{0}}\mathbf{E}=\widehat{\mathcal{H}}_{\mathrm{AV}}(0,\hat{\Gamma}_{\mathrm{CP}}^{\mathbf{E}},0)+\delta\widehat{\mathcal{H}}^{*}_{\mathrm{CP},1}(0,\hat{\Gamma}_{\mathrm{CP}}^{\mathbf{E}},0,0),\qquad\hat{\Gamma}_{\mathrm{CP}}^{\mathbf{E}}=L^{-1}\Gamma_{\mathrm{CP}}^{\mathbf{E}},$

with $\delta=\rho\alpha^{3}L^{4}$,

$\displaystyle\widehat{\mathcal{H}}_{\mathrm{AV}}(0,\hat{\Gamma},0)$
$\displaystyle={\frac{16L^{6}}{\rho_{0}}\mathcal{H}_{\mathrm{AV}}(0,L\hat{\Gamma},0)}=(1-12\hat{\Gamma}-12\hat{\Gamma}^{2})\left(1+{\delta}\frac{U_{2}^{0,0}}{2}\right),$
$\displaystyle\widehat{\mathcal{H}}_{\mathrm{CP},1}^{{*}}(0,\hat{\Gamma},0,0)$
$\displaystyle={\frac{16L^{2}}{\rho_{1}}\left(\mathcal{H}_{\mathrm{CP},1}(0,L\hat{\Gamma},0,0)-\mathcal{H}_{\mathrm{AV},1}(0,L\hat{\Gamma},0)\right)}$
$\displaystyle=-2U_{2}^{1,0}\sqrt{(1-2\hat{\Gamma})(3+2\hat{\Gamma})}(2\hat{\Gamma}+1)-\frac{U_{2}^{2,0}}{2}(1-2\hat{\Gamma})(3+2\hat{\Gamma}).$

Let $\hat{\Gamma}_{\mathrm{AV}}^{\mathbf{E}}$ be such that

$\frac{16L^{6}}{\rho_{0}}\mathbf{E}=\widehat{\mathcal{H}}_{\mathrm{AV}}(0,\hat{\Gamma}_{\mathrm{AV}}^{\mathbf{E}},0).$

One has, from definition (3.10) of $\hat{a}_{0},\hat{c}_{0}$, that

$\hat{\Gamma}_{\mathrm{CP}}^{\mathbf{E}}-\hat{\Gamma}_{\mathrm{AV}}^{\mathbf{E}}=-\delta\frac{\widehat{\mathcal{H}}^{*}_{\mathrm{CP},1}(0,\hat{\Gamma}_{\mathrm{CP}}^{\mathbf{E}},0,0)}{\partial_{\hat{\Gamma}}\widehat{\mathcal{H}}_{\mathrm{AV}}(0,\hat{\Gamma}_{\mathrm{AV}}^{\mathbf{E}},0)}+\mathcal{O}(\delta^{2})=-\delta\hat{c}_{0}\big{(}\hat{\Gamma}_{\mathrm{CP}}^{\mathbf{E}};\delta\big{)}+\mathcal{O}(\delta^{2})$

and the proof is finished.

###  6.2 Perturbative analysis.

We fix $L\in[L_{\min},L_{\max}]$ and we consider the Hamiltonian
$\mathcal{H}_{\mathrm{CP}}$ in (2.16). Let
$(\eta,\Gamma,\xi,h)=(0,\Gamma(t;L),0,h(t;L))$ be a periodic orbit satisfying
$\Gamma(0;L)=\Gamma_{0}$ and $h(0;L)=0$ (see Theorem 3.7). We recall that in
(3.2) we have introduced the parameter $\delta=\alpha^{3}L^{4}\rho$.

####  6.2.1 Proof of Proposition 3.12

The proof of Proposition 3.12 relies on a simple perturbative analysis. In
order to make our analysis independent on $L$, we perform the change

$\hat{\Gamma}=\Gamma/L,\qquad s=3\rho_{0}t/(4L^{7})$

to the differential equation in Theorem 3.7 and we obtain the new system

$\displaystyle\hat{h}^{\prime}=$ $\displaystyle\,a_{0}(\hat{\Gamma};L)+\delta
a_{1}(\hat{\Gamma},\hat{h};L)$ $\displaystyle:=$
$\displaystyle-1-2\hat{\Gamma}-\alpha^{3}L^{4}\frac{\rho}{2}U_{2}^{0,0}(1+2\hat{\Gamma})+\delta\frac{1}{6}\Big{(}4U_{2}^{1,0}\frac{1-4\hat{\Gamma}-4\hat{\Gamma}^{2}}{\sqrt{(1-2\hat{\Gamma})(3+2\hat{\Gamma})}}\cos\hat{h}$
$\displaystyle-U_{2}^{2,0}(1+2\hat{\Gamma})\cos(2\hat{h})\Big{)}$ (6.4)
$\displaystyle\hat{\Gamma}^{\prime}=$
$\displaystyle-\delta\frac{1}{12}\Big{(}2U_{2}^{1,0}\sqrt{(1-2\hat{\Gamma})(3+2\hat{\Gamma})}(2\hat{\Gamma}+1)\sin\hat{h}+U_{2}^{2,0}(1-2\hat{\Gamma})(3+2\hat{\Gamma})\sin(2\hat{h})\Big{)}.$

By the analyticity of system (6.2.1) with respect to $\delta$, we have that

$\hat{\Gamma}(s)=\hat{\Gamma}_{0}+\delta\hat{\Gamma}_{1}(s)+\delta^{2}\hat{\Gamma}_{2}(s,\delta),\qquad\hat{h}(s)=\hat{h}_{0}(s)+\delta\hat{h}_{1}(s)+\delta^{2}\hat{h}_{2}(s,\delta)$

with $\hat{\Gamma}_{0}\in\mathbb{R}$ the initial condition and

$\displaystyle\hat{h}^{\prime}_{0}=$ $\displaystyle
a_{0}(\hat{\Gamma}_{0};L),$ $\displaystyle\hat{\Gamma}^{\prime}_{1}=$
$\displaystyle-\frac{1}{12}\Big{(}2U_{2}^{1,0}\sqrt{(1-2\hat{\Gamma}_{0})(3+2\hat{\Gamma}_{0})}(2\hat{\Gamma}_{0}+1)\sin\hat{h}_{0}(s)$
$\displaystyle+U_{2}^{2,0}(1-2\hat{\Gamma}_{0})(3+2\hat{\Gamma}_{0})\sin(2\hat{h}_{0}(s))\Big{)},$
$\displaystyle\hat{h}^{\prime}_{1}=$
$\displaystyle\partial_{\hat{\Gamma}}a_{0}(\hat{\Gamma}_{0};L)\hat{\Gamma}_{1}(s)+a_{1}(\hat{\Gamma}_{0},\hat{h}_{0}(s);L).$

Then, imposing that $\hat{h}(0)=0$ and $\hat{\Gamma}(0)=\hat{\Gamma}_{0}$,
that is $\hat{\Gamma}_{1}(0)=0$, we have that

$\displaystyle h_{0}(s)=$ $\displaystyle a_{0}(\hat{\Gamma}_{0};L)s$
$\displaystyle\hat{\Gamma}_{1}(s)=$ $\displaystyle
c_{0}(\hat{\Gamma}_{0};L)+c_{1}(\hat{\Gamma}_{0};L)\cos\big{(}a_{0}(\hat{\Gamma}_{0};L)s\big{)}+c_{2}(\hat{\Gamma}_{0};L)\cos\big{(}2a_{0}(\hat{\Gamma}_{0};L)s\big{)}$
$\displaystyle h_{1}(s)=$
$\displaystyle\partial_{\hat{\Gamma}}a_{0}(\hat{\Gamma}_{0};L)c_{0}(\hat{\Gamma}_{0};L)s+d_{1}(\hat{\Gamma}_{0};L)\sin\big{(}a_{0}(\hat{\Gamma}_{0};L)s\big{)}+d_{2}(\hat{\Gamma}_{0};L)\sin\big{(}2a_{0}(\hat{\Gamma}_{0};L)s\big{)}$

where it is straightforward to check that the constants
$a_{0},c_{0,1,2},d_{1,2}$ correspond to
$\hat{a}_{0},\hat{c}_{0,1,2},\hat{d}_{1,2}$, the ones defined in (3.10) and
(3.11), taking into account that the dependence with respect to the parameter
$L$ comes from the fact that

$a_{0}(\hat{\Gamma}_{0};L)=\hat{a}_{0}\left(\hat{\Gamma}_{0};\rho\alpha^{3}L^{4}\right).$

With respect to the period of the periodic orbit, if
$\hat{h}(\widehat{\mathcal{T}}(\hat{\Gamma}_{0};\delta))=2\pi$, that is
$\widehat{\mathcal{T}}(\hat{\Gamma}_{0};\delta)$ is the period of the periodic
orbit, then

$\widehat{\mathcal{T}}(\hat{\Gamma}_{0};\delta)=\frac{2\pi+\mathcal{O}(\delta^{2})}{|\hat{a}_{0}(\hat{\Gamma}_{0};\delta)+\delta\partial_{\hat{\Gamma}}\hat{a}_{0}(\hat{\Gamma}_{0};\delta)\hat{c}_{0}(\hat{\Gamma}_{0};\delta)|}=\frac{2\pi+\mathcal{O}(\delta)}{|\hat{a}_{0}(\hat{\Gamma}_{0};\delta)|}.$

Undoing the change of coordinates we have that
$\mathcal{T}(\Gamma_{0};L)=\widehat{\mathcal{T}}(\Gamma_{0}L^{-1};\delta){\frac{4L^{7}}{3\rho_{0}}}$
and

$\Gamma(t;L)=L\hat{\Gamma}\left({\frac{3\rho_{0}}{4L^{7}}}t;\rho\alpha^{3}L^{4}\right),\qquad
h(t;L)=h\left({\frac{3\rho_{0}}{4L^{7}}}t;\rho\alpha^{3}L^{4}\right).$

This completes the proof of Proposition 3.12 .

####  6.2.2 Proof of Theorem 3.15

Along this section, whenever there is no danger of confusion we will omit the
dependence on $L\in[L_{\min},L_{\max}]$.

It is convenient to introduce the (small) parameter $\epsilon=\alpha^{3}$ and
$X(\xi,\eta,\Gamma,h;\epsilon)$333With this new order in the variables the
linearized vector field around the periodic because has block structure., the
vector field associated to
$\mathcal{H}_{\epsilon}:=\mathcal{H}_{\mathrm{AV}}+\epsilon(\mathcal{H}_{\mathrm{CP},1}-\mathcal{H}_{\mathrm{AV},1})$.

The first step is to characterize the variational equation of $X$ around the
periodic orbit, $(0,0,\Gamma(t),h(t))$. Note that

$DX(0,0,\Gamma(t),h(t);\epsilon)=\left(\begin{array}[]{cc}\partial_{\xi,\eta}X_{\xi,\eta}(0,0,\Gamma(t),h(t);\epsilon)&0_{2\times
2}\\\ 0_{2\times
2}&\partial_{\hat{\Gamma},h}X_{\Gamma,h}(0,0,\Gamma(t),h(t);\epsilon)\end{array}\right).$

Then, if $\Phi(t;\epsilon)$ is the fundamental matrix of the linearized system

$\dot{z}=DX(0,0,\Gamma(t),h(t);\epsilon)z$

satisfying $\Phi(0;\epsilon)=\mathrm{Id}$, it has also the block form

$\Phi(t;\epsilon)=\left(\begin{array}[]{cc}\Psi(t;\epsilon)&0_{2\times 2}\\\
0_{2\times
2}&\widetilde{\Psi}(t;\epsilon)\end{array}\right),\qquad\Psi(0;\epsilon)=\widetilde{\Psi}(0;\epsilon)=\mathrm{Id}.$
(6.5)

Since the system is Hamiltonian, $\mathrm{det}\Phi(t;\epsilon)=1$ for all $t$
and therefore we also have that
$\mathrm{det}\Psi(t;\epsilon)=\mathrm{det}\widetilde{\Psi}(t;\epsilon)=1$.

The character of the periodic orbit $(0,0,\Gamma(t),h(t))$, is determined by
the eigenvalues of the monodromy matrix $\Psi(\mathcal{T};\epsilon)$, with
$\mathcal{T}=\mathcal{T}(\Gamma(0);L)$ the period of the periodic orbit. To
study the monodromy matrix, we first notice that

$\partial_{\xi,\eta}X_{\xi,\eta}(0,0,\Gamma(t),h(t);\epsilon)=\partial_{\xi,\eta}X_{\mathrm{AV}}(0,0,\Gamma(t))+\epsilon\frac{\rho_{1}}{L^{2}}\tilde{B}(t)$

where $X_{\mathrm{AV}}(\xi,\eta,\Gamma)$ is the vector field associated to the
$h$-averaged system $\mathcal{H}_{\mathrm{AV}}$ and, using formulas (2.19) and
(A.1) of $\mathcal{H}_{\mathrm{AV},1}$ and $\mathcal{H}_{\mathrm{CP},1}$
respectively, it is not difficult to see that the matrix $\tilde{B}$ has the
form

$\tilde{B}(t)=\left(\begin{array}[]{cc}\sum_{j=1}^{3}\tilde{b}_{11}^{j}(\Gamma(t),L)\sin
jh(t)&\sum_{j=1}^{3}\tilde{b}_{12}^{j}(\Gamma(t),L)\cos jh(t)\\\
\sum_{j=1}^{3}\tilde{b}_{21}^{j}(\Gamma(t),L)\cos
jh(t)&-\sum_{j=1}^{3}\tilde{b}_{11}^{j}(\Gamma(t),L)\sin
jh(t)\end{array}\right).$

The elements $\tilde{b}_{ij}$ depend on $L,\Gamma$ through the functions
$\mathcal{D}_{l,k}(0,\Gamma,0;L)$ (see Table 6). For future computations, we
note that $\tilde{b}_{ik}^{j}$ are linear functions on
$\mathcal{D}_{i,k}(0,\Gamma,0;L)$, namely

$\tilde{b}_{ik}^{j}(\Gamma,L)=\mathbf{b}_{ik}^{j}\big{(}\\{\mathcal{D}_{i,k}(0,\Gamma,0;L)\\}_{i=0,1,2,\,k=0,2},L\\{\mathcal{D}_{i,1}(0,\Gamma,0;L)\\}_{i=1,2}\big{)},$
(6.6)

with $\mathbf{b}_{ik}^{j}$ linear functions. Therefore the matrix
$\Psi(t;\epsilon)$ (see (6.5)) is the fundamental matrix of the two
dimensional linear system

$\dot{\zeta}=\left[\partial_{\xi,\eta}X_{\mathrm{AV}}(0,0,\Gamma(t))+\epsilon\frac{\rho_{1}}{L^{2}}\tilde{B}(t)\right]\zeta,\;\;\text{such
that }\;\;\Psi(0,\epsilon)=\mathrm{Id},\;\mathrm{det}\,\Psi(t,\epsilon)=1.$
(6.7)

######  Lemma 6.2.

Let $(0,0,\Gamma(t),h(t))$ be a $\mathcal{T}$-periodic orbit of the coplanar
Hamiltonian $\mathcal{H}_{\mathrm{CP}}$ with $\Gamma(0)=\Gamma_{0},h(0)=0$.

Consider the new variable $\hat{\Gamma}=\Gamma/L$, the new independent time
$s=3\rho_{0}t/(4L^{7})$, the initial condition $\hat{\Gamma}_{0}=\Gamma_{0}/L$
and the parameters $\delta$ and $\sigma_{0}$ defined by

$\delta=\epsilon
L^{4}\rho=\alpha^{3}L^{4}\frac{\rho_{1}}{\rho_{0}},\qquad\sigma_{0}=\sigma_{0}(\hat{\Gamma}_{0};\delta)=\hat{\Gamma}_{0}+\delta\hat{c}_{0}(\hat{\Gamma}_{0};\delta)$

with $\hat{c}_{0}$ introduced in (3.10). Then

  1. 1.

In the new variables, the variational equation (6.7) becomes

$\zeta^{\prime}=\left[\mathbf{X}(\sigma_{0})+\delta\widehat{B}(s)+\mathcal{O}(\delta^{2})\right]\zeta,\qquad\int_{0}^{\widehat{\mathcal{T}}}\widehat{B}(u)\,du=\mathcal{O}(\delta),$
(6.8)

where $\mathbf{X}(\sigma_{0})$ is defined as

$\displaystyle\mathbf{X}(\sigma_{0}):=$
$\displaystyle\frac{4L^{7}}{3\rho_{0}}\partial_{\xi,\eta}X_{\mathrm{AV}}(0,0,L\sigma_{0})$
(6.9) $\displaystyle=$
$\displaystyle\frac{1}{4}\left(\begin{array}[]{cc}0&-(\mathbf{a}(\sigma_{0})+\delta\mathbf{c}_{-}(\sigma_{0}))\\\
\mathbf{a}(\sigma_{0})+\delta\mathbf{c}_{+}(\sigma_{0})&0\end{array}\right)$

where $\mathbf{a}(\sigma)=1-16\sigma-20\sigma^{2}$ has been defined previously
in (4.9) and $\mathbf{c}_{\pm}$ have been introduced in Remark (4.1):

$\mathbf{c}_{\pm}(\sigma)=\frac{U_{2}^{0,0}}{2}\mathbf{a}(\sigma)\pm\frac{5U_{2}^{1,0}}{3}\mathbf{b}(\sigma),\qquad\mathbf{b}(\sigma)=(2\sigma+3)\sqrt{3-4\sigma-4\sigma^{2}}.$
(6.10)

  2. 2.

The fundamental matrix $\widehat{\Psi}(s,\delta)$ of (6.8) satisfies that

$\widehat{\Psi}(s,\delta)=\Psi\left(\frac{3\rho_{0}}{4L^{7}}s,\frac{\delta\rho_{0}}{\rho_{1}L^{4}}\right),\qquad\widehat{\Psi}(0,\delta)=\mathrm{Id},\qquad\mathrm{det}\widehat{\Psi}(s,\delta)=1.$

######  Remark 6.3.

The change of variables, presented in the previous lemma, allows to keep our
perturbative analysis uniform in $L\in[L_{\min},L_{\max}]$.

###### Proof.

Performing the change $\hat{\Gamma}=\Gamma/L$ and $s=4L^{7}t/(3\rho_{0})$,
system (6.7) becomes

$\zeta^{\prime}=\frac{4L^{7}}{3\rho_{0}}\partial_{\xi,\eta}X_{\mathrm{AV}}(0,0,L\hat{\Gamma}(s))\zeta+\delta\frac{4L}{3}\tilde{B}\left(\frac{4L^{7}}{3\rho_{0}}s\right)\zeta.$
(6.11)

Using Table 6 for the explicit expressions of $\mathcal{D}_{l,k}$, we obtain
that

$\displaystyle\mathcal{D}_{l,k}(0,L\hat{\Gamma},0;L)$
$\displaystyle=\frac{1}{L}\mathbf{d}_{l,k}(\hat{\Gamma}),$ $\displaystyle
l=0,1,2,\,k=0,2$ $\displaystyle\mathcal{D}_{l,1}(0,L\hat{\Gamma},0;L)$
$\displaystyle=\frac{1}{L^{2}}\mathbf{d}_{l,1}(\hat{\Gamma}),$ $\displaystyle
l=0,1,2,$

for smooth functions $\mathbf{d}_{l,k}$. Then, by formula (6.6) of
$\tilde{b}_{ij}$ (we recall that the functions $\mathbf{b}_{ij}$ are linear),
we conclude that the equation (6.11) can be rewritten as

$\zeta^{\prime}=\frac{4L^{7}}{3\rho_{0}}\partial_{\xi,\eta}X_{\mathrm{AV}}(0,0,L\hat{\Gamma}(s))\zeta+\delta\mathcal{B}(s)\zeta$
(6.12)

with

$\mathcal{B}(s)=\left(\begin{array}[]{cc}\sum_{j=1}^{3}{b}_{11}^{j}(\hat{\Gamma}(s))\sin
j\hat{h}(s)&\sum_{j=1}^{3}{b}_{12}^{j}(\hat{\Gamma}(s))\cos j\hat{h}(s)\\\
\sum_{j=1}^{3}{b}_{21}^{j}(\hat{\Gamma}(s))\cos
j\hat{h}(s)&-\sum_{j=1}^{3}b_{11}^{j}(\hat{\Gamma}(s))\sin
j\hat{h}(s)\end{array}\right).$

As a consequence, the fundamental matrix $\widehat{\Psi}(s,\delta)$ of (6.12)
satisfies that

$\widehat{\Psi}(s,\delta)=\Psi\left(\frac{4L^{7}}{3\rho_{0}}s,\frac{\delta\rho_{0}}{\rho_{1}L^{4}}\right),\qquad\widehat{\Psi}(0,\delta)=\mathrm{Id},\qquad\mathrm{det}\widehat{\Psi}(s,\delta)=1.$

This proves the second item in the lemma, provided that system (6.12)
coincides with (6.8).

Now we study $\partial_{\xi,\eta}X_{\mathrm{AV}}$. By (4.6),

$\partial_{\xi,\eta}X_{\mathrm{AV}}(0,0,L\hat{\Gamma})=\frac{\rho_{0}}{8L^{10}}\left(\begin{array}[]{cc}0&-X^{\eta}(L\hat{\Gamma};L)\\\
X^{\xi}(L\hat{\Gamma};L)&0\end{array}\right).$

Then, using the expressions (4.8) for $X^{\xi,\eta}(L\hat{\Gamma};L)$ and the
formulas for $d_{0},d_{1}$ in (4.2), we conclude that

$\frac{4L^{7}}{3\rho_{0}}\partial_{\xi,\eta}X_{\mathrm{AV}}(0,0,L\hat{\Gamma})=\frac{1}{4}\left(\begin{array}[]{cc}0&-(\mathbf{a}(L\hat{\Gamma})+\delta\mathbf{c}_{-}(L\hat{\Gamma}))\\\
\mathbf{a}(L\hat{\Gamma})+\delta\mathbf{c}_{+}(L\hat{\Gamma})&0\end{array}\right)$
(6.13)

with $\mathbf{a}(\sigma)=1-16\sigma-20\sigma^{2}$ (and
$\mathbf{b}(\sigma)=(2\sigma+3)\sqrt{3-4\sigma-4\sigma^{2}}$) defined in (4.9)
and

$\mathbf{c}_{\pm}(\sigma)=\frac{U_{2}^{0,0}}{2}\mathbf{a}(\sigma)\pm\frac{5U_{2}^{1,0}}{3}\mathbf{b}(\sigma).$

On the other hand, since Proposition 3.12 implies that
$\hat{\Gamma}(s)=\hat{\Gamma}_{0}+\delta\hat{c}_{0}(\hat{\Gamma}_{0})+\delta
b(s)+\mathcal{O}(\delta^{2})$, with

$\int_{0}^{\widehat{\mathcal{T}}}b(s)\,ds=\int_{0}^{\frac{2\pi}{|\hat{a}_{0}|}}b(s)\,ds+\mathcal{O}(\delta)=\mathcal{O}(\delta),$

we have that the variational equation for $\widehat{\Psi}$ in (6.12) is of the
form

$\zeta^{\prime}=\left[\frac{4L^{7}}{3\rho_{0}}\partial_{\xi,\eta}X_{\mathrm{AV}}(0,0,L(\hat{\Gamma}_{0}+\delta\hat{c}_{0}(\hat{\Gamma}_{0})))+\delta\widehat{B}(s)+\mathcal{O}(\delta^{2})\right]\zeta,\qquad\int_{0}^{\widehat{\mathcal{T}}}\widehat{B}(s)\,dt=\mathcal{O}(\delta)$

where
$\partial_{\xi,\eta}X_{\mathrm{AV}}(0,0,L(\hat{\Gamma}_{0}+\delta\hat{c}_{0}(\hat{\Gamma}_{0})))$
is a constant matrix (see (6.13)) studied in Section 4.2, and the
$\mathcal{O}(\delta^{k})$ are uniform in $\hat{\Gamma}_{0},L$. Taking
$\sigma_{0}=\sigma_{0}(\hat{\Gamma}_{0};\delta)=\hat{\Gamma}_{0}+\delta\hat{c}_{0}(\hat{\Gamma}_{0})$
and using again (6.13) the result holds true.

∎

We write the fundamental matrix $\widehat{\Psi}$ as

$\widehat{\Psi}(s;\delta)=\widehat{\Psi}_{0}(s;\sigma_{0})+\delta\widehat{\Psi}_{1}(s;\sigma_{0})+\delta^{2}\widehat{\Psi}_{2}(s;\sigma_{0},\delta),$

where, $\widehat{\Psi}_{0}(s;\sigma_{0}),\widehat{\Psi}_{1}(s;\sigma_{0})$
depend on $\delta$ through the parameter $\sigma_{0}$ and $\widehat{\Psi}_{2}$
is bounded.

The matrix $\widehat{\Psi}_{0}(s;\sigma_{0})$ is the fundamental matrix of the
constant coefficients linear system

$\zeta^{\prime}=\mathbf{X}(\sigma_{0})\zeta,\qquad\text{with}\qquad\widehat{\Psi}_{0}(0;\sigma_{0})=\mathrm{Id},$

with eigenvalues $e^{\pm\hat{\lambda}_{0}(\sigma_{0})}$ satisfying, from
expression (6.9) of $\mathbf{X}(\sigma_{0})$,

$|\hat{\lambda}_{0}(\sigma_{0})|=\frac{1}{4}\sqrt{\big{|}\mathbf{a}(\sigma_{0})+\delta\mathbf{c}_{+}(\sigma_{0})\big{|}\big{|}\mathbf{a}(\sigma_{0})+\delta\mathbf{c}_{-}(\sigma_{0})\big{|}}.$
(6.14)

The matrix $\widehat{\Psi}_{1}(s;\sigma_{0})$ satisfies

$\widehat{\Psi}_{1}^{\prime}(s;\sigma_{0})=\mathbf{X}(\sigma_{0})\widehat{\Psi}_{1}(s;\sigma_{0})+\widehat{B}(s)\widehat{\Psi}_{0}(s;\sigma_{0}),\qquad\widehat{\Psi}_{1}(0;\sigma_{0})=0$

and therefore

$\widehat{\Psi}_{1}(s;\sigma_{0})=\widehat{\Psi}_{0}(s;\sigma_{0})\int_{0}^{s}(\widehat{\Psi}_{0}(u;\sigma_{0})^{-1}\widehat{B}(u)\widehat{\Psi}_{0}(u;\sigma_{0})du.$
(6.15)

Whenever there is no danger of confusion, we will omit the dependence on
$\hat{\Gamma}_{0},\sigma_{0},\delta$. In particular, we write
$\widehat{\mathcal{T}}=\widehat{\mathcal{T}}(\hat{\Gamma}_{0},\delta)$,
$\widehat{\Psi}_{0}=\widehat{\Psi}_{0}(\widehat{\mathcal{T}},\sigma_{0})$,
$\widehat{\Psi}_{1}=\widehat{\Psi}_{1}(\widehat{\mathcal{T}},\sigma_{0})$ and
$\widehat{\Psi}=\widehat{\Psi}(\widehat{\mathcal{T}},\delta)$. In addition, we
will denote by $M$ a generic constant, independent on
$\hat{\Gamma}_{0},\delta$ that, along the proof, can change its value.

In order to prove Theorem 3.15, we need to compare the eigenvalues of
$\widehat{\Psi}$ with the ones of $\widehat{\Psi}_{0}$. To do so, we first
compare the trace of $\widehat{\Psi}$ with the one of $\widehat{\Psi}_{0}$.

######  Lemma 6.4.

Let $\hat{\lambda}_{0}=\hat{\lambda}_{0}(\sigma_{0})$ be such that the
eigenvalues of $\mathbf{X}(\sigma_{0})$ are $\pm\hat{\lambda}_{0}$. Then,
there exists $\delta_{0}>0$ such that, if $\delta\in[0,\delta_{0}]$,

$\mathrm{tr}(\widehat{\Psi})=\mathrm{tr}(\widehat{\Psi}_{0})+\mathcal{O}(\mathbf{L}\delta^{2})+\delta^{3}T(\delta)$

with $T(\delta)$ bounded on $[0,\delta_{0}]$ and

$\mathbf{L}:=\min_{k\in\mathbb{N}}\\{|\widehat{\mathcal{T}}|\hat{\lambda}_{0}|-\pi
k|\\}.$

######  Remark 6.5.

The quantity $\mathbf{L}$ is well defined and
$0\leq\mathbf{L}\leq\frac{\pi}{2}$. Indeed, for a given
$\widehat{\mathcal{T}}|\hat{\lambda}_{0}|$, the integer part
$k^{\prime}=\left[\frac{\widehat{\mathcal{T}}|\hat{\lambda}_{0}|}{\pi}\right]\in\mathbb{N}$
satisfies $\widehat{\mathcal{T}}|\hat{\lambda}_{0}|-k^{\prime}\pi\in[0,\pi)$.
Therefore if

$\widehat{\mathcal{T}}|\hat{\lambda}_{0}|-k^{\prime}\pi\in\left[0,\frac{\pi}{2}\right],\qquad\text{implies
}\qquad\mathbf{L}=\widehat{\mathcal{T}}|\hat{\lambda}_{0}|-k^{\prime}\pi\leq\frac{\pi}{2}.$

On the contrary, when
$\widehat{\mathcal{T}}|\hat{\lambda}_{0}|-k^{\prime}\pi\in\left[\frac{\pi}{2},\pi\right)$,
$\mathbf{L}=|\widehat{\mathcal{T}}|\hat{\lambda}_{0}|-(k^{\prime}+1)\pi|\leq\frac{\pi}{2}$.

###### Proof.

Let $\kappa>1$ big enough (its value will be determined later). We distinguish
two cases, namely $\mathbf{L}\geq\kappa\delta$ and
$\mathbf{L}\leq\kappa\delta$.

When $|\widehat{\mathcal{T}}\hat{\lambda}_{0}|\leq\kappa\delta$, from (6.14)
and Remark 3.14, we deduce that $|\mathbf{a}(\sigma_{0})|\leq M\kappa\delta$.
Then, since $\mathbf{a}(\sigma)=1-16\sigma-20\sigma^{2}$, its unique positive
zero is $\frac{m}{2}$, with $m$ defined in (3.4) and therefore the condition
$|\widehat{\mathcal{T}}\hat{\lambda}_{0}|\leq\kappa\delta$ implies that

$\left|\sigma_{0}-\frac{m}{2}\right|\leq M\kappa\delta$

for some constant $M$. Therefore, using that

$\mathbf{X}(\sigma_{0})=\mathbf{X}\left(\frac{m}{2}\right)+\mathcal{O}(\delta)=\mathcal{O}(\delta)$

we obtain that

$\widehat{\Psi}_{0}(\widehat{\mathcal{T}};\sigma_{0})=\widehat{\Psi}_{0}\left(\widehat{\mathcal{T}};\frac{m}{2}\right)+\mathcal{O}(\delta)=\mathrm{Id}+\mathcal{O}(\delta).$

If $|\widehat{\mathcal{T}}|\hat{\lambda}_{0}|-k\pi|\leq\kappa\delta$ for some
$k\in\mathbb{N}\backslash\\{0\\}$, by Theorem 3.2 the eigenvalues of
$\mathbf{X}(\sigma_{0})$ are of the form $\pm i|\hat{\lambda}_{0}|$, if
$\delta$ is small enough. Therefore, if $P$ is such that
$P\mathbf{X}(\sigma_{0})P^{-1}$ is a diagonal matrix,

$\displaystyle\widehat{\Psi}_{0}$
$\displaystyle=P\left(\begin{array}[]{cc}e^{i|\hat{\lambda}_{0}|\widehat{\mathcal{T}}}&0\\\
0&e^{-i|\hat{\lambda}_{0}|\widehat{\mathcal{T}}}\end{array}\right)P^{-1}=P\left(\begin{array}[]{cc}(-1)^{k}+\mathcal{O}(\delta)&0\\\
0&(-1)^{k}+\mathcal{O}(\delta)\end{array}\right)P^{-1}$
$\displaystyle=(-1)^{k}\mathrm{Id}+\mathcal{O}(\delta).$

Therefore, if $k$ is such that
$\mathbf{L}=|\widehat{\mathcal{T}}|\hat{\lambda}_{0}|-k\pi|$,
$\widehat{\Psi}_{0}=(-1)^{k}\mathrm{Id}+\mathcal{O}(\delta)$. Then, from
(6.15) and (6.8),

$\displaystyle\widehat{\Psi}_{1}$
$\displaystyle=((-1)^{k}\mathrm{Id}+\mathcal{O}(\delta))\int_{0}^{\widehat{\mathcal{T}}}((-1)^{k}\mathrm{Id}+\mathcal{O}(\delta))^{-1}\widehat{B}(u)((-1)^{k}\mathrm{Id}+\mathcal{O}(\delta))\,du$
$\displaystyle=((-1)^{k}\mathrm{Id}+\mathcal{O}(\delta))\left(\int_{0}^{\widehat{\mathcal{T}}}(-1)^{k}\widehat{B}(u)\,du+\mathcal{O}(\delta)\right)=\mathcal{O}(\delta).$

Therefore $\widehat{\Psi}=\widehat{\Psi}_{0}+\mathcal{O}(\delta^{2})$. We
write

$\widehat{\Psi}=\left(\begin{array}[]{cc}a&b\\\
c&d\end{array}\right)=\left(\begin{array}[]{cc}a_{0}+\delta^{2}a_{1}+\mathcal{O}(\delta^{3})&b_{0}+\delta^{2}b_{1}+\mathcal{O}(\delta^{3})\\\
c_{0}+\delta^{2}c_{1}+\mathcal{O}(\delta^{3})&d_{0}+\delta^{2}d_{1}+\mathcal{O}(\delta^{3})\end{array}\right)$

where the elements $a_{j},b_{j},c_{j},d_{j}$ are the corresponding ones for
$\widehat{\Psi}_{j}$, $j=0,1$. Imposing that
$\mathrm{det}(\widehat{\Psi})=\mathrm{det}(\widehat{\Psi}_{0})=1$, we have

$\delta^{2}(a_{0}d_{1}+a_{1}d_{0}-b_{0}c_{1}-b_{1}c_{0})+\mathcal{O}(\delta^{3})=0.$

Then, since $a_{0},d_{0}=1+\mathcal{O}(\delta)$ (in fact
$1+\mathcal{O}(\delta^{2}))$ and $b_{0},c_{0}=\mathcal{O}(\delta)$, we obtain
that $d_{1}+a_{1}=\mathcal{O}(\delta)$. As a consequence

$\mathrm{tr}(\widehat{\Psi})=\mathrm{tr}(\widehat{\Psi}_{0})+\mathcal{O}(\delta^{3})$

and, because we are assuming that $\mathbf{L}\leq\kappa\delta$, the lemma is
proven for this case.

Now we deal with the case
$|\widehat{\mathcal{T}}\hat{\lambda}_{0}|\geq\mathbf{L}\geq\kappa\delta$ with
$\kappa$ a large enough constant. If $\delta$ is small enough, by Theorem 3.2,
the eigenvalues of $\mathbf{X}(\sigma_{0})$ are of the form $\pm
i|\hat{\lambda}_{0}|$ with $|\lambda_{0}|$ in (6.14) and moreover, again by
Remark 3.14,

$\left|\sigma_{0}-\frac{m}{2}\right|\geq M\kappa\delta$

for some constant $M$. As a consequence, recalling that
$\mathbf{a}\left(\frac{m}{2}\right)=0$ and
$\partial_{\sigma}\mathbf{a}(\sigma)<0$, for $\sigma>0$, we have that

$|\mathbf{a}(\sigma_{0})|\geq\left|\mathbf{a}\left(\frac{m}{2}\pm\kappa\delta\right)\right|\geq
M\kappa\delta.$

In addition, it is clear from the definition of $\mathbf{c}_{\pm}$ in the
second item of Lemma 6.2 that, $|\mathbf{c}_{\pm}(\sigma_{0})|\leq M$. Then,

$\frac{\mathbf{a}(\sigma_{0})+\delta\mathbf{c}_{\pm}(\sigma_{0})}{\mathbf{a}(\sigma_{0})+\delta\mathbf{c}_{\mp}(\sigma_{0})}\geq\frac{1-\delta\frac{|\mathbf{c}_{\pm}(\sigma_{0})|}{|\mathbf{a}(\sigma_{0})}|}{1+\delta\frac{|\mathbf{c}_{\mp}(\sigma_{0})|}{|\mathbf{a}(\sigma_{0})|}}\geq
1-\frac{M}{\kappa},\qquad\frac{\mathbf{a}(\sigma_{0})+\delta\mathbf{c}_{\pm}(\sigma_{0})}{\mathbf{a}(\sigma_{0})+\delta\mathbf{c}_{\mp}(\sigma_{0})}\leq
1+\frac{M^{\prime}}{\kappa}$

for some constants $M,M^{\prime}$ and $\kappa$ large enough. This trivial fact
implies that, $P(\sigma_{0};\delta)$, the constant linear change of variables
that transforms $\mathbf{X}(\sigma_{0})$ in its diagonal form satisfies that
$\|P(\sigma_{0};\delta)\|,\|P^{-1}(\sigma_{0};\delta)\|$ are bounded uniformly
in $\delta$. After this (constant) linear change of variables,
$\widetilde{\Psi}=P\widehat{\Psi}P^{-1}=\widetilde{\Psi}_{0}+\delta\widetilde{\Psi}_{1}+\mathcal{O}(\delta^{2})$
with $\widetilde{\Psi}_{0}$

$\widetilde{\Psi}_{0}=\left(\begin{array}[]{cc}e^{i|\hat{\lambda}_{0}|\widehat{\mathcal{T}}}&0\\\
0&e^{-i|\hat{\lambda}_{0}|\widehat{\mathcal{T}}}\end{array}\right)=(-1)^{k}\left(\begin{array}[]{cc}e^{i\mathbf{L}}&0\\\
0&e^{-i\mathbf{L}}\end{array}\right)$

and
$\widetilde{\Psi}_{1}=\widetilde{\Psi}_{1}(\widetilde{\mathcal{T}};\sigma_{0})$
with

$\widetilde{\Psi}_{1}(s;\sigma_{0})=\widetilde{\Psi}_{0}(s;\sigma_{0})\int_{0}^{s}(\widetilde{\Psi}_{0}(u;\sigma_{0}))^{-1}P^{-1}(\sigma_{0},\delta)\widehat{B}(u)P(\sigma_{0},\delta)\widetilde{\Psi}_{0}(u;\sigma_{0})du.$

Using the bound for $\widehat{B}$ in (6.8) and the fact that $P,P^{-1}$ are
uniformly bounded with respect to $\delta$, we obtain that

$\int_{0}^{\widehat{\mathcal{T}}}P^{-1}(\sigma_{0};\delta)\widehat{B}(u)P(\sigma_{0};\delta)\,du=P^{-1}(\sigma_{0};\delta)\left[\int_{0}^{\widehat{\mathcal{T}}}\widehat{B}(u)\,du\right]P(\sigma_{0};\delta)=\mathcal{O}(\delta).$

Then, one easily checks that

$\displaystyle\widetilde{\Psi}$
$\displaystyle=\widetilde{\Psi}_{0}+\delta\widetilde{\Psi}_{1}+\mathcal{O}(\delta^{2})$
$\displaystyle=(-1)^{k}\left(\begin{array}[]{cc}e^{i\mathbf{L}}&0\\\
0&e^{-i\mathbf{L}}\end{array}\right)\left[\mathrm{Id}+\delta\left(\begin{array}[]{cc}0&{b}_{12}(\hat{\Gamma}_{0};\delta)\\\
b_{21}(\hat{\Gamma}_{0};\delta)&0\end{array}\right)\right]+\mathcal{O}(\delta^{2})$

with, writing $P^{-1}\widehat{B}(u)P=\big{(}\hat{b}_{ij}(u))_{i,j}$,

$b_{12}(\hat{\Gamma}_{0};\delta)=\int_{0}^{\widehat{\mathcal{T}}}\hat{b}_{12}(u)e^{-2i\mathbf{L}u}\,du,\qquad
b_{21}(\hat{\Gamma}_{0};\delta)=\int_{0}^{\widehat{\mathcal{T}}}\hat{b}_{12}(u)e^{2i\mathbf{L}u}\,du.$

Therefore, using (6.8),

$\widetilde{\Psi}=\widetilde{\Psi}_{0}\left[\mathrm{Id}+\mathbf{L}\delta\left(\begin{array}[]{cc}0&\beta_{1}(\hat{\Gamma}_{0})\\\
\beta_{2}(\hat{\Gamma}_{0})&0\end{array}\right)+\mathcal{O}(\delta^{2})\right]$

for some constants $\beta_{1},\beta_{2}$. We emphasize that, when
$\mathbf{L}\geq M$, namely is not small when $\delta$ is small, this
expression still holds true. We write now

$\widetilde{\Psi}=(-1)^{k}\left(\begin{array}[]{cc}e^{i\mathbf{L}}+a_{2}\delta^{2}&\beta_{1}\mathbf{L}\delta+b_{2}\delta^{2}\\\
\beta_{2}\mathbf{L}\delta+c_{2}\delta^{2}&e^{-i\mathbf{L}}+d_{2}\delta^{2}\end{array}\right)$

Therefore $\mathrm{det}(\widetilde{\Psi})=1$ implies that

$1=1+\delta^{2}(e^{i\mathbf{L}}d_{2}+e^{-i\mathbf{L}}a_{2})-\beta_{1}\beta_{2}\mathbf{L}^{2}\delta^{2}+\mathcal{O}(\delta^{3}\mathbf{L})+\mathcal{O}(\delta^{4}).$

Then, since $e^{\pm i\mathbf{L}}=1+\mathcal{O}(\mathbf{L})$

$a_{2}+d_{2}=\mathcal{O}(\mathbf{L})+\mathcal{O}(\delta\mathbf{L})+\mathcal{O}(\delta^{2})=\mathcal{O}(\mathbf{L}),$

where we have used that $\mathbf{L}\geq\kappa\delta$. As a consequence,

$\mathrm{tr}(\widehat{\Psi})=\mathrm{tr}(\widetilde{\Psi})=\mathrm{tr}(\widetilde{\Psi}_{0})+\mathcal{O}(\mathbf{L}\delta^{2})=\mathrm{tr}(\widehat{\Psi}_{0})+\mathcal{O}(\mathbf{L}\delta^{2})$

and the lemma is proven. ∎

We denote by $\mu=e^{\hat{\lambda}\widehat{\mathcal{T}}}$ one eigenvalue of
$\widehat{\Psi}$. Denoting $\tau=\mathrm{tr}(\widehat{\Psi})$ and using that
$\mathrm{det}(\widehat{\Psi})=1$, we can assume that

$\mu=\frac{1}{2}\tau+\frac{1}{2}\sqrt{\tau^{2}-4}.$

Let $\mu_{0}=e^{\hat{\lambda}_{0}\widehat{\mathcal{T}}}$ be such that

$\mu_{0}=\frac{1}{2}\tau_{0}+\frac{1}{2}\sqrt{\tau_{0}^{2}-4}.$

with $\tau_{0}=\mathrm{tr}(\widehat{\Psi}_{0})$. Then

$2|\mu-\mu_{0}|\leq|\tau-\tau_{0}|+\left|\sqrt{\tau^{2}-4}-\sqrt{\tau_{0}^{2}-4}\right|.$

By Lemma 6.4, $|\tau-\tau_{0}|\leq M\delta^{2}$, therefore, to prove Theorem
3.15 we have to analyze

$\Delta\mu:=\left|\sqrt{\tau^{2}-4}-\sqrt{\tau_{0}^{2}-4}\right|.$

######  Lemma 6.6.

Take $C_{1},C_{2}>0$ two constants. There exist $\delta_{0}>0$ and a constant
$C_{*}>0$ such that, for $\delta\in[0,\delta_{0}]$, we have that

  * •

When $\mathbf{L}\geq C_{1}\delta$, then $\Delta\mu\leq C_{*}\delta^{2}.$

  * •

If $C_{2}\delta^{3/2}\leq\mathbf{L}\leq C_{1}\delta$, then $\Delta\mu\leq
C_{*}\delta^{3}\mathbf{L}^{-1}\leq C_{*}\delta^{3/2}C_{2}^{-1}.$

  * •

Otherwise, that is if $\mathbf{L}\leq C_{2}\delta^{3/2}$, then $\Delta\mu\leq
C_{*}\delta^{3/2}$.

As a consequence $\Delta\mu\leq M\delta^{3/2}$ for some positive constant $M$.

###### Proof.

We recall that, by Remark 6.5, $\mathbf{L}\leq\frac{\pi}{2}$. Assume first
that $\mathbf{L}\geq\kappa\delta$ for some $\kappa$ big enough. Then by
Theorem 3.2, the eigenvalues of $\mathbf{X}(\hat{\Gamma}_{0})$ in (6.9) are of
the form $\pm i|\hat{\lambda}_{0}|$ with $|\hat{\lambda}_{0}|$ in (6.14).
Then,
$\tau_{0}=e^{i|\hat{\lambda}_{0}|\widehat{\mathcal{T}}}+e^{-i|\hat{\lambda}_{0}|\widehat{\mathcal{T}}}$
satisfies

$|\tau_{0}^{2}-4|=4-4\cos^{2}\big{(}\widehat{\mathcal{T}}\hat{\lambda}_{0}\big{)}=4\sin^{2}\big{(}\widehat{\mathcal{T}}\hat{\lambda}_{0}\big{)}=4\sin^{2}\mathbf{L}.$

In view of Lemma 6.4, $\tau=\tau_{0}+\mathbf{L}\delta^{2}\tau_{1}$ with
$|\tau_{1}|\leq M$. Then

$\displaystyle\Delta\mu$
$\displaystyle=\big{|}\sqrt{\tau_{0}^{2}-4}\big{|}\left|1-\sqrt{1+\frac{2\delta^{2}\mathbf{L}\tau_{0}\tau_{1}+\delta^{4}\mathbf{L}^{2}\tau_{1}^{2}}{\tau_{0}^{2}-4}}\right|\leq\frac{2\delta^{2}\mathbf{L}|\tau_{0}||\tau_{1}|}{\sqrt{\tau_{0}^{2}-4}}+M\delta^{2}|\tau_{1}|$
$\displaystyle\leq M\delta^{2}|\tau_{1}|.$

In the last bound we have used that $\sqrt{\tau_{0}^{2}-4}=2\sin\mathbf{L}\geq
M\mathbf{L}$, provided $0\leq\mathbf{L}\leq\frac{\pi}{2}$.

Now assume that $\mathbf{L}\leq\kappa\delta$. By Lemma 6.4,
$\tau=\tau_{0}+\delta^{3}\tau_{1}$ with $|\tau_{1}|\leq M$. We notice that,
using that either $\tau_{0}=\pm(e^{\mathbf{L}}+e^{-\mathbf{L}})$ or
$\tau_{0}=\pm(e^{i\mathbf{L}}+e^{-i\mathbf{L}})$. Then

$|\tau_{0}^{2}-4|=2\mathbf{L}^{2}+\mathcal{O}(\mathbf{L}^{4}).$

When $\mathbf{L}^{2}\geq M\delta^{3}$, for some suitable constant

$\displaystyle\Delta\mu$
$\displaystyle=|\sqrt{\tau_{0}^{2}-4}|\left|1-\sqrt{1+\frac{2\delta^{3}\tau_{0}\tau_{1}+\delta^{6}\tau_{1}^{2}}{\tau_{0}^{2}-4}}\right|\leq
M\delta^{3}\frac{|\tau_{1}|}{\sqrt{|\tau_{0}^{2}-4|}}$ $\displaystyle\leq
M\delta^{3}\mathbf{L}^{-1}|\tau_{1}|.$

On the other hand, when $\mathbf{L}^{2}\leq M\delta^{3}$, clearly,

$\Delta\mu\leq
2\sqrt{|\tau_{0}^{2}-4|}+\delta^{3/2}\sqrt{2|\tau_{0}||\tau_{1}|}+\delta^{3}|\tau_{1}|\leq
M\delta^{3/2}|\tau_{1}|^{1/2}.$

∎

Theorem 3.15 is a straightforward consequence of Lemma 6.6 just taking into
account that, using definitions (3.12), (3.16) and (3.15) of $\mathcal{T}$,
$\mathcal{T}^{(0)}$ and $\Gamma_{\mathrm{CP}}^{(0)}$ we have that, by Taylor’s
theorem,

$\widehat{\mathcal{T}}^{(0)}(\hat{\Gamma}_{\mathrm{CP}}^{(0)}(\hat{\Gamma}_{0};\delta);\delta)=\frac{4L^{7}}{3\rho_{0}}{\mathcal{T}}^{(0)}(\Gamma_{\mathrm{CP}}^{(0)}(\Gamma_{0};L);L)=\widehat{\mathcal{T}}(\hat{\Gamma}_{0};\delta)+\mathcal{O}(\delta^{2}).$

Therefore,

$e^{i|\hat{\lambda}_{0}|\widehat{\mathcal{T}}^{(0)}}=e^{i|\hat{\lambda}_{0}|\widehat{\mathcal{T}}}+\mathcal{O}(\delta^{2}).$

####  6.2.3 Proof of Lemma 3.18

To prove this result, by expressions (6.14) and (3.16) of
$|\hat{\lambda}_{0}|$ and $\mu_{\mathrm{CP}}^{(0)}$ respectively, we need to
solve the equation

$\frac{2\pi}{\big{|}\hat{a}_{0}(\sigma;\delta)\big{|}}\frac{1}{4}\sqrt{\big{|}\mathbf{a}(\sigma)+\delta\mathbf{c}_{+}(\sigma)\big{|}\big{|}\mathbf{a}(\sigma)+\delta\mathbf{c}_{-}(\sigma)\big{|}}=k\pi,$

with $\hat{a}_{0}$ defined in (3.10), for $k=1,2$. Using definition (6.10) of
$\mathbf{c}_{\pm}$, we rewrite it as

$H_{k}(\sigma;\delta)=\mathbf{a}^{2}(\sigma)+\delta\mathbf{a}^{2}(\sigma)U_{2}^{0,0}+\delta^{2}\mathbf{c}_{+}(\sigma)\mathbf{c}_{-}(\sigma)-4k^{2}\hat{a}_{0}^{2}(\sigma;\delta).$

Let $\hat{\Gamma}_{k}^{(0)}$ be the values such that this equation holds true
for $\delta=0$ (see Table 1 for the exact values of $\Gamma_{k}^{(0)}$). One
can easily check that
$\partial_{\sigma}H_{k}(\hat{\Gamma}_{k}^{(0)};\delta)\neq 0$ and therefore,
by the implicit function theorem, there exists $\hat{\Gamma}_{k}(\delta)$ such
that $H_{k}(\hat{\Gamma}_{k}(\delta);\delta)=0$ and the lemma holds true. In
addition a simply computation proves that

$\partial_{\delta}\hat{\Gamma}_{k}(0)=-\frac{\partial_{\delta}H_{k}(\hat{\Gamma}_{k}^{(0)};0)}{\partial_{\sigma}H_{k}(\hat{\Gamma}_{k}^{(0)};0)}=0$

and therefore
$\hat{\Gamma}_{k}(\delta)=\hat{\Gamma}^{(0)}_{k}+\mathcal{O}(\delta^{2})$.

####  6.2.4 Proof of Theorem 3.19

The proof of Theorem 3.19 is a nontrivial consequence of Lemma 6.6. Fix
$\delta_{0}>0$ small enough (to be determined later) and consider
$\delta\in[0,\delta_{0}]$. By Theorem 3.2 and Lemma 3.18, there exist
$\sigma_{\pm}$ and $\sigma_{k}$, $k=1,2$ such that

$\hat{\lambda}_{0}(\sigma_{\pm})=0,\qquad\widehat{\mathcal{T}}^{(0)}(\sigma_{k};\delta)|\hat{\lambda}_{0}(\sigma_{k})|=\frac{2\pi|\hat{\lambda}_{0}(\sigma_{k})|}{|\hat{a}_{0}(\sigma_{k};\delta)|}=k\pi,\;k=1,2$

with $\hat{a}_{0}$ defined in (3.10). Notice also that

$\sigma_{\pm}=\frac{m}{2}\pm\delta
C_{0}+\mathcal{O}(\delta^{2}),\qquad\sigma_{k}=\hat{\Gamma}_{k}^{(0)}+C_{k}\delta+\mathcal{O}(\delta^{2}).$

Therefore, since $1\leq\nu\leq 2$, we can assume that the parameter $\sigma$
satisfies, taking $\delta$ small enough,

$|\sigma-\sigma_{\pm}|\geq C\delta^{\nu},\qquad|\sigma-\sigma_{k}|\geq
C\delta^{\frac{1+\nu}{2}}.$

We now observe that, using the definition of $\mathbf{c}_{\pm}$ in (6.10)

$16|\hat{\lambda}_{0}(\sigma)|^{2}=\mathbf{a}^{2}(\sigma)(1+\delta
U_{2}^{0,0})+\mathcal{O}(\delta^{2})$ (6.16)

with $\mathbf{a}=1-16\sigma-20\sigma^{2}$. In addition

$\widehat{\mathcal{T}}^{(0)}(\sigma;\delta):=\frac{2\pi}{|\hat{a}_{0}(\sigma;\delta)|}=\widehat{\mathcal{T}}(\sigma;\delta)+\mathcal{O}(\delta^{2}).$
(6.17)

We distinguish the cases $|\mathbf{a}(\sigma)|\geq\delta|\log\delta|$ and
$|\mathbf{a}(\sigma)|<\delta|\log\delta|$. We recall that
$\mathbf{a}\left(\frac{m}{2}\right)=0$. Then, the cases considered are
equivalent to $\left|\sigma-\frac{m}{2}\right|\geq\kappa\delta|\log\delta|$
and $\left|\sigma-\frac{m}{2}\right|<\kappa\delta|\log\delta|$ for some
constant $\kappa$.

We start with $\left|\sigma-\frac{m}{2}\right|\geq\kappa\delta|\log\delta|$.
From (6.16) we have that

$4|\hat{\lambda}_{0}(\sigma)|=|\mathbf{a}(\sigma)|\big{(}1+\mathcal{O}(|\log\delta|^{-1})\big{)}.$
(6.18)

  * •

If $\sigma\in\left(0,\frac{m}{2}-\kappa\delta|\log\delta|\right)$, by
Proposition 3.3, the function $|\hat{\lambda}_{0}|$ is decreasing. Therefore,
taking into account (6.17) and definition 3.10 of $\hat{a}_{0}$, if $\delta$
is small enough,

$\frac{1}{2}\left|\hat{\lambda}_{0}\left(\frac{m}{2}-\delta|\log\delta|\right)\right|\widehat{\mathcal{T}}\left(\frac{m}{2}-\delta|\log\delta|;\delta\right)\leq|\hat{\lambda}_{0}(\sigma)|\widehat{\mathcal{T}}(\sigma;\delta)\leq\frac{3}{2}|\hat{\lambda}_{0}(0)|\widehat{\mathcal{T}}(0;\delta)$

Then
$|\hat{\lambda}_{0}(\sigma)|\widehat{\mathcal{T}}(\sigma;\delta)|\leq\frac{3\pi}{4}(1+M\delta)$
and $|\hat{\lambda}_{0}(\sigma)|\widehat{\mathcal{T}}(\sigma;\delta)\geq
M\delta|\log\delta|$. Those bounds imply that, for all $k$

$\big{|}|\hat{\lambda}_{0}(\sigma)|\widehat{\mathcal{T}}(\sigma;\delta)-k\pi\big{|}\geq
M\delta|\log\delta|$

and, in conclusion,

$\mathbf{L}\geq M\delta|\log\delta|.$

By Lemma 6.6, $|\mu-\mu_{0}|\leq M\delta^{2}$ and the result follows for this
case.

  * •

If $\sigma>\frac{m}{2}+\delta|\log\delta|$, we notice that, using (6.18) and
(6.17),

$\displaystyle\widehat{\mathcal{T}}(\sigma;\delta)|\hat{\lambda}_{0}(\sigma)|-k\pi=$
$\displaystyle\widehat{\mathcal{T}}(\sigma;\delta)|\hat{\lambda}_{0}(\sigma)|-\widehat{\mathcal{T}}(\sigma_{k};\delta)|\hat{\lambda}_{0}(\sigma_{k})|$
$\displaystyle=$
$\displaystyle\left[\partial_{\sigma}\left(\frac{1}{4}\frac{|\mathbf{a}(s)|}{|\hat{a}_{0}(s;0)|}\right)+\mathcal{O}(|\log\delta|^{-1})\right]|\sigma-\sigma_{k}|$
$\displaystyle\widehat{\mathcal{T}}(\sigma;\delta)|\hat{\lambda}_{0}(\sigma)|=$
$\displaystyle\left[\partial_{\sigma}\left(\frac{1}{4}\frac{|\mathbf{a}(t)|}{|\hat{a}_{0}(t;0)|}\right)+\mathcal{O}(|\log\delta|^{-1})\right]|\sigma-\sigma_{+}|,$

with $s\in\overline{\sigma,\sigma_{k}}$, $t\in\overline{\sigma,\sigma_{\pm}}$
and $\hat{a}_{0}$ defined in (3.10). Since
$\sigma,\sigma_{k},\sigma_{+}>\frac{m}{2}$ (see Lemma 3.18),

$\partial_{\sigma}\left(\frac{1}{4}\frac{|\mathbf{a}(s)|}{|\hat{a}_{0}(s;0)|}\right)=\frac{1}{4|\hat{a}_{0}(s;0)|^{2}}(40s^{2}+40s+18)\geq
M.$

Therefore $\mathbf{L}\geq
M\min\left\\{\delta|\log\delta|,\delta^{\frac{1+\nu}{2}}\right\\}=M\delta^{\frac{1+\nu}{2}}\geq
M\delta^{\frac{3}{2}}$. Therefore, Lemma 6.6 assures that

$|\mu-\mu_{0}|\leq M\delta^{3}\mathbf{L}^{-1}\leq M\delta^{\frac{5-\nu}{2}}.$

Now we deal with the case
$\left|\sigma-\frac{m}{2}\right|\leq\kappa\delta|\log\delta|$, which in
particular implies that $|\hat{\lambda}_{0}(\sigma)|\leq M\delta|\log\delta|$.
Then, by continuity, in this case
$\mathbf{L}=\widehat{\mathcal{T}}|\hat{\lambda}_{0}(\sigma)|$, if $\delta$ is
small enough. We first observe that, $\sigma_{\pm}=\Gamma_{\pm}L^{-1}$ with
$\Gamma_{\pm}$ defined in Theorem 3.2 and $\sigma_{\pm}=\frac{m}{2}\pm
C_{0}\delta+\mathcal{O}(\delta^{2})$ for some constant $C_{0}$ (see (3.18)).
Therefore $|\sigma-\sigma_{\pm}|\leq\kappa\delta|\log\delta|$ for a
(different) constant $\kappa$.

Assume now that
$C\delta^{\nu}\leq|\sigma-\sigma_{-}|\leq|\sigma-\sigma_{+}|\leq\kappa\delta|\log\delta|$.
We have that (see Theorem 3.2)

$\mathbf{a}(\sigma_{-})+\delta\mathbf{c}_{-}(\sigma_{-})=0.$

Taking $\Delta\sigma=\sigma-\sigma_{-}=\mathcal{O}(\delta|\log\delta|)$, we
have that

$\displaystyle|\hat{\lambda}_{0}(\sigma)|^{2}=$
$\displaystyle\left|\partial_{\sigma}\mathbf{a}(\sigma_{-})\Delta\sigma+\Delta\sigma\mathcal{O}(\delta)\right|\left|\partial_{\sigma}\mathbf{a}(\sigma_{-})\Delta\sigma+\delta\big{(}\mathbf{c}_{+}(\sigma_{-})-\mathbf{c}_{-}(\sigma_{-})\big{)}+\Delta\sigma\mathcal{O}(\delta)\right|$
$\displaystyle=$
$\displaystyle\left|\left[\partial_{\sigma}\mathbf{a}(\sigma_{-})\right]^{2}(\Delta\sigma)^{2}+\partial_{\sigma}\mathbf{a}(\sigma_{-})\big{(}\mathbf{c}_{+}(\sigma_{-})-\mathbf{c}_{-}(\sigma_{-})\big{)}\Delta\sigma\delta+\Delta\sigma\mathcal{O}(\delta^{2}|\log\delta|^{2})\right|.$

We note that $|\partial_{\sigma}\mathbf{a}(\sigma)|=16+40\sigma>16$ and,
therefore, we can write $|\hat{\lambda}_{0}(\sigma)|^{2}$ as

$|\hat{\lambda}_{0}(\sigma)|^{2}=\left[\partial_{\sigma}\mathbf{a}(\sigma_{-})\right]^{2}\Lambda$

where $\Lambda$, using formula (6.10) for $\mathbf{c}_{\pm}$, formula (3.18)
of $C_{0}$ and that $\sigma_{-}=\frac{m}{2}+\mathcal{O}(\delta)$, is given by

$\Lambda=\left|(\Delta\sigma)^{2}-2\delta\Delta\sigma
C_{0}+\Delta\sigma\mathcal{O}(\delta^{2}|\log\delta|^{2})\right|=|\Delta\sigma|\left|\Delta\sigma-2\delta
C_{0}+\mathcal{O}(\delta^{2}|\log\delta|^{2})\right|.$

Since $|\sigma-\sigma_{-}|\leq|\sigma-\sigma_{+}|$ and
$\sigma_{\pm}=\frac{m}{2}\pm
C_{0}\delta+\mathcal{O}(\delta^{2}|\log\delta|^{2})$, $\Delta\sigma<2\delta
C_{0}$ and then, using that $|\Delta\sigma|\geq C\delta^{\nu}$ with $C<C_{0}$,
for $\delta$ small enough

$\Lambda=|\Delta\sigma|(2\delta
C_{0}-\Delta\sigma++\mathcal{O}(\delta^{2}|\log\delta|^{2}))\geq
C\delta^{\nu}(2\delta
C_{0}-C\delta^{\nu}+\mathcal{O}(\delta^{2}|\log\delta|^{2}))\geq
M\delta^{1+\nu}$

for some constant $M>0$. Therefore, since $1+\nu\leq 3$,
$\mathbf{L}=\widehat{\mathcal{T}}|\hat{\lambda}_{0}|\geq
M\delta^{\frac{1+\nu}{2}}\geq M\delta^{3/2}$ and by Lemma 6.6,

$|\mu-\mu_{0}|\leq C_{*}\delta^{\frac{5-\nu}{2}}.$ (6.19)

An analogous analysis can be done if
$C\delta^{\nu}\leq|\sigma-\sigma_{+}|\leq|\sigma-\sigma_{-}|\leq\kappa\delta|\log\delta|$
and, as a consequence, (6.19) holds true for
$C\delta^{\nu}\leq|\sigma-\sigma_{\pm}|\leq\kappa\delta|\log\delta|$.

##  Appendix A Expression of the Hamiltonian

This appendix is devoted to obtain the expressions for the coplanar and
averaged Hamiltonians introduced in Section 2.2.

##### Slow-fast Delaunay coordinates $(y,x)$.

Analogously to the definition of the coplanar Hamiltonian
$\mathcal{H}_{\mathrm{CP},1}$ in Poincaré coordinates (see (2.17)), we define

${\mathrm{H}}_{\mathrm{CP},1}(y,\Gamma,x,h)={\mathrm{H}}_{1}(y,\Gamma,x,h,\Omega_{\mathrm{M}};0).$

Recall that, by Remark 2.1, ${\mathrm{H}}_{1}$ is independent of
$\Omega_{\mathrm{M}}$ when $i_{\mathrm{M}}=0$. Then, by (2.11) and (2.5), one
has that

$\displaystyle{\mathrm{H}}_{\mathrm{CP},1}(y,\Gamma,x,h)$
$\displaystyle=\frac{\rho_{1}}{L^{2}}\sum_{m=0}^{2}\sum_{p=0}^{2}\hat{c}_{m}U_{2}^{m,0}D_{m,p}(y,\Gamma)\cos\Big{(}\psi_{m,p,0}(x,h)\Big{)},$

where $D_{m,p}={\tilde{D}}_{m,p}\circ\Upsilon^{\mathrm{Del}}$,
$\psi_{m,p,0}=\tilde{\psi}_{m,p,0}\circ\Upsilon^{\mathrm{Del}}$, with
${\tilde{D}}_{m,p}$ and $\tilde{\psi}_{m.p,0}$ defined in (2.6) and (2.8)
respectively, and

$\displaystyle\hat{c}_{0}=c_{0,s}=\frac{1}{2},\qquad\hat{c}_{1}=c_{1,s}=\frac{1}{3},\qquad\hat{c}_{2}=c_{2,s}=-\frac{1}{12}.$

Notice that these last definitions are possible since, for
$s,m\in\\{0,1,2\\}$, the constants $c_{m,s}$ as defined in (2.9) do not depend
on $s$. In addition, $U_{2}^{m,0}(\epsilon)$ is the Giacaglia function given
in Table 4.

$m$ | $U_{2}^{m,0}(\epsilon)$ | $\simeq$
---|---|---
0 | $1-6C^{2}+6C^{4}$ | 0.762646
1 | $-3CS^{-1}\left(2C^{4}-3C^{2}+1\right)$ | 0.547442
2 | $6C^{2}S^{-2}\left(C^{2}-1\right)^{2}$ | 0.237353
Table 4: The Giacaglia function $U_{2}^{m,0}(\epsilon)$ for the moon
perturbation (see [8]) where $C=\cos\frac{\epsilon}{2}$ and
$S=\sin\frac{\epsilon}{2}$ and its value for $\epsilon=23.44$˚.

Applying the slow-fast change of coordinates, one obtains that

$\displaystyle\psi_{m,p,0}(x,h)$ $\displaystyle=(1-p)x-(1-p-m)h.$

See Table 5 for the explicit expressions of $\psi_{m,p,0}(x,h)$ and
$D_{m,p}(y,\Gamma)$.

$m$ | $p$ | $D_{m,p}(y,\Gamma)$ | $\psi_{m,p,0}$
---|---|---|---
0 | 0 | $-\frac{15}{64}(Ly)^{-2}(y-\Gamma)(3y+\Gamma)\left(L-2y\right)\left(L+2y\right)$ | $x-h$
0 | 1 | $\phantom{-}\frac{1}{32}(Ly)^{-2}(y^{2}-6y\Gamma-3\Gamma^{2})\left(5L^{2}-12y^{2}\right)$ | $0$
0 | 2 | $-\frac{15}{64}(Ly)^{-2}(y-\Gamma)(3y+\Gamma)\left(L-2y\right)\left(L+2y\right)$ | $-x-h$
1 | 0 | $\phantom{-}\frac{15}{32}(Ly)^{-2}\sqrt{(y-\Gamma)(3y+\Gamma)}\left(3y+\Gamma\right)\left(L-2y\right)\left(L+2y\right)$ | $x$
1 | 1 | $-\frac{3}{16}(Ly)^{-2}\sqrt{(y-\Gamma)(3y+\Gamma)}(\Gamma+y)\left(5L^{2}-12y^{2}\right)$ | $h$
1 | 2 | $-\frac{15}{32}(Ly)^{-2}\sqrt{(y-\Gamma)(3y+\Gamma)}\left(y-\Gamma\right)\left(L-2y\right)\left(L+2y\right)$ | $-x+2h$
2 | 0 | $\phantom{-}\frac{15}{32}(Ly)^{-2}\left(3y+\Gamma\right)^{2}\left(L-2y\right)\left(L+2y\right)$ | $x+h$
2 | 1 | $\phantom{-}\frac{3}{16}(Ly)^{-2}(y-\Gamma)(3y+\Gamma)\left(5L^{2}-12y^{2}\right)$ | $2h$
2 | 2 | $\phantom{-}\frac{15}{32}(Ly)^{-2}\left(y-\Gamma\right)^{2}\left(L-2y\right)\left(L+2y\right)$ | $-x+3h$
Table 5: Computation of the functions $(D_{m,p})_{m,p\in\\{0,1,2\\}}$ and
$(\psi_{m,p,0})_{m,p\in\\{0,1,2\\}}$ for the prograde case.

Now, we consider the $h$-average of ${\mathrm{H}}_{\mathrm{CP},1}$,

${\mathrm{H}}_{\mathrm{AV},1}(y,\Gamma,x)=\frac{1}{2\pi}\int_{0}^{2\pi}{\mathrm{H}}_{\mathrm{CP},1}(y,\Gamma,x,h){\mathrm{d}}h$

The definition of $\psi_{m,p,0}(x,h)$ implies that the only terms independent
of $h$ are for the couples $(m,p)=(0,1)$ and $(m,p)=(1,0)$ (see Table 5). As a
consequence, we obtain the expression of the averaged Hamiltonian presented in
(2.20).

##### Poincaré coordinates $(\eta,\xi)$.

Taking into account that the Poincaré change of coordinates satisfies that

$\displaystyle y=\frac{2L-\xi^{2}-\eta^{2}}{4},\qquad(L-2y)\cos
x=\frac{\xi^{2}-\eta^{2}}{2},\qquad(L-2y)\sin x=\xi\eta,$

one has that

$\displaystyle\mathcal{H}_{\mathrm{CP},1}=$
$\displaystyle\,\frac{\rho_{1}}{L^{2}}\sum_{m=0}^{2}\hat{c}_{m}U_{2}^{m,0}\cdot$
(A.1)
$\displaystyle\Big{[}\mathcal{D}_{m,0}(\eta,\Gamma,\xi)\left(\frac{\xi^{2}-\eta^{2}}{2}\cos\big{(}(1-m)h\big{)}+\xi\eta\sin\big{(}(1-m)h\big{)}\right)$
$\displaystyle+\mathcal{D}_{m,1}(\eta,\Gamma,\xi)\left(8L^{2}+12L(\xi^{2}+\eta^{2})-3(\xi^{2}+\eta^{2})^{2}\right)\cos(mh)$
$\displaystyle+\mathcal{D}_{m,2}(\eta,\Gamma,\xi)\left(\frac{\xi^{2}-\eta^{2}}{2}\cos\big{(}(1+m)h\big{)}+\xi\eta\sin\big{(}(1+m)h\big{)}\right)\Big{]},$

where the functions $(\mathcal{D}_{m,p})_{m,p\in\\{0,1,2\\}}$ are given in
Table 6.

$\begin{array}[]{ccl}\hline\cr
m&p&\lx@intercol\hfil\mathcal{D}_{m,p}(\eta,\Gamma,\xi)\hfil\lx@intercol\\\
\hline\cr
0&0&-\frac{15}{128}{L^{-2}(2L-M)^{-2}}{(2L-M-4\Gamma)(6L-3M+4\Gamma)(4L-M)}\\\
\hline\cr
0&1&\phantom{-}\frac{1}{128}{L^{-2}(2L-M)^{-2}}((2L-M)^{2}-24(2L-M)\Gamma-48\Gamma^{2})\\\
\hline\cr
0&2&-\frac{15}{128}{L^{-2}(2L-M)^{-2}}{(2L-M-4\Gamma)(6L-3M+4\Gamma)(4L-M)}\\\
\hline\cr
1&0&\phantom{-}\frac{15}{64}{L^{-2}(2L-M)^{-2}}{\sqrt{(2L-M-4\Gamma)}(6L-3M+4\Gamma)^{3/2}(4L-M)}\\\
\hline\cr
1&1&-\frac{3}{64}{L^{-2}(2L-M)^{-2}}{\sqrt{(2L-M-4\Gamma)(6L-3M+4\Gamma)}(2L-M+4\Gamma)}\\\
\hline\cr
1&2&-\frac{15}{64}{L^{-2}(2L-M)^{-2}}{(2L-M-4\Gamma)^{3/2}\sqrt{6L-3M+4\Gamma}(4L-M)}\\\
\hline\cr
2&0&\phantom{-}\frac{15}{64}{L^{-2}(2L-M)^{-2}}{(6L-3M+4\Gamma)^{2}(4L-M)}\\\
\hline\cr
2&1&\phantom{-}\frac{3}{64}{L^{-2}(2L-M)^{-2}}{(2L-M-4\Gamma)(6L-3M+4\Gamma)}\\\
\hline\cr
2&2&\phantom{-}\frac{15}{64}{L^{-2}(2L-M)^{-2}}{(2L-M-4\Gamma)^{2}(4L-M)}\\\
\hline\cr\end{array}$ Table 6: Computation of the functions
$(\mathcal{D}_{m,p})_{m,p\in\\{0,1,2\\}}$ with $M:=\xi^{2}+\eta^{2}$.

Let us now consider the $h$-averaged Hamiltonian $\mathcal{H}_{\mathrm{AV},1}$
as defined in (2.18). Analogously to ${\mathrm{H}}_{\mathrm{AV},1}$, one sees
that the only terms independent of $h$ are for the couples $(m,p)=(0,1)$ and
$(1,0)$. As a consequence, one obtains the expression in (2.19).

## References

  * [1] E. M. Alessi, I. Baldomá, M. Giralt, and M. Guardia.  On the Arnold diffusion mechanism in Medium Earth Orbit.  Preprint, available at https://arxiv.org/abs/2403.01283, 2024.
  * [2] E. M. Alessi, A. Rossi, G. B. Valsecchi, L. Anselmo, C. Pardini, C. Colombo, H. G. Lewis, J. Daquin, F. Deleflie, M. Vasile, F. Zuiani, and K. Merz.  Effectiveness of GNSS disposal strategies.  Acta Astronautica, 99:292–302, 2014.
  * [3] S. Breiter.  Lunisolar resonances revisited.  Celest. Mec. Dyn. Astron., 81:81–91, 2001.
  * [4] C. Chao.  Applied Orbit Perturbation and Maintenance.  American Institute of Aeronautics and Astronautics/Aerospace Press, Chao, C., ‘Applied Orbit Perturbation and Maintenance’, American Institute of Aeronautics and Astronautics/Aerospace Press, Reston, Virginia/El Segundo, California, 2005.
  * [5] G. Cook.  Luni-solar perturbations of the orbit of an earth satellite.  J. Roy. Astron. Soc., 6:271–291, 1962.
  * [6] J. Daquin, E. Legnaro, I. Gkolias, and C. Efthymiopoulos.  A deep dive into the $2g+h$ resonance: separatrices, manifolds and phase space structure of navigation satellites.  Celestial Mechanics and Dynamical Astronomy, 134:6, 2022.
  * [7] J. Daquin, A. J. Rosengren, E. M. Alessi, F. Deleflie, G. B. Valsecchi, and A. Rossi.  The dynamical structure of the MEO region: long-term evolution stability, chaos, and transport.  Celestial Mechanics and Dynamical Astronomy, 124:335–366, 2016\.
  * [8] G. E. O. Giacaglia.  Lunar perturbations of artificial satellites of the earth.  Celestial Mechanics and Dynamical Astronomy, 9:239–267, 1974.
  * [9] I. Gkolias, J. Daquin, D. K. Skoulidou, K. Tsiganis, and C. Efthymiopoulos.  Chaotic transport of navigation satellites.  Chaos, 29:101106, 2019.
  * [10] S. Hughes.  Earth satellite orbits with resonant lunisolar perturbations. i. resonances dependent only on inclination.  Proceedings of the Royal Society of London. Series A, Mathematical and Physical Sciences, 372:243–264, 1980.
  * [11] S. Hughes.  Earth satellite orbits with resonant lunisolar perturbations i. resonances dependent only on inclination.  Proceedings of the Royal Society of London. A. Mathematical and Physical Sciences, 372(1749):243–264, 1980.
  * [12] A. B. Jenkin and R. A. Gick.  Dilution of disposal orbit collision for the Medium Earth Orbit constellation.  In Proceedings of the 4th European Conference on Space Debris, ESA SP-587, pages 309–314. ESA/ESOC, 2005.
  * [13] W. Kaula.  Theory of Satellite Geodesy: Applications of Satellites to Geodesy.  Blaisdell Publi. CO., 1966.
  * [14] W. M. Kaula.  Development of the lunar and solar disturbing functions for a close satellite.  The Astronomical Journal, 67:300–303, 1962.
  * [15] E. Legnaro and C. Efthymiopoulos.  A detailed dynamical model for inclination-only dependent lunisolar resonances. Effect on the “eccentricity growth” mechanism.  Advances in Space Research, 72:2460–2480, 2023.
  * [16] E. Legnaro, C. Efthymiopoulos, and M. Harsoula.  Semi-analytical estimates for the chaotic diffusion in the Second Fundamental Model of Resonance. Application to Earth’s navigation satellites.  Physica D, 456:133946, 2023.
  * [17] P. Musen.  On the long-period lunar and solar effects on the motion of an artificial satellite, 2.  J. Geophys. Res., 66:2797–2805, 1961.
  * [18] M. Pellegrino, D. Scheeres, and B. J. Streetman.  The feasibility of targeting chaotic regions in the GNSS regime.  The Journal of the Astronautical Sciences, 68:553–584, 2021.
  * [19] E. Perozzi, A. E. Roy, B. A. Steves, and G. B. Valsecchi.  Significant high number commensurabilities in the main lunar problem. I: The Saros as a near periodicity of the Moon’s orbit.  Celestial Mechanics and Dynamical Astronomy, 52:241–261, 1991.
  * [20] J. Radtke, R. Domínguez-González, S. K. Flegel, N. Sánchez-Ortiz, and K. Merz.  Impact of eccentricity build-up and graveyard disposal strategies on MEO navigation constellations.  Advances in Space Research, 56:2626–2644, 2015.
  * [21] A. J. Rosengren, E. M. Alessi, A. Rossi, and G. B. Valsecchi.  Chaos in navigation satellite orbits caused by the perturbed motion of the Moon.  Monthly Notices of the Royal Astronomical Society, 449:3522–3526, 2015.
  * [22] A. Rossi.  Resonant dynamics of medium earth orbits: space debris issues.  Celestial Mechanics and Dynamical Astronomy, 100:267–286, 2008\.
  * [23] A. E. Roy.  The use of the saros in lunar dynamical studies.  The Moon, 7:6–13, 1973.
  * [24] V. G. Szebehely.  Theory of orbits : the restricted problem of three bodies.  Academic Press, New York [etc.], 1967.


# Stochastic foundations of $g$-subdiffusion process

Tadeusz Kosztołowicz<EMAIL_ADDRESS>Institute of Physics,
Jan Kochanowski University,
Uniwersytecka 7, 25-406 Kielce, Poland  Aldona Dutkiewicz<EMAIL_ADDRESS>Faculty of Mathematics and Computer Science,
Adam Mickiewicz University, Uniwersytetu Poznańskiego 4, 61-614 Poznań, Poland

###### Abstract

Recently, in the paper: T. Kosztołowicz and A. Dutkiewicz, Phys. Rev. E 104,
014118 (2021) the $g$–subdiffusion equation with fractional Caputo time
derivative with respect to another function $g$ has been considered. This
equation offers new possibilities for modelling diffusion such as a process in
which a type of diffusion evolves continuously over time. However, the
equation has not been derived from a stochastic model and the stochastic
interpretation of $g$–subdiffusion has been unknown. In this paper we show
stochastic foundations of this process. We derive the equation by means of a
modified Continuous Time Random Walk model. Interpretation of the
$g$–subdiffusion process is also discussed.

Introduction. Subdiffusion occurs in media in which the movement of diffusing
molecules is very difficult due to the complex internal structure of the
medium. Within the Continuous Time Random Walk (CTRW) model, a distribution of
time between particle jumps $\psi$ has a heavy tail for subdiffusion,
$\psi(t)\sim 1/t^{1+\alpha}$, $0<\alpha<1$ mk ; ks ; barkai2000 . This model
leads to the “ordinary” subdiffusion equation with the fractional order Caputo
derivative. Recently, a more general subdiffusion equation with the Caputo
derivative with respect to another function $g$ has been considered kd , see
also Ref. sz ; we call it the $g$–subdiffusion equation which describes the
$g$–subdiffusion process. As shown in Ref. kd , this equation describes a
process in which a type of diffusion can change over time. Unfortunately,
$g$–subdiffusion has not had a stochastic interpretation yet. We show how to
derive the $g$–subdiffusion equation by means of modified CTRW model and we
discuss the interpretation of this process.

“Ordinary” subdiffusion equation. The fractional subdiffusion equation with
“ordinary” Caputo derivative of the order $\alpha\in(0,1)$ is kd

$\frac{{}^{C}\partial^{\alpha}P(x,t)}{\partial
t^{\alpha}}=D\frac{\partial^{2}P(x,t)}{\partial x^{2}},$ (1)

where the Caputo fractional derivative is defined for $0<\alpha<1$ as

$\frac{{}^{C}d^{\alpha}}{dt^{\alpha}}f(t)=\frac{1}{\Gamma(1-\alpha)}\int_{0}^{t}(t-u)^{-\alpha}f^{\prime}(u)du,$
(2)

$\alpha$ is a subdiffusion parameter and $D$ is a generalized diffusion
coefficient measured in the units of ${\rm m^{2}/s^{\alpha}}$. To solve the
equation the Laplace transform $\mathcal{L}$ can be used,

$\mathcal{L}[f(t)](s)=\int_{0}^{\infty}{\rm e}^{-st}f(t)dt.$ (3)

Due to the relation

$\mathcal{L}\left[\frac{{}^{C}d^{\alpha}f(t)}{dt^{\alpha}}\right](s)=s^{\alpha}\mathcal{L}[f(t)](s)-s^{\alpha-1}f(0),$
(4)

where $0<\alpha\leq 1$, we get

$\displaystyle s^{\alpha}\mathcal{L}[P(x,t)](s)-s^{\alpha-1}P(x,0)$ (5)
$\displaystyle=D\frac{\partial^{2}\mathcal{L}[P(x,t)](s)}{\partial x^{2}}.$

$G$–subdiffusion equation. In this paper functions describing $g$–subdiffusion
are denoted by tilda. The $g$-subdiffusion equation reads

$\frac{{}^{C}\partial^{\alpha}_{g}\tilde{P}(x,t)}{\partial
t^{\alpha}}=D\frac{\partial^{2}\tilde{P}(x,t)}{\partial x^{2}},$ (6)

where $0<\alpha<1$, the Caputo derivative with respect to another function $g$
is defined as abd

$\frac{{}^{C}d^{\alpha}_{g}f(t)}{dt^{\alpha}}=\frac{1}{\Gamma(1-\alpha)}\int_{0}^{t}(g(t)-g(u))^{-\alpha}f^{\prime}(u)du,$
(7)

the function $g$ fulfils the conditions $g(0)=0$, $g(\infty)=\infty$, and
$g^{\prime}(t)>0$ for $t>0$, its values are given in a time unit. When
$g(t)=t$, the $g$-Caputo fractional derivative takes a form of the “ordinary”
Caputo derivative. To solve Eq. (6) the $g$–Laplace transform can be used,
this transform is defined as jarad

$\mathcal{L}_{g}[\tilde{f}(t)](s)=\int_{0}^{\infty}{\rm
e}^{-sg(t)}\tilde{f}(t)g^{\prime}(t)dt.$ (8)

Due to the property jarad

$\mathcal{L}_{g}\left[\frac{{}^{C}d^{\alpha}_{g}}{dt^{\alpha}}\tilde{f}(t)\right](s)=s^{\alpha}\mathcal{L}_{g}[\tilde{f}(t)](s)-s^{\alpha-1}\tilde{f}(0),$
(9)

the procedure of solving Eq. (6) is similar to the procedure of solving
“ordinary” subdiffusion equation by means of the “ordinary” Laplace transform
method. In terms of the $g$-Laplace transform the $g$-subdiffusion equation is

$\displaystyle
s^{\alpha}\mathcal{L}_{g}[\tilde{P}(x,t)](s)-s^{\alpha-1}\tilde{P}(x,0)$ (10)
$\displaystyle=D\frac{\partial^{2}\mathcal{L}_{g}[P(x,t)](s)}{\partial
x^{2}}.$

Using the $g$-Laplace transform to Eq. (6) yields Eq. (10) in the same form as
Eq. (5).

Model of particle random walk. To derive the subdiffusion equation we use a
simple model of a particle random walk along a one–dimensional homogeneous
lattice. Usually, in the CTRW model both a particle jump length and waiting
time for a particle jump are random variables. In our considerations, we
assume that the jump length distribution $\lambda$ has the form
$\lambda(x)=\frac{1}{2}[\delta(x-\epsilon)+\delta(x+\epsilon)]$, where
$\delta$ is the delta Dirac function. Only the choice of a particle jump
direction is random, its length $\epsilon$ is a parameter. We start with the
particle random walk model in which the particle positions and time are
discrete. Next, we move to continuous variables. Random walk with discrete
time $n$ is described by the equation
$P_{n+1}(m)=\frac{1}{2}P_{n}(m+1)+\frac{1}{2}P_{n}(m-1)$, where $P_{n}(m)$ is
a probability that a diffusing particle is at the position $m$ after $n$-th
step. Let the initial particle position be $m=0$. Moving from discrete $m$ to
continuous $x$ spatial variable we assume $x=m\epsilon$ and
$P_{n}(x)=P_{n}(m)/\epsilon$, where $\epsilon$ is a distance between discrete
sites. The above equations and the relation
$[P_{n}(x+\epsilon)+P_{n}(x-\epsilon)-2P_{n}(x)]/\epsilon^{2}=\partial^{2}P_{n}(x)/\partial
x^{2}$, $\epsilon\rightarrow 0$, provide the following equation in the limit
of small $\epsilon$

$P_{n+1}(x)-P_{n}(x)=\epsilon^{2}\frac{\partial^{2}P_{n}(x)}{\partial x^{2}}.$
(11)

To move from discrete to continuous time we use the formula mw

$P(x,t)=\sum_{n=0}^{\infty}Q_{n}(t)P_{n}(x),$ (12)

where $Q_{n}(t)$ is the probability that a diffusing particle takes $n$ step
in the time interval $(0,t)$. The function $Q_{n}$ is determined differently
for the “ordinary” subdiffusion and $g$–subdiffusion. In the following, we
find the rule for determining the functions $Q_{n}$ and the explicit form of
the functions $\psi$ for both processes. These functions, together with Eqs.
(11) and (12), provide “ordinary” subdiffusion and $g$–subdiffusion equations.

The case of “ordinary” subdiffusion. In this case the function $Q_{n}$ is a
convolution of $n$ distributions $\psi$ of a waiting time for a particle to
jump and a function $U$ which is the probability that a particle does not
change its position after $n$-th step,

$Q_{n}(t)=(\underbrace{\psi\ast\psi\ast\ldots\ast\psi}_{n\;times}\ast U)(t),$
(13)

where the convolution is defined as

$(f\ast h)(t)=\int_{0}^{t}f(u)h(t-u)du.$ (14)

The “ordinary” Laplace transform has the following property that makes the
transform useful in determining the function $Q_{n}$

$\mathcal{L}[(f\ast h)(t)](s)=\mathcal{L}[f(t)](s)\mathcal{L}[h(t)](s).$ (15)

From Eqs. (12), (13), and (15) we have

$\mathcal{L}[P(x,t)](s)=\mathcal{L}[U(t)](s)\sum_{n=0}^{\infty}\mathcal{L}^{n}[\psi(t)](s)P_{n}(x).$
(16)

Combining Eqs. (11), (12), and (16) we get

$\displaystyle\frac{2(1-\mathcal{L}[\psi(t)](s))}{\epsilon^{2}\mathcal{L}[\psi(t)](s)}\mathcal{L}[P(x,t)](s)$
(17)
$\displaystyle-\frac{2\mathcal{L}[U(t)](s)}{\epsilon^{2}\mathcal{L}[\psi(t)](s)}P(x,0)=\frac{\partial^{2}\mathcal{L}[P(x,t)](s)}{\partial
x^{2}}.$

Eq. (17) coincides with Eq. (5) only if

$\displaystyle\frac{1-\mathcal{L}[\psi(t)](s)]}{\mathcal{L}[\psi(t)](s)}=\frac{\epsilon^{2}s^{\alpha}}{2D},\;\frac{\mathcal{L}[U(t)](s)}{\mathcal{L}[\psi(t)](s)}=\frac{\epsilon^{2}s^{\alpha-1}}{2D}.$

The solutions to the above equations are

$\mathcal{L}[\psi(t)](s)=\frac{1}{1+\frac{\epsilon^{2}s^{\alpha}}{2D}},$ (18)

and

$\mathcal{L}[U(t)](s)=\frac{\epsilon^{2}s^{\alpha-1}}{2D\left(1+\frac{\epsilon^{2}s^{\alpha}}{2D}\right)}=\frac{1-\mathcal{L}[\psi(t)](s)}{s}.$
(19)

Due to the relations

$\mathcal{L}[1](s)=\frac{1}{s}\;,\;\mathcal{L}\left[\int_{0}^{t}f(u)du\right](s)=\frac{\mathcal{L}[f(t)](s)}{s},$
(20)

we get

$U(t)=1-\int_{0}^{t}\psi(u)du.$ (21)

In order to find the function $\psi$ we use the relation

$\displaystyle\mathcal{L}^{-1}[s^{\nu}{\rm
e}^{-as^{\beta}}](t)=\frac{1}{t^{1+\nu}}\sum_{k=0}^{\infty}\frac{1}{k!\Gamma(-\nu-\beta
k)}\left(-\frac{a}{t^{\beta}}\right)^{k}$ (22) $\displaystyle\equiv
f_{\nu,\beta}(t:a),$

where $a,\beta>0$, $\Gamma$ is the Gamma-Euler function. The function
$f_{\nu,\beta}$ is the Wright function and the special case of the H-Fox
function. To find the inverse Laplace transform of Eq. (18) first we calculate
the inverse Laplace transform of the function ${\rm e}^{-as^{\beta}}/(1+\tau
s^{\alpha})$, where $\tau=\epsilon^{2}/2D$ and $a,\beta>0$, using the formula
$1/(1+u)=\sum_{n=0}^{\infty}u^{n}$ when $|u|<1$. We get

$\displaystyle\mathcal{L}\left[\frac{{\rm e}^{-as^{\beta}}}{1+\tau
s^{\alpha}}\right](s)=$
$\displaystyle\left\\{\begin{array}[]{c}\frac{1}{\tau}\sum_{n=0}^{\infty}\left(-\frac{1}{\tau}\right)^{n}s^{-(n+1)\alpha}{\rm
e}^{-as^{\beta}},\;s>\frac{1}{\tau^{1/\alpha}},\\\ \\\
\sum_{n=0}^{\infty}(-\tau)^{n}s^{n\alpha}{\rm
e}^{-as^{\beta}},\;s<\frac{1}{\tau^{1/\alpha}},\end{array}\right.$ (26)

Next, we take the limit of $a\rightarrow 0$. From Eqs. (22), (Stochastic
foundations of $g$-subdiffusion process), and the relations
$f_{\nu,\beta}(t;0)=1/\Gamma(-\nu)t^{1+\nu}$, $1/\Gamma(0)=0$ we get

$\displaystyle\psi(t)=\left\\{\begin{array}[]{c}\frac{1}{\tau}\sum_{n=0}^{\infty}\left(-\frac{1}{\tau}\right)^{n}\frac{t^{(n+1)\alpha-1}}{\Gamma((n+1)\alpha)},\;t<\tau^{1/\alpha},\\\
\\\
\sum_{n=0}^{\infty}\left(-\tau\right)^{n+1}\frac{t^{-(n+1)\alpha-1}}{\Gamma(-(n+1)\alpha)},\;t>\tau^{1/\alpha},\end{array}\right.$
(30)

We have $\psi(t)\approx\alpha\tau/\Gamma(1-\alpha)t^{1+\alpha}$ in the limit
of $t\rightarrow\infty$. The function $\psi$ was already derived using the
relation $\mathcal{L}^{-1}[1/(1+\tau
s^{\alpha})]=t^{\alpha-1}E_{\alpha,\alpha}(-t^{\alpha}/\tau)$, where
$E_{\alpha,\alpha}(z)=\sum_{n=0}^{\infty}z^{n}/\Gamma(\alpha(n+1))$ is the
two–parameter Mittag–Leffler function, see for example Ref. m . Then, the
function $\psi$ corresponds to Eq. (30) but for the case of
$t<\tau^{1/\alpha}$ only.

The case of $g$–subdiffusion. To get Eq. (10) we use the $g$–Laplace
transform. This transform has the following property jarad

$\mathcal{L}_{g}[(f\ast_{g}h)(t)](s)=\mathcal{L}_{g}[f(t)](s)\mathcal{L}_{g}[h(t)](s),$
(31)

where the $g$–convolution is defined as

$(f\ast_{g}h)(t)=\int_{0}^{t}f(u)h(g^{-1}(g(t)-g(u))g^{\prime}(u)du.$ (32)

We involve the $g$–convolution in the CTRW model. Then, the procedure for
deriving the $g$–subdiffusion equation using the $g$–Laplace transform is
analogous to the procedure for deriving the ”ordinary” subdiffusion equation
using the “ordinary” Laplace transform. Assuming

$\tilde{P}(x,t)=\sum_{n=0}^{\infty}\tilde{Q}_{n}(t)P_{n}(x).$ (33)

and

$\tilde{Q}_{n}(t)=(\underbrace{\tilde{\psi}\ast_{g}\tilde{\psi}\ast_{g}\ldots\ast_{g}\tilde{\psi}}_{n\;times}\ast_{g}\tilde{U})(t),$
(34)

from Eqs. (31), (33), and (34) we obtain

$\mathcal{L}_{g}[\tilde{P}(x,t)](s)=\sum_{n=0}^{\infty}\mathcal{L}_{g}[\tilde{U}(t)](s)\mathcal{L}^{n}_{g}[\tilde{\psi}(t)](s)P_{n}(x).$
(35)

From Eqs. (11) and (35) we get

$\displaystyle\frac{1-\mathcal{L}_{g}[\tilde{\psi}(t)](s)}{\epsilon^{2}\mathcal{L}_{g}[\tilde{\psi}(t)](s)}\mathcal{L}_{g}[\tilde{P}(x,t)](s)$
(36)
$\displaystyle-\frac{\mathcal{L}_{g}[\tilde{U}(t)](s)}{\epsilon^{2}\mathcal{L}_{g}[\tilde{\psi}(t)](s)}\tilde{P}(x,0)=\frac{\partial^{2}\mathcal{L}[\tilde{P}(x,t)](s)}{\partial
x^{2}}.$

Eq. (36) is consistent with Eq. (10) only when

$\mathcal{L}_{g}[\tilde{\psi}(t)](s)=\frac{1}{1+\frac{\epsilon^{2}s^{\alpha}}{2D}}$
(37)

and

$\mathcal{L}_{g}[\tilde{U}(t)](s)=\frac{\epsilon^{2}s^{\alpha-1}}{2D\left(1+\frac{\epsilon^{2}s^{\alpha}}{2D}\right)}.$
(38)

Comparing Eqs. (37) and (38) with Eqs. (18) and (19), respectively, we get

$\mathcal{L}_{g}[\tilde{\psi}(t)](s)=\mathcal{L}[\psi(t)](s),$ (39)
$\mathcal{L}_{g}[\tilde{U}(t)](s)=\mathcal{L}[U(t)](s).$ (40)

From the relation

$\mathcal{L}_{g}[\tilde{f}(t)](s)=\mathcal{L}[\tilde{f}(g^{-1}(t))](s),$ (41)

we get the following rule kd

$\mathcal{L}_{g}[\tilde{f}(t)](s)=\mathcal{L}[f(t)](s)\Leftrightarrow\tilde{f}(t)=f(g(t)).$
(42)

Due to Eq. (42), from Eqs. (39) and (40) we obtain

$\tilde{\psi}(t)=\psi(g(t)),$ (43)

and

$\tilde{U}(t)=U(g(t)).$ (44)

Eqs. (30) and (43) provide

$\displaystyle\tilde{\psi}(t)=\left\\{\begin{array}[]{c}\frac{1}{\tau}\sum_{n=0}^{\infty}\left(-\frac{1}{\tau}\right)^{n}\frac{g^{(n+1)\alpha-1}(t)}{\Gamma((n+1)\alpha)},\;t<g^{-1}(\tau^{1/\alpha}),\\\
\\\
\sum_{n=0}^{\infty}\left(-\tau\right)^{n+1}\frac{g^{-(n+1)\alpha-1}(t)}{\Gamma(-(n+1)\alpha)},\;t>g^{-1}(\tau^{1/\alpha}),\end{array}\right.$
(48)

We get $\psi(t)\approx\alpha\tau/\Gamma(1-\alpha)g^{1+\alpha}(t)$ when
$t\rightarrow\infty$.

We link the $g$–convolution with the “ordinary” convolution. Let
$\tilde{f}(t)=f(g(t))$ and $\tilde{h}(t)=h(g(t))$. After simple calculation we
get

$(\tilde{f}\ast_{g}\tilde{h})(t)=(f\ast h)(g(t)).$ (49)

From Eqs. (33), (34), and (49) we have

$\tilde{P}(x,t)=\sum_{n=0}^{\infty}Q_{n}(g(t))P_{n}(x).$ (50)

Comparing Eqs. (50) and (12) we obtain

$\tilde{P}(x,t)=P(x,g(t)).$ (51)

Interpretation. The $g$–subdiffusion process is associated to “ordinary”
subdiffusion controlled by the same parameter $\alpha$. The waiting time for a
particle jump in the $g$–subdiffusion process is controlled by the functions
$\psi$ and $g$. A particle jump that would occur with some probability after
time $t$ in ”ordinary” subdiffusion process will occur with the same
probability after time $\tilde{t}=g^{-1}(t)$ in the $g$–subdiffusion process.
If $g(t)<t$ we have $t<\tilde{t}$, subdiffusion is then slowed down. When
$g(t)>t$ subdiffusion is accelerated. An example of $g$–subdiffusion is the
diffusion of molecules in a medium consisting of matrix in which there are
narrow channels. If channels have a complicated geometric structure and
diffusing molecules do not interact with the matrix, then ”ordinary”
subdiffusion controlled by the parameter $\alpha$ occurs. When temporary
penetration of a molecule into the matrix is possible then the molecule
”disappears” from channels and may diffuse further upon returning to a
channel. In this case, ”ordinary” subdiffusion is slowed down. Such a process
occurs in a vessel filled with alginate beads immersed in water in which
colistin antibiotic diffuses kdlwa . When the matrix provides the diffusing
molecules with additional energy, subdiffusion can be accelerated.

Final remarks. We have shown that the $g$–subdiffusion equation can be derived
by means of the modified CTRW model (we call it the $g$–CTRW model). In the
$g$-CTRW model we use the $g$-convolution and the $g$-Laplace transform
instead of the “ordinary” convolution and the“ordinary” Laplace transform,
respectively, which are used in the ”ordinary” CTRW model.

We note that the condition $\mathcal{L}_{g}[\tilde{\psi}(t)](0)=1$ does not
guarantee that the function $\tilde{\psi}$ is normalized. Therefore,
$\tilde{\psi}$ is not a probability distribution. Thus, it seems that the
$g$-CTRW model is merely a mathematical procedure. However, this model can be
interpreted as an ”ordinary” CTRW model in which the time scale is controlled
by the function $g(t)$, see Eqs. (43)–(51). The key issue for the
$g$–subdiffusion process is determining the parameter $\alpha$ and the
function $g$. An example of their determination from empirical data is shown
in kdlwa .

In practice, the transformations made in deriving the $g$–subdiffusion
equation within the $g$-CTRW model are the same as in deriving the ”ordinary”
subdiffusion equation using ”ordinary” CTRW. Within the ”ordinary” CTRW,
subdiffusion-reaction equations kl  as well as the Green’s functions and
membrane boundary conditions for a system in which a thin membrane separates
different subdiffusive media tk2019  have been derived. Within the $g$–CTRW
model the same procedures can also be used to derive $g$–subdiffusion-reaction
equations, Green’s functions, and boundary conditions at the membrane for the
processes described by $g$–subdiffusion equations.

We suppose that the $g$–subdiffusion model can be used to describe diffusion
of antibiotics in a biofilm. Biofilm usually has a gel structure. When the
antibiotic does not interact with bacteria, the “ordinary” antibiotic
subdiffusion in the biofilm is expected. However, bacteria in the biofilm have
different defense mechanisms against the action of the antibiotic. These
mechanisms may hinder or even facilitate antibiotic subdiffusion, see Ref. km
and the references cited therein. Thus, the application of the
$g$–subdiffusion equation to describe this process may be effective.

## References

  * (1) R. Metzler and J. Klafter, Phys. Rep. 339, 1 (2000).
  * (2) J. Klafter and I.M. Sokolov, First step in random walks. From tools to applications, (Oxford UP, New York, 2011).
  * (3) E. Barkai, R. Metzler, and J. Klafter, Phys. Rev. E 61, 132 (2000).
  * (4) T. Kosztołowicz and A. Dutkiewicz, Phys. Rev. E 104, 014118 (2021).
  * (5) B. Samet and Y. Zhou, RACSAM 113, 2887 (2019); R. Garra, A. Giusti, and F. Mainardi, Ricerche Mat. 67, 899 (2018).
  * (6) W. Abdelhedi, Comp. Appl. Math. 40, 53 (2021); R. Almeida, Commun. Nonlinear Sci. Numer. Simul. 44, 460 (2017); A.A. Kilbas, H.M. Srivastava, and J.J. Trujillo, Theory and Applications of Fractional Differential Equations (North-Holland Mathematics Studies, 204, Elsevier, Amsterdam, 2006); J.V.D.C. Sousa and E.C. de Oliveira, Commun. Nonlinear Sci. Numer. Simul. 60, 72 (2018).
  * (7) F. Jarad and T. Abdeljawad, Discrete And Continuous Dynamical Systems-S 13, 709 (2020); F. Jarad, T. Abdeljawad, S. Rashid, and Z. Hammouch, Adv. Differ. Equ. 2020, 303 (2020).
  * (8) E.W. Montroll and G.H. Weiss, J. Math. Phys. 6, 167 (1965).
  * (9) T. Kosztołowicz, J. Phys. A 37, 10779 (2004).
  * (10) F. Mainardi, WSEAS Trans. Math. 19, 74 (2020); Entropy 22, 1359 (2020).
  * (11) T. Kosztołowicz, A. Dutkiewicz, K. Lewandowska, S. Wa̧sik, and M. Arabski, arXiv: cond-mat. 2107.02419 (2021).
  * (12) T. Kosztołowicz and K. Lewandowska, Phys. Rev. E 90, 032136 (2014).
  * (13) T. Kosztołowicz, Phys. Rev. E 99, 022127 (2019); Int. J. Heat Mass Transf. 111, 1322 (2017).
  * (14) T. Kosztołowicz and R. Metzler, Phys. Rev. E 102, 032408 (2020); T. Kosztołowicz, R. Metzler, S. Wa̧sik, and M. Arabski, PLoS One 15, e0243003 (2020).