$\left\|v-v^{\prime}\right\|_{g^{-1}}\leq\delta_{v}$, then

$\left|-H(x^{\prime},v^{\prime})+H(x,v)\right|\leq\Theta\left(\left|f(x)-f(x^{\prime})\right|+\left(\delta_{x}\left\|v\right\|_{g^{-1}}^{2}+\delta_{v}^{2}+\delta_{v}\left\|v\right\|_{g^{-1}}\right)+n\delta_{x}\right),$

where the Hamiltonian is
$H(x,v)=f(x)+\frac{1}{2}v^{\top}g(x)^{-1}v+\frac{1}{2}\log\det g(x)$.

###### Proof.

We consider each term separately. For the second term,

$\displaystyle\left|\frac{1}{2}v^{\top}g^{-1}v-\frac{1}{2}v^{\prime\top}g^{\prime-1}v\right|$
$\displaystyle\leq\frac{1}{2}\underbrace{\left|v^{\top}g^{-1}v-v^{\top}g^{\prime-1}v\right|}_{F}+\frac{1}{2}\underbrace{\left|v^{\top}g^{\prime-1}v-v^{\prime\top}g^{\prime-1}v^{\prime}\right|}_{S}.$

For $F$, we have $F\leq O(\delta_{x})\left\|v\right\|_{g^{-1}}^{2}$ by Lemma
26-3. For $S$, it follows that

$\displaystyle S$
$\displaystyle=\left|\left\|v\right\|_{g^{\prime-1}}^{2}-\left\|v^{\prime}\right\|_{g^{\prime-1}}^{2}\right|\leq\left\|v-v^{\prime}\right\|_{g^{\prime-1}}\left\|v+v^{\prime}\right\|_{g^{\prime-1}}$
$\displaystyle\leq
O(1)\left\|v-v^{\prime}\right\|_{g^{-1}}(\left\|v\right\|_{g^{-1}}+\left\|v^{\prime}\right\|_{g^{-1}})$
$\displaystyle\leq O(1)\delta_{v}(\delta_{v}+\left\|v\right\|_{g^{-1}}).$

Therefore, the second term is bounded by
$O(1)\left(\delta_{x}\left\|v\right\|_{g^{-1}}^{2}+\delta_{v}^{2}+\delta_{v}\left\|v\right\|_{g^{-1}}\right)$.

For the third term,

$\left|\frac{1}{2}\left(\log\det g(x)-\log\det
g(x^{\prime})\right)\right|=\frac{1}{2}\left|\log\det
g^{\prime-\frac{1}{2}}gg^{\prime-\frac{1}{2}}\right|\leq O(n\delta_{x}),$

where the inequality follows from Lemma 26 and the fact that the determinant
is the product of eigenvalues. ∎

Since the ideal RHMC preserves the Hamiltonian along its Hamiltonian curve,
$H(x,v)=H(x_{h},v_{h})$. Hence, we can obtain a lower bound on the acceptance
probability by computing either

$\min\left(1,\frac{e^{-H(\bar{x}_{h},\bar{v}_{h})}}{e^{-H(x,v)}}\right)\
\text{or}\
\min\left(1,\frac{e^{-H(\bar{x}_{h},\bar{v}_{h})}}{e^{-H(x_{h},v_{h})}}\right).$

######  Lemma 32.

Let $(x_{h},v_{h})$ and $(\bar{x},\bar{v})$ be the points obtained by the
ideal RHMC and discretized RHMC with a sensitive numerical integrator starting
at $x\in\mathbb{\mathcal{M}}_{\rho}$ with $v\in V_{\text{good}}^{x}$. If the
step size $h$ satisfies

$h^{2}\leq\frac{10^{-10}}{n+\sqrt{M_{1}}+M_{2}},\>h^{2}C_{x}(x,v)\leq\frac{10^{-10}}{n+\sqrt{M_{1}}+\sqrt{M_{1}^{*}}},\,h^{2}C_{v}(x,v)\leq\frac{10^{-10}}{\sqrt{n+\sqrt{M_{1}}}},$

then the rejection probability of the Metropolis filter is bounded by
$10^{-3}$.

###### Proof.

We use the first condition on the step size to obtain
$\left\|v_{h}\right\|_{g^{-1}}=O\left(h\left(n+\sqrt{M_{1}}\right)\right)=O\left(\sqrt{n+\sqrt{M_{1}}}\right)$
by Proposition 28-1. Then the claims follows from

$\displaystyle
d_{\text{TV}}(\overline{\mathcal{P}}_{x}^{\prime},\overline{\mathcal{P}}_{x})$
$\displaystyle\leq
10^{5}\left(\delta_{x}\sqrt{M_{1}^{*}}+\delta_{x}\left(n+\left\|v_{h}\right\|_{g^{-1}}^{2}\right)+\delta_{v}\left(\delta_{v}+\left\|v_{h}\right\|_{g^{-1}}\right)\right)$
$\displaystyle\leq
10^{6}h^{2}\left(C_{x}(x,v)\left(n+\sqrt{M_{1}}+\sqrt{M_{1}^{*}}\right)+C_{v}(x,v)\left(C_{v}(x,v)h^{2}+\sqrt{n+\sqrt{M_{1}}}\right)\right)$
$\displaystyle\leq 10^{-4}+10^{-20}+10^{-4}\leq 10^{-3},$

where we used the second and third step-size conditions in the last
inequality. ∎

Putting three main parts together, we obtain the result on the mixing rate of
RHMC discretized by a sensitive numerical integrator.

###### Proof of Theorem 24.

By Lemma 14, 29 and 32, we have
$d_{\text{TV}}(\overline{\mathcal{P}}_{x},\overline{\mathcal{P}}_{y})\leq\frac{9}{10}$
if $d_{\phi}(x,y)\leq h$ for $x,y\in\mathbb{\mathcal{M}}_{\rho}$. Then the
claim follows by reproducing the proof of Proposition 15. ∎

##  7 Numerical integrators

We examine two numerical integrators commonly used in practice, the implicit
midpoint method (IMM) in Section 7.1 and the generalized Leapfrog method (LM)
in Section 7.2. To this end, we bound $C_{x}$ and $C_{v}$, the second-
orderness parameters, and then find a condition on step size for the
sensitivity. We note that these integrators are symplectic (so measure-
preserving) and time-reversible (see [HHIL06]).

###  7.1 Implicit midpoint method

For an initial condition $(x,v)$ and step size $h$, the implicit midpoint
method attempts to find the solution $(x^{\prime},v^{\prime})$ for the
following implicit equation:

$x^{\prime}=x+h\frac{\partial H}{\partial
v}\left(\frac{x+x^{\prime}}{2},\frac{v+v^{\prime}}{2}\right),\
v^{\prime}=v-h\frac{\partial H}{\partial
x}\left(\frac{x+x^{\prime}}{2},\frac{v+v^{\prime}}{2}\right).$

In general, these implicit equations require several iterations so that an
initial guess for this equation converges to the fixed point
$(x^{\prime},v^{\prime})$.

In this section, we consider a variant of IMM in Algorithm 2 instead. It has
computational benefits over the original IMM, since iterations for finding the
fixed point of the integrator run with a simpler Hamiltonian
$H_{2}(x,v)=\frac{1}{2}v^{\top}g(x)^{-1}v$ instead of $H=H_{1}+H_{2}$. We then
prove that if for $x\in\mathbb{\mathcal{M}}_{\rho}$ and $v\in
V_{\text{good}}^{x}$ step size $h$ satisfies
$h^{2}\left(n+\sqrt{M_{1}}\right)\leq 10^{-10}$, then IMM is second-order with
$C_{x}(x,v)=O\left(n+\sqrt{M_{1}}\right)$ and
$C_{v}(x,v)=O\left(\sqrt{n+\sqrt{M_{1}}}\left(n+\sqrt{M_{1}}+M_{2}^{*}\right)\right)$.
Moreover, if the step size $h$ satisfies
$h^{2}\leq\min\left(\frac{10^{-10}}{\left(n+\sqrt{M_{1}}\right)^{2}},\frac{10^{-5}}{\sqrt{n}\bar{R}_{1}}\right)$
in addition to the step-size conditions in Proposition 28, then IMM is
sensitive at $x\in\mathbb{\mathcal{M}}_{\rho}$ and $v\in V_{\text{good}}^{x}$.

Input: Initial point $x$, velocity $v$, step size $h$

// Step 1: Solve $\frac{dx}{dt}=\frac{\partial H_{1}(x,v)}{\partial
v},\frac{dv}{dt}=-\frac{\partial H_{1}(x,v)}{\partial x}$

Set $x_{\frac{1}{3}}\leftarrow x$ and $v_{\frac{1}{3}}\leftarrow
v-\frac{h}{2}\frac{\partial H_{1}(x,v)}{\partial x}$.

// Step 2: Solve $\frac{dx}{dt}=\frac{\partial H_{2}(x,v)}{\partial
v},\frac{dv}{dt}=-\frac{\partial H_{2}(x,v)}{\partial x}$ (Implicit)

Find $(x_{{\frac{2}{3}}},v_{{\frac{2}{3}}})$ such that

$\displaystyle x_{{\frac{2}{3}}}$
$\displaystyle=x_{\frac{1}{3}}+h\frac{\partial H_{2}}{\partial
v}\left(\frac{x_{\frac{1}{3}}+x_{{\frac{2}{3}}}}{2},\frac{v_{\frac{1}{3}}+v_{{\frac{2}{3}}}}{2}\right),$
$\displaystyle v_{{\frac{2}{3}}}$
$\displaystyle=v_{\frac{1}{3}}-h\frac{\partial H_{2}}{\partial
x}\left(\frac{x_{\frac{1}{3}}+x_{{\frac{2}{3}}}}{2},\frac{v_{\frac{1}{3}}+v_{{\frac{2}{3}}}}{2}\right).$

// Step 3: Solve $\frac{dx}{dt}=\frac{\partial H_{1}(x,v)}{\partial
v},\frac{dv}{dt}=-\frac{\partial H_{1}(x,v)}{\partial x}$

Set $x_{1}\leftarrow x_{\frac{2}{3}}$ and $v_{1}\leftarrow
v_{\frac{2}{3}}-\frac{h}{2}\frac{\partial H_{1}}{\partial
x}\left(x_{\frac{2}{3}},v_{\frac{2}{3}}\right)$.

Output: $x_{1,}v_{1}$

Algorithm 2 Implicit Midpoint Method

####  7.1.1 Second-order

######  Lemma 33.

For $x\in\mathbb{\mathcal{M}}_{\rho}$ and $v\in V_{\text{good}}^{x}$, let
$g=g(x)$ and $h$ step size of IMM with $h^{2}\left(n+\sqrt{M_{1}}\right)\leq
10^{-10}$. Let $(\overline{x},\overline{v})$ be the point obtained from RHMC
discretized by IMM with the step size $h$ and initial condition $(x,v)$.

  1. 1.

$\left\|x-\overline{x}\right\|_{g}=O\left(h\sqrt{n}+h^{2}\left(n+\sqrt{M_{1}}\right)\right)$
and $\left\|v-\overline{v}\right\|_{g^{-1}}=O\left(h\left(\left\|\nabla
f(\bar{x})\right\|_{g^{-1}}+n+\sqrt{M_{1}}\right)\right)$.

  2. 2.

$C_{x}(x,v)=O\left(n+\sqrt{M_{1}}\right)$.

  3. 3.

$C_{v}(x,v)=O\left(\sqrt{n+\sqrt{M_{1}}}\left(n+\sqrt{M_{1}}+M_{2}^{*}\right)\right)$
.

##### Proof of 1.

Let $\bar{x}=x_{{\frac{2}{3}}}=\overline{T}_{x}(v)$ and
$v_{\frac{1}{3}},v_{\frac{2}{3}},\bar{v}$ be the velocity points obtained when
starting with $(x,v)$. Then $x_{\mathrm{mid}}$ and $v_{\mathrm{mid}}$ satisfy

$\displaystyle x_{{\frac{2}{3}}}$
$\displaystyle=x_{\frac{1}{3}}+hg_{\mathrm{mid}}^{-1}v_{\mathrm{mid}},$ (7.1)
$\displaystyle v_{{\frac{2}{3}}}$
$\displaystyle=v_{\frac{1}{3}}+\frac{h}{2}Dg_{\mathrm{mid}}\left[g_{\mathrm{mid}}^{-1}v_{\mathrm{mid}},g_{\mathrm{mid}}^{-1}v_{\mathrm{mid}}\right],$

where
$x_{\mathrm{mid}}=\frac{x_{\frac{1}{3}}+x_{{\frac{2}{3}}}}{2},v_{\mathrm{mid}}=\frac{v_{\frac{1}{3}}+v_{{\frac{2}{3}}}}{2}$
and $g_{\mathrm{mid}}=g\left(x_{\mathrm{mid}}\right)$. Since
$\left\|\bar{x}-x\right\|_{g_{\mathrm{mid}}}=\left\|x_{{\frac{2}{3}}}-x\right\|_{g_{\mathrm{mid}}}\to
0$ as $h\to 0$, we can take $h_{0}>0$ such that
$\left\|x_{{\frac{2}{3}}}-x\right\|_{g_{\mathrm{mid}}}\leq\frac{1}{1000}$ for
$h\leq h_{0}$ with the equality held at $h=h_{0}$, or
$\left\|x_{\frac{2}{3}}-x\right\|_{g_{\mathrm{mid}}}\leq\frac{1}{1000}$ for
any $h>0$.

We start with the former. By adding $v_{\frac{1}{3}}$ to the second line of
(7.1) and dividing by $2$, we have

$h_{0}v_{\mathrm{mid}}=h_{0}v_{\frac{1}{3}}+\frac{h_{0}^{2}}{4}Dg_{\mathrm{mid}}\left[g_{\mathrm{mid}}^{-1}v_{\mathrm{mid}},g_{\mathrm{mid}}^{-1}v_{\mathrm{mid}}\right].$
(7.2)

As
$\left\|x_{{\frac{2}{3}}}-x\right\|_{g_{\mathrm{mid}}}=h_{0}\left\|v_{\mathrm{mid}}\right\|_{g_{\mathrm{mid}}^{-1}}$
from the first line of (7.1), taking the $g_{\mathrm{mid}}^{-1}$-norm on both
sides of (7.2) and using Lemma 27 yield

$\displaystyle\frac{1}{1000}$
$\displaystyle=h_{0}\left\|v_{\mathrm{mid}}\right\|_{g_{\mathrm{mid}}^{-1}}$
$\displaystyle\leq
h_{0}\left\|v_{\frac{1}{3}}\right\|_{g_{\mathrm{mid}}^{-1}}+\frac{h_{0}^{2}}{2}\left\|v_{\mathrm{mid}}\right\|_{g_{\mathrm{mid}}^{-1}}^{2}\leq
h_{0}\left\|v_{\frac{1}{3}}\right\|_{g_{\mathrm{mid}}^{-1}}+\frac{1}{2000},$

and we obtain $\frac{1}{2000}\leq
h_{0}\left\|v_{\frac{1}{3}}\right\|_{g_{\mathrm{mid}}^{-1}}$. Recall that
$\left\|\bar{x}-x\right\|_{g_{\mathrm{mid}}}\leq 1/1000$ for $h\leq h_{0}$, so
we can swap the local norms between $x_{\mathrm{mid}}$ and $x$ due to Lemma
26, losing a multiplicative constant like $1.001$. Using Lemma 27 on the first
step of the numerical integrator,

$\displaystyle\left\|v_{\frac{1}{3}}\right\|_{g_{\mathrm{mid}}^{-1}}$
$\displaystyle\leq 1.001\left\|v_{\frac{1}{3}}\right\|_{g^{-1}}\leq
1.001\left(\left\|v\right\|_{g^{-1}}+h_{0}\left(n+\sqrt{M_{1}}\right)\right)$
(7.3) $\displaystyle\leq 200\sqrt{n}+2h_{0}\left(n+\sqrt{M_{1}}\right).$

Due to $1/2000\leq
h_{0}\left\|v_{\frac{1}{3}}\right\|_{g_{\mathrm{mid}}^{-1}}$, it follows that

$\displaystyle\frac{1}{2000}$ $\displaystyle\leq
h_{0}\left\|v_{\frac{1}{3}}\right\|_{g_{\mathrm{mid}}^{-1}}\leq
200\sqrt{n}h_{0}+2h_{0}^{2}\left(n+\sqrt{M_{1}}\right),$

and solving this for $h_{0}$ we have
$h_{0}\geq\frac{1}{10^{4}\sqrt{n+\sqrt{M_{1}}}}$. For the case of
$\left\|x_{\frac{2}{3}}-x\right\|_{g_{\mathrm{mid}}}\leq\frac{1}{1000}$ for
any $h>0$, we can simply think of $h_{0}$ as $\infty$.

Now for $h\leq h_{0}$, we can obtain from (7.1)

$\displaystyle\left\|x_{{\frac{2}{3}}}-x\right\|_{g_{\mathrm{mid}}}$
$\displaystyle\leq h\left\|v_{\mathrm{mid}}\right\|_{g_{\mathrm{mid}}^{-1}},$
(7.4)
$\displaystyle\left\|v_{{\frac{2}{3}}}-v_{\frac{1}{3}}\right\|_{g_{\mathrm{mid}}^{-1}}$
$\displaystyle\leq
h\left\|v_{\mathrm{mid}}\right\|_{g_{\mathrm{mid}}^{-1}}^{2}.$

Using this and
$h\left\|v_{\mathrm{mid}}\right\|_{g_{\mathrm{mid}}^{-1}}=\left\|x_{{\frac{2}{3}}}-x\right\|_{g_{\mathrm{mid}}}\leq\frac{1}{1000}$,
we can bound $\left\|v_{\mathrm{mid}}\right\|_{g_{\mathrm{mid}}^{-1}}$ by

$\displaystyle\left\|v_{\mathrm{mid}}\right\|_{g_{\mathrm{mid}}^{-1}}$
$\displaystyle=\left\|v_{\frac{1}{3}}+\frac{v_{{\frac{2}{3}}}-v_{\frac{1}{3}}}{2}\right\|_{g_{\mathrm{mid}}^{-1}}\leq\left\|v_{\frac{1}{3}}\right\|_{g_{\mathrm{mid}}^{-1}}+\frac{1}{2}h\left\|v_{\mathrm{mid}}\right\|_{g_{\mathrm{mid}}^{-1}}^{2}$
$\displaystyle\leq\left\|v_{\frac{1}{3}}\right\|_{g_{\mathrm{mid}}^{-1}}+\frac{1}{2000}\left\|v_{\mathrm{mid}}\right\|_{g_{\mathrm{mid}}^{-1}},$

and thus

$\displaystyle\left\|v_{\mathrm{mid}}\right\|_{g_{\mathrm{mid}}^{-1}}\leq\frac{2000}{1999}\left\|v_{\frac{1}{3}}\right\|_{g_{\mathrm{mid}}^{-1}}$
$\displaystyle\underset{\eqref{eq:LMIMM-
ineq1}}{\leq}200\sqrt{n}+2h\left(n+\sqrt{M_{1}}\right).$ (7.5)

Putting this back to (7.4) for step size
$h\leq\frac{1}{10^{5}\sqrt{n+\sqrt{M_{1}}}}$, we have

$\displaystyle\left\|x_{{\frac{2}{3}}}-x\right\|_{g_{\mathrm{mid}}}$
$\displaystyle\leq 200h\sqrt{n}+2h^{2}\left(n+\sqrt{M_{1}}\right),$
$\displaystyle\left\|v_{{\frac{2}{3}}}-v_{\frac{1}{3}}\right\|_{g_{\mathrm{mid}}^{-1}}$
$\displaystyle\leq 125h\left(n+\sqrt{M_{1}}\right).$

Hence by substituting the step size into above and switching local norms
properly, we have $\left\|x_{{\frac{2}{3}}}-x\right\|_{g}\leq 10^{-8}$.

By applying Lemma 27 to
$\overline{v}=v_{{\frac{2}{3}}}-\frac{h}{2}\frac{\partial H_{1}}{\partial
x}\left(x_{{\frac{2}{3}}},v_{{\frac{2}{3}}}\right)$ in the third step, we also
have

$\displaystyle\left\|\bar{v}-v\right\|_{g^{-1}}$
$\displaystyle\leq\left\|\bar{v}-v_{\frac{2}{3}}\right\|_{g^{-1}}+\left\|v_{\frac{2}{3}}-v_{\frac{1}{3}}\right\|_{g^{-1}}+\left\|v_{\frac{1}{3}}-v\right\|_{g^{-1}}$
$\displaystyle\leq 1.001h\left(\left\|\nabla
f(x_{{\frac{2}{3}}})\right\|_{g^{-1}}+n+125\left(n+\sqrt{M_{1}}\right)+\left(n+\sqrt{M_{1}}\right)\right)$
$\displaystyle\leq 200h\left(\left\|\nabla
f(x_{{\frac{2}{3}}})\right\|_{g^{-1}}+n+\sqrt{M_{1}}\right).$

In conclusion,

$\displaystyle\left\|\bar{x}-x\right\|_{g}$ $\displaystyle\leq
200h\sqrt{n}+3h^{2}\left(n+\sqrt{M_{1}}\right),$
$\displaystyle\left\|\bar{v}-v\right\|_{g^{-1}}$ $\displaystyle\leq
200h\left(\left\|\nabla f(\bar{x})\right\|_{g^{-1}}+n+\sqrt{M_{1}}\right).$

##### Proof of 2.

For $t\in[0,h]$, let $(x_{t},v_{t})$ be the Hamiltonian curve of the ideal
RHMC at time $t$ starting from $(x,v)$. Recall that for $g_{t}=g(x_{t})$

$\displaystyle T_{x}(v)$ $\displaystyle=x+\int_{0}^{h}\frac{\partial
H}{\partial v}(x_{t},v_{t})dt=x+\int_{0}^{h}g_{t}^{-1}v_{t}dt,$
$\displaystyle\overline{T}_{x}(v)$
$\displaystyle=\overline{x}=x+hg_{\mathrm{mid}}^{-1}v_{\mathrm{mid}}.$

Thus,

$\displaystyle\left\|T_{x}(v)-\overline{T}_{x}(v)\right\|_{g_{\mathrm{mid}}}$
$\displaystyle=\left\|\left(x+\int_{0}^{h}g_{t}^{-1}v_{t}dt\right)-\left(x+hg_{\mathrm{mid}}^{-1}v_{\mathrm{mid}}\right)\right\|_{g_{\mathrm{mid}}}$
$\displaystyle=\left\|\int_{0}^{h}\left(g_{t}^{-1}v_{t}-g_{\mathrm{mid}}^{-1}v_{\mathrm{mid}}\right)dt\right\|_{g_{\mathrm{mid}}}$
$\displaystyle\leq
h\max_{t\in[0,h]}\left\|g_{t}^{-1}v_{t}-g_{\mathrm{mid}}^{-1}v_{\mathrm{mid}}\right\|_{g_{\mathrm{mid}}}$

By Lemma 26-11,

$\displaystyle\left\|g_{t}^{-1}v_{t}-g_{\mathrm{mid}}^{-1}v_{\mathrm{mid}}\right\|_{g_{\mathrm{mid}}}$
$\displaystyle\lesssim\left\|v_{t}-v_{\mathrm{mid}}\right\|_{g_{\mathrm{mid}}^{-1}}+\left\|x_{t}-x_{\mathrm{mid}}\right\|_{g_{\mathrm{mid}}}\left\|v_{t}\right\|_{g_{\mathrm{mid}}^{-1}}$
$\displaystyle\leq\frac{1}{2}\left(\left\|v_{t}-v_{\frac{1}{3}}\right\|_{g_{\mathrm{mid}}^{-1}}+\left\|v_{t}-v_{{\frac{2}{3}}}\right\|_{g_{\mathrm{mid}}^{-1}}\right)$
$\displaystyle\quad+\frac{1}{2}\left(\left\|x_{t}-x\right\|_{g_{\mathrm{mid}}}+\left\|x_{t}-\overline{x}\right\|_{g_{\mathrm{mid}}}\right)\left\|v_{t}\right\|_{g_{\mathrm{mid}}^{-1}}$
$\displaystyle\lesssim\left(\left\|v-v_{t}\right\|_{g^{-1}}+\left\|v-v_{\frac{1}{3}}\right\|_{g^{-1}}+\left\|v_{t}-v\right\|_{g^{-1}}+\left\|v-v_{{\frac{2}{3}}}\right\|_{g^{-1}}\right)$
$\displaystyle\quad+\left(\left\|x_{t}-x\right\|_{g}+\left\|x_{t}-x\right\|_{g}+\left\|\overline{x}-x\right\|_{g}\right)\left(\left\|v_{t}-v\right\|_{g^{-1}}+\left\|v\right\|_{g^{-1}}\right)$
$\displaystyle\lesssim\left(\left\|v-v_{t}\right\|_{g^{-1}}+\left\|v-v_{\frac{1}{3}}\right\|_{g^{-1}}+\left\|v_{\frac{1}{3}}-v_{{\frac{2}{3}}}\right\|_{g^{-1}}\right)$
$\displaystyle\quad+\left(\left\|x_{t}-x\right\|_{g}+\left\|\overline{x}-x\right\|_{g}\right)\left(\left\|v_{t}-v\right\|_{g^{-1}}+\left\|v\right\|_{g^{-1}}\right).$

Using our bounds on $\left\|\bar{x}-x\right\|_{g},\left\|x_{t}-x\right\|_{g}$
and
$\left\|v\right\|_{g^{-1}},\left\|v_{t}-v\right\|_{g^{-1}},\left\|v-v_{\frac{1}{3}}\right\|_{g^{-1}},\left\|v_{\frac{1}{3}}-v_{{\frac{2}{3}}}\right\|_{g^{-1}}$,
we conclude that
$\max_{t\in[0,h]}\left\|g_{t}^{-1}v_{t}-g_{\mathrm{mid}}^{-1}v_{\mathrm{mid}}\right\|_{g_{\mathrm{mid}}}\leq
10^{4}h\left(n+\sqrt{M_{1}}\right)$, and thus

$\left\|T_{x}(v)-\overline{T}_{x}(v)\right\|_{g}\leq
10^{4}h^{2}\left(n+\sqrt{M_{1}}\right).$

##### Proof of 3.

From the algorithm,

$\displaystyle v_{h}$ $\displaystyle=v-\int_{0}^{h}\frac{\partial H}{\partial
x}(x_{t},v_{t})dt$
$\displaystyle=v-\int_{0}^{h}\frac{\mathcal{\partial}H_{1}}{\mathcal{\partial}x}(x_{t},v_{t})dt-\int_{0}^{h}\frac{\mathcal{\partial}H_{2}}{\mathcal{\partial}x}(x_{t},v_{t})dt,$
$\displaystyle\bar{v}$
$\displaystyle=v_{{\frac{2}{3}}}-\frac{h}{2}\frac{\mathcal{\partial}H_{1}}{\mathcal{\partial}x}(x_{{\frac{2}{3}}},v_{{\frac{2}{3}}})=v_{\frac{1}{3}}-h\frac{\mathcal{\partial}H_{2}}{\mathcal{\partial}x}(x_{\mathrm{mid}},v_{\mathrm{mid}})-\frac{h}{2}\frac{\mathcal{\partial}H_{1}}{\mathcal{\partial}x}(x_{{\frac{2}{3}}},v_{{\frac{2}{3}}})$
$\displaystyle=v-\frac{h}{2}\left(\frac{\mathcal{\partial}H_{1}}{\mathcal{\partial}x}(x,v)+\frac{\mathcal{\partial}H_{1}}{\mathcal{\partial}x}(x_{{\frac{2}{3}}},v_{{\frac{2}{3}}})\right)-h\frac{\mathcal{\partial}H_{2}}{\mathcal{\partial}x}(x_{\mathrm{mid}},v_{\mathrm{mid}}).$

Thus,

$\displaystyle\left\|v_{h}-\bar{v}\right\|_{g_{\mathrm{mid}}^{-1}}$ (7.6)
$\displaystyle=\left\|\int_{0}^{h}\left(\frac{\mathcal{\partial}H_{1}}{\mathcal{\partial}x}(x_{t},v_{t})-\frac{1}{2}\left(\frac{\mathcal{\partial}H_{1}}{\mathcal{\partial}x}(x,v)+\frac{\mathcal{\partial}H_{1}}{\mathcal{\partial}x}(x_{{\frac{2}{3}}},v_{{\frac{2}{3}}})\right)\right)dt+\int_{0}^{h}\left(\frac{\mathcal{\partial}H_{2}}{\mathcal{\partial}x}(x_{t},v_{t})-\frac{\mathcal{\partial}H_{2}}{\mathcal{\partial}x}(x_{\mathrm{mid}},v_{\mathrm{mid}})\right)dt\right\|_{g_{\mathrm{mid}}^{-1}}$
$\displaystyle\leq
h\underbrace{\max_{t\in[0,h]}\left\|\frac{\mathcal{\partial}H_{1}}{\mathcal{\partial}x}(x_{t},v_{t})-\frac{1}{2}\left(\frac{\mathcal{\partial}H_{1}}{\mathcal{\partial}x}(x,v)+\frac{\mathcal{\partial}H_{1}}{\mathcal{\partial}x}(x_{{\frac{2}{3}}},v_{{\frac{2}{3}}})\right)\right\|_{g_{\mathrm{mid}}^{-1}}}_{F}$
$\displaystyle\quad+h\underbrace{\max_{t\in[0,h]}\left\|\frac{\mathcal{\partial}H_{2}}{\mathcal{\partial}x}(x_{t},v_{t})-\frac{\mathcal{\partial}H_{2}}{\mathcal{\partial}x}(x_{\mathrm{mid}},v_{\mathrm{mid}})\right\|_{g_{\mathrm{mid}}^{-1}}}_{S}.$

We separately bound $F$ and $S$. For $F$, by following the proof of
Proposition 26-12, we have

$\displaystyle F$ $\displaystyle\lesssim
n\left(\left\|v_{t}-v\right\|_{g^{-1}}+\left\|v_{t}-v_{{\frac{2}{3}}}\right\|_{g^{-1}}\right)+M_{2}^{*}\left(\left\|x_{h}-x\right\|_{g}+\left\|x_{h}-\bar{x}\right\|_{g}\right).$

Using our bounds on
$\left\|x_{h}-\bar{x}\right\|_{g},\left\|x_{h}-x\right\|_{g}$ and
$\left\|v_{t}-v\right\|_{g^{-1}},\left\|v-v_{{\frac{2}{3}}}\right\|_{g^{-1}}$,
we obtain

$F\lesssim
h\left(n+\sqrt{M_{1}}\right)^{3/2}+h\sqrt{n+\sqrt{M_{1}}}M_{2}^{*}.$

We remark that the smoothness of $f$ guarantees that $M_{2}^{*}$ is bounded by
some constant for all sufficiently small $h$.

Similarly for $S$, we have that for
$\delta_{v}=\left\|v_{t}-v_{\mathrm{mid}}\right\|_{g_{\mathrm{mid}}^{-1}}$ and
$\delta_{x}=\left\|x_{t}-x_{\mathrm{mid}}\right\|_{g_{\mathrm{mid}}}$

$S\lesssim\max_{t\in[0,h]}\left(\delta_{v}+\delta_{x}\left\|v_{\mathrm{mid}}\right\|_{g_{\mathrm{mid}}^{-1}}\right)\left(\left\|v_{\mathrm{mid}}\right\|_{g_{\mathrm{mid}}^{-1}}+\left\|v_{t}\right\|_{g_{\mathrm{mid}}^{-1}}\right).$

Using our bounds on $\left\|\bar{x}-x\right\|_{g},\left\|x_{t}-x\right\|_{g}$
and
$\left\|v\right\|_{g^{-1}},\left\|v_{t}-v\right\|_{g^{-1}},\left\|v-v_{\frac{1}{3}}\right\|_{g^{-1}},\left\|v_{\frac{1}{3}}-v_{{\frac{2}{3}}}\right\|_{g^{-1}}$,
it follows that

$S\lesssim h\left(n+\sqrt{M_{1}}\right)^{3/2}.$

Substituting the bounds on $F$ and $S$ into (7.6) we can conclude that

$\left\|v_{h}-\bar{v}\right\|_{g^{-1}}\leq
10^{10}h^{2}\sqrt{n+\sqrt{M_{1}}}\left(n+\sqrt{M_{1}}+M_{2}^{*}\right).$

####  7.1.2 Sensitivity

We use Proposition 30 to show that IMM is sensitive. Here we assume that
$\log\det g(x)$ is convex in $\mathbb{\mathcal{M}}$, which is the case for the
logarithmic barriers of polytopes.

######  Lemma 34.

For $x\in\mathbb{\mathcal{M}}_{\rho}$ and $v\in V_{\text{good}}^{x}$, if
$\log\det g(x)$ is convex in $x$, then IMM is sensitive at $(x,v)$ for step
size $h$ satisfying
$h^{2}\leq\min\left(\frac{10^{-10}}{\left(n+\sqrt{M_{1}}\right)^{2}},\frac{10^{-5}}{\sqrt{n}\bar{R}_{1}}\right)$
and the step-size conditions in Proposition 28.

###### Proof.

Let $\overline{F}(x,v)=(\bar{x},\bar{v})$ and $\overline{T}_{x}(v)=\bar{x}$.
We lower bound $\left|D\overline{T}_{x}(v)\right|.$ Recall that one iteration
of IMM consists of three steps with input $(x,v)$ and output $(x_{1},v_{1})$,
as described in the following diagram:

$\left(x,v\right)\overset{X}{\longrightarrow}\left(x_{\frac{1}{3}},v_{\frac{1}{3}}\right)\overset{Y}{\longrightarrow}\left(x_{{\frac{2}{3}}},v_{{\frac{2}{3}}}\right)\overset{Z}{\longrightarrow}\left(x_{1},v_{1}\right),$

where each of the maps $X,Y$ and $Z$ is defined by

$\displaystyle X\left(x,v\right)$
$\displaystyle=\left(x,v-\frac{h}{2}\frac{\partial H_{1}(x,v)}{\partial
x}\right),$ $\displaystyle Y\left(x_{\frac{1}{3}},v_{\frac{1}{3}}\right)$
$\displaystyle=\left(x_{\frac{1}{3}}+h\frac{\partial
H_{2}(x_{\mathrm{mid}},v_{\mathrm{mid}})}{\partial
v},v_{\frac{1}{3}}-h\frac{\partial
H_{2}(x_{\mathrm{mid}},v_{\mathrm{mid}})}{\partial x}\right),$ $\displaystyle
Z\left(x_{{\frac{2}{3}}},v_{{\frac{2}{3}}}\right)$
$\displaystyle=\left(x_{{\frac{2}{3}}},v_{{\frac{2}{3}}}-\frac{h}{2}\frac{\partial
H_{1}(x_{{\frac{2}{3}}},v_{{\frac{2}{3}}})}{\partial x}\right),$

for $x_{\mathrm{mid}}=\frac{x_{\frac{1}{3}}+x_{{\frac{2}{3}}}}{2}$ and
$v_{\mathrm{mid}}=\frac{v_{\frac{1}{3}}+v_{{\frac{2}{3}}}}{2}$. Due to
$\overline{T}_{x}(v)=\pi_{x}\circ(Z\circ Y\circ X)(x,v)$, it follows that
$D\overline{T}_{x}(v)$ is the upper-right $n\times n$ submatrix of $D(Z\circ
Y\circ X)(x,v)$. From direct computation, we have

$\displaystyle DX(x,v)$ $\displaystyle=\left[\begin{array}[]{cc}I&0\\\
*&I\end{array}\right],$ $\displaystyle DY(x_{\frac{1}{3}},v_{\frac{1}{3}})$
$\displaystyle=\left[\begin{array}[]{cc}P&Q\\\ R&S\end{array}\right],$
$\displaystyle DZ(x_{{\frac{2}{3}}},v_{{\frac{2}{3}}})$
$\displaystyle=\left[\begin{array}[]{cc}I&0\\\ *&I\end{array}\right],$

and due to $D(Z\circ Y\circ X)=DZ\cdot DY\cdot DX$ we have
$D\overline{T}_{x}(v)=Q$. Thus, it suffices to focus on the second step only
(i.e., the map $Y$).

Now let us represent the map $Y$ in a compact way. With two symbols

$r=\left[\begin{array}[]{c}x\\\ v\end{array}\right]\in\mathbb{R}^{2n}\
\text{and}\ J=\left[\begin{array}[]{cc}0&I_{n}\\\
-I_{n}&0\end{array}\right]\in\mathbb{R}^{2n\times 2n},$

the second step can be rewritten as

$r_{{\frac{2}{3}}}=r_{\frac{1}{3}}+hJ\nabla_{(x,v)}H_{2}(r_{\mathrm{mid}}),$

where $r_{*}=\left[\begin{array}[]{c}x_{*}\\\ v_{*}\end{array}\right]$ for
$*\in\left\\{\frac{1}{3},{\frac{2}{3},\mathrm{mid}}\right\\}$ and
$H_{2}(x,v)=\frac{1}{2}v^{\top}g(x)^{-1}v$. Differentiating both sides by
$r_{\frac{1}{3}}$,

$\frac{\partial
r_{{\frac{2}{3}}}}{\mathcal{\partial}r_{\frac{1}{3}}}=I_{2n}+hJ\nabla^{2}H_{2}(r_{\mathrm{mid}})\left(\frac{1}{2}I_{2n}+\frac{1}{2}\frac{\partial
r_{{\frac{2}{3}}}}{\mathcal{\partial}r_{\frac{1}{3}}}\right).$

As $DY\left(r_{\frac{1}{3}}\right)=\frac{\partial
r_{{\frac{2}{3}}}}{\mathcal{\partial}r_{\frac{1}{3}}}$, we have that

$\left(I_{2n}-\frac{h}{2}J\nabla^{2}H_{2}(r_{\mathrm{mid}})\right)DY\left(r_{\frac{1}{3}}\right)=I_{2n}+\frac{h}{2}J\nabla^{2}H_{2}(r_{\mathrm{mid}}).$

For $G(x)\stackrel{{\scriptstyle\mathrm{{\scriptscriptstyle
def}}}}{{=}}\left[\begin{array}[]{cc}g(x)^{\frac{1}{2}}&0\\\
0&g(x)^{-\frac{1}{2}}\end{array}\right]$, we have

$\displaystyle G(x_{\mathrm{mid}})$
$\displaystyle\left(I_{2n}-\frac{h}{2}J\nabla^{2}H_{2}(r_{\mathrm{mid}})\right)G(x_{\mathrm{mid}})^{-1}G(x_{\mathrm{mid}})DY\left(r_{\frac{1}{3}}\right)G(x_{\mathrm{mid}})^{-1}$
$\displaystyle=G(x_{\mathrm{mid}})\left(I_{2n}+\frac{h}{2}J\nabla^{2}H_{2}(r_{\mathrm{mid}})\right)G(x_{\mathrm{mid}})^{-1}$

and

$\displaystyle\left(I_{2n}-\frac{h}{2}G(x_{\mathrm{mid}})J\nabla^{2}H_{2}(r_{\mathrm{mid}})G(x_{\mathrm{mid}})^{-1}\right)G(x_{\mathrm{mid}})DY\left(r_{\frac{1}{3}}\right)G(x_{\mathrm{mid}})^{-1}$
$\displaystyle=I_{2n}+\frac{h}{2}G(x_{\mathrm{mid}})J\nabla^{2}H_{2}(r_{\mathrm{mid}})G(x_{\mathrm{mid}})^{-1}.$
(7.7)

Let us look into the term $B\stackrel{{\scriptstyle\mathrm{{\scriptscriptstyle
def}}}}{{=}}G(x_{\mathrm{mid}})J\nabla^{2}H_{2}(r_{\mathrm{mid}})G(x_{\mathrm{mid}})^{-1}$.
By direct computation, for block matrices $B_{1},B_{2},B_{3},B_{4}$ of size
$n\times n$ we have

$\displaystyle B$
$\displaystyle\stackrel{{\scriptstyle\mathrm{{\scriptscriptstyle
def}}}}{{=}}\left[\begin{array}[]{cc}B_{1}&B_{2}\\\
B_{3}&B_{4}\end{array}\right]$
$\displaystyle=\left[\begin{array}[]{cc}g(x_{\mathrm{mid}})^{\frac{1}{2}}\\\
&g(x_{\mathrm{mid}})^{-\frac{1}{2}}\end{array}\right]\left[\begin{array}[]{cc}\frac{\mathcal{\partial}^{2}H_{2}}{\mathcal{\partial}v\mathcal{\partial}x}(r_{\mathrm{mid}})&\frac{\mathcal{\partial}^{2}H_{2}}{\mathcal{\partial}v^{2}}(r_{\mathrm{mid}})\\\
-\frac{\mathcal{\partial}^{2}H_{2}}{\mathcal{\partial}x^{2}}(r_{\mathrm{mid}})&-\left(\frac{\mathcal{\partial}^{2}H_{2}}{\mathcal{\partial}x\mathcal{\partial}v}(r_{\mathrm{mid}})\right)^{\top}\end{array}\right]\left[\begin{array}[]{cc}g(x_{\mathrm{mid}})^{-\frac{1}{2}}\\\
&g(x_{\mathrm{mid}})^{\frac{1}{2}}\end{array}\right]$

and thus

$\displaystyle B_{1}$
$\displaystyle=g(x_{\mathrm{mid}})^{-\frac{1}{2}}Dg(x_{\mathrm{mid}})\left[g(x_{\mathrm{mid}})^{-1}v_{\mathrm{mid}}\right]g(x_{\mathrm{mid}})^{-\frac{1}{2}},$
(7.8) $\displaystyle B_{2}$ $\displaystyle=I_{n},$ $\displaystyle B_{3}$
$\displaystyle=g(x_{\mathrm{mid}})^{-\frac{1}{2}}\bigg{(}-v_{\mathrm{mid}}^{\top}g(x_{\mathrm{mid}})^{-1}Dg(x_{\mathrm{mid}})g(x_{\mathrm{mid}})^{-1}Dg(x_{\mathrm{mid}})g(x_{\mathrm{mid}})^{-1}v_{\mathrm{mid}}$
$\displaystyle\leavevmode\nobreak\ \leavevmode\nobreak\
+\frac{1}{2}v_{\mathrm{mid}}^{\top}g(x_{\mathrm{mid}})^{-1}D^{2}g(x_{\mathrm{mid}})g(x_{\mathrm{mid}})^{-1}v_{\mathrm{mid}}\bigg{)}g(x_{\mathrm{mid}})^{-\frac{1}{2}}$
$\displaystyle=-B_{1}^{2}+\frac{1}{2}g(x_{\mathrm{mid}})^{-\frac{1}{2}}D^{2}g(x_{\mathrm{mid}})\left[g(x_{\mathrm{mid}})^{-1}v_{\mathrm{mid}},g(x_{\mathrm{mid}})^{-1}v_{\mathrm{mid}}\right]g(x_{\mathrm{mid}})^{-\frac{1}{2}},$
$\displaystyle B_{4}$ $\displaystyle=-B_{1}^{\top}.$

Now we bound the operator norm of $B_{i}$ for each $i\in[4]$ as follows.

$\displaystyle\left\|B_{1}\right\|$ $\displaystyle=\left\|B_{4}\right\|$
$\displaystyle=\max_{p,q:\left\|p\right\|_{2},\left\|q\right\|_{2}\leq
1}p^{\top}g(x_{\mathrm{mid}})^{-\frac{1}{2}}Dg(x_{\mathrm{mid}})\left[g(x_{\mathrm{mid}})^{-1}v_{\mathrm{mid}}\right]g(x_{\mathrm{mid}})^{-\frac{1}{2}}q$
$\displaystyle=\max_{p,q}Dg(y_{\mathrm{mid}})\left[g(x_{\mathrm{mid}})^{-1}v_{\mathrm{mid}},g(x_{\mathrm{mid}})^{-\frac{1}{2}}p,g(x_{\mathrm{mid}})^{-\frac{1}{2}}q\right]$
$\displaystyle\leq
2\max_{p,q}\left\|g(x_{\mathrm{mid}})^{-1}v_{\mathrm{mid}}\right\|_{g(x_{\mathrm{mid}})}\left\|g(x_{\mathrm{mid}})^{-\frac{1}{2}}p\right\|_{g(x_{\mathrm{mid}})}\left\|g(x_{\mathrm{mid}})^{-\frac{1}{2}}q\right\|_{g(x_{\mathrm{mid}})}$
$\displaystyle=2\max_{p,q}\left\|g(x_{\mathrm{mid}})^{-1}v_{\mathrm{mid}}\right\|_{g(x_{\mathrm{mid}})}\left\|p\right\|_{2}\left\|q\right\|_{2}\leq
2\left\|v_{\mathrm{mid}}\right\|_{g_{\mathrm{mid}}^{-1}}$ $\displaystyle\leq
O\left(\sqrt{n+\sqrt{M_{1}}}\right),$

where we used (7.5) guaranteed by the condition of
$h^{2}\left(n+\sqrt{M_{1}}\right)\leq 10^{-10}$ ($\small{1}⃝$ in Proposition
28). For $B_{2}$ and $B_{3}$, we have

$\displaystyle\left\|B_{2}\right\|$ $\displaystyle=1,$
$\displaystyle\left\|B_{3}\right\|$
$\displaystyle\leq\left\|B_{1}\right\|^{2}+\frac{1}{2}\max_{p,q:\left\|p\right\|_{2},\left\|q\right\|_{2}\leq
1}D^{2}g(x_{\mathrm{mid}})\left[g(x_{\mathrm{mid}})^{-1}v_{\mathrm{mid}},g(x_{\mathrm{mid}})^{-1}v_{\mathrm{mid}},g(x_{\mathrm{mid}})^{-\frac{1}{2}}p,g(x_{\mathrm{mid}})^{-\frac{1}{2}}q\right]$
$\displaystyle\leq
O\left(n+\sqrt{M_{1}}\right)+3\max_{p,q}\left\|g(x_{\mathrm{mid}})^{-1}v_{\mathrm{mid}}\right\|_{g(x_{\mathrm{mid}})}^{2}\left\|g(x_{\mathrm{mid}})^{-\frac{1}{2}}p\right\|_{g(x_{\mathrm{mid}})}\left\|g(x_{\mathrm{mid}})^{-\frac{1}{2}}q\right\|_{g(x_{\mathrm{mid}})}$
$\displaystyle=O\left(n+\sqrt{M_{1}}\right)+O\left(n+\sqrt{M_{1}}\right)\max_{p,q}\left\|p\right\|_{2}\left\|q\right\|_{2}$
$\displaystyle=O\left(n+\sqrt{M_{1}}\right),$

where the second inequality for $\left\|B_{3}\right\|$ follows from the highly
self-concordance of $\phi$. Due to
$\left\|B\right\|\leq\sum_{i=1}^{4}\left\|B_{i}\right\|$, we have
$\left\|\frac{h}{2}B\right\|=O\left(h\left(n+\sqrt{M_{1}}\right)\right)$.
Hence, the condition of $h^{2}\left(n+\sqrt{M_{1}}\right)^{2}\leq 10^{-10}$
ensures that the inverse of $I_{2n}-\frac{h}{2}B$ exists, and it can be
written as a series of matrices,

$\left(I_{2n}-\frac{h}{2}B\right)^{-1}=\sum_{i=0}^{\infty}(hB/2)^{i}.$

By substituting this series into (7.7),

$\displaystyle
G(x_{\mathrm{mid}})DY\left(r_{\frac{1}{3}}\right)G(x_{\mathrm{mid}})^{-1}$
$\displaystyle=\sum_{i=0}^{\infty}(hB/2)^{i}\left(I_{2n}+\frac{h}{2}B\right)=\sum_{i=0}^{\infty}(hB/2)^{i}+\sum_{i=1}^{\infty}(hB/2)^{i}$
$\displaystyle=I_{2n}+2\sum_{i=1}^{\infty}(hB/2)^{i}.$

By multiplying $\left[\begin{array}[]{cc}I_{n}&0\end{array}\right]^{\top}$ to
the left and $\left[\begin{array}[]{c}0\\\ I_{n}\end{array}\right]$ to the
right on both sides,

$g(x_{\mathrm{mid}})^{\frac{1}{2}}D\overline{T}_{x}(v)g(x_{\mathrm{mid}})^{\frac{1}{2}}=2\sum_{i=1}^{\infty}\left[\begin{array}[]{cc}I_{n}&0\end{array}\right]^{\top}(hB/2)^{i}\left[\begin{array}[]{c}0\\\
I_{n}\end{array}\right].$

From (7.8), $B$ is of the form

$B=\left[\begin{array}[]{cc}C&I_{n}\\\ -C^{2}+R&-C\end{array}\right],$

where $C\in\mathbb{R}^{n\times n}$ is symmetric and $R\in\mathbb{R}^{n\times
n}$, and thus by Lemma 58

$g(x_{\mathrm{mid}})^{\frac{1}{2}}D\overline{T}_{x}(v)g(x_{\mathrm{mid}})^{\frac{1}{2}}=h\sum_{i=0}^{\infty}(h^{2}R/2)^{i},$

where
$R=\frac{1}{2}g(x_{\mathrm{mid}})^{-\frac{1}{2}}D^{2}g(x_{\mathrm{mid}})\left[g(x_{\mathrm{mid}})^{-1}v_{\mathrm{mid}},g(x_{\mathrm{mid}})^{-1}v_{\mathrm{mid}}\right]g(x_{\mathrm{mid}})^{-\frac{1}{2}}$.
Thus for $E\stackrel{{\scriptstyle\mathrm{{\scriptscriptstyle
def}}}}{{=}}\sum_{i=1}^{\infty}(h^{2}R/2)^{i}$

$\frac{1}{h}g(x_{\mathrm{mid}})^{\frac{1}{2}}D\overline{T}_{x}(v)g(x_{\mathrm{mid}})^{\frac{1}{2}}=I+E.$
(7.9)

We now bound its operator norm, trace and Frobenius norm. It is easy to see
that

$\displaystyle\left\|E\right\|_{2}$ $\displaystyle\lesssim\sum_{i\geq
1}\left(\frac{h^{2}}{2}\left(n+\sqrt{M_{1}}\right)\right)^{i},$
$\displaystyle\mathrm{Tr}(E)$ $\displaystyle\lesssim\sum_{i\geq
1}\left(\frac{h^{2}}{2}\mathrm{Tr}(R)\right)^{i}\leq\sum_{i\geq
1}\left(\frac{h^{2}}{2}n\left(n+\sqrt{M_{1}}\right)\right),$
$\displaystyle\left\|E\right\|_{F}$ $\displaystyle\lesssim\sum_{i\geq
1}\left(\frac{h^{2}}{2}\sqrt{n}\left(n+\sqrt{M}\right)\right)^{i},$

where we used the following estimations

$\displaystyle\left\|R\right\|_{2}$ $\displaystyle\leq
O\left(n+\sqrt{M_{1}}\right)$ $\displaystyle\mathrm{Tr}(R)$
$\displaystyle=\frac{1}{2}\mathbb{E}_{p\sim\mathcal{N}(0,I)}p^{\top}g(x_{\mathrm{mid}})^{-\frac{1}{2}}D^{2}g(x_{\mathrm{mid}})[g(x_{\mathrm{mid}})^{-1}v_{\mathrm{mid}},g^{-1}v_{\mathrm{mid}}]g(x_{\mathrm{mid}})^{-\frac{1}{2}}p$
$\displaystyle\leq\mathbb{E}D^{2}g(x_{\mathrm{mid}})\left[g(x_{\mathrm{mid}})^{-1}v_{\mathrm{mid}},g(x_{\mathrm{mid}})^{-1}v_{\mathrm{mid}},g(x_{\mathrm{mid}})^{-\frac{1}{2}}p,g(x_{\mathrm{mid}})^{-\frac{1}{2}}p\right]$
$\displaystyle\leq\mathbb{E}\left\|v_{\mathrm{mid}}\right\|_{g_{\mathrm{mid}}^{-1}}^{2}\left\|p\right\|_{2}^{2}=O\left(n\left(n+\sqrt{M_{1}}\right)\right),$
$\displaystyle\left\|R\right\|_{F}$
$\displaystyle\leq\sqrt{n}\left\|R\right\|_{2}=O\left(\sqrt{n}\left(n+\sqrt{M}\right)\right).$

Therefore, the step-size condition of
$h^{2}\left(n+\sqrt{M_{1}}\right)^{2}\leq 10^{-10}$ ensures that these three
quantities can be made smaller than $10^{-8}$. Applying Lemma 60 to (7.9), we
have

$\displaystyle
e^{\mathrm{Tr}(E)}e^{-\left\|E\right\|_{F}^{2}}\leq\left|\frac{1}{h}g(x_{\mathrm{mid}})^{\frac{1}{2}}D\overline{T}_{x}(v)g(x_{\mathrm{mid}})^{\frac{1}{2}}\right|$
$\displaystyle\leq e^{\mathrm{Tr}(E)}e^{\left\|E\right\|_{F}^{2}},$

and thus

$\left|D\overline{T}_{x}(v)\right|\geq(1-10^{-6})\cdot\frac{h^{n}}{\left|g(x_{\mathrm{mid}})\right|}.$

Since $\log\det g(x)=\log\det\nabla^{2}\phi(x)$ is convex in $x$, it follows
that

$\displaystyle\log\left|g(x_{\mathrm{mid}})\right|$
$\displaystyle=\log\left|g\left(\frac{x_{\frac{1}{3}}+x_{{\frac{2}{3}}}}{2}\right)\right|\leq\frac{1}{2}\left(\log\left|g(x_{\frac{1}{3}})\right|+\log\left|g(x_{{\frac{2}{3}}})\right|\right)=\frac{1}{2}\log\left|g(x_{\frac{1}{3}})g(x_{{\frac{2}{3}}})\right|$
$\displaystyle=\log\sqrt{\left|g(x)\right|\left|g(\bar{x})\right|},$

and thus
$\left|D\overline{T}_{x}(v)\right|\geq\frac{(1-10^{-6})h^{n}}{\sqrt{\left|g(x)\right|\left|g(\bar{x})\right|}}$.
Due to the step-size conditions of $h^{2}\sqrt{n}\bar{R}_{1}\leq 10^{-5}$ and
in Proposition 28, we can use Proposition 30 to conclude that
$\overline{T}_{x}(v)$ is sensitive at $(x,v)$. ∎

###  7.2 Generalized Leapfrog method (Störmer–Verlet)

We now analyze the generalized Leapfrog method (Algorithm 3), which is
symplectic and reversible in the Riemannian settings. In a similar way we
analyzed IMM, we show that if step size $h$ satisfies
$h^{2}\left(n+\sqrt{M_{1}}\right)\leq 10^{-10}$, then LM is second-order for
$x\in\mathbb{\mathcal{M}}_{\rho}$ and $v\in V_{\text{good}}^{x}$ with
$C_{x}(x,v)=O\left(n+\sqrt{M_{1}}\right)$ and
$C_{v}(x,v)=O\left(\sqrt{n+\sqrt{M_{1}}}\left(n+\sqrt{M_{1}}+M_{2}^{*}\right)\right)$.
Next, if the step size $h$ satisfies
$h^{2}\leq\min\left(\frac{10^{-20}}{n^{2}(n+\sqrt{M_{1}})},\frac{10^{-5}}{\sqrt{n}\bar{R}_{1}}\right)$
in addition to the step-size conditions in Proposition 28, then LM is
sensitive at $x\in\mathbb{\mathcal{M}}_{\rho}$ and $v\in V_{\text{good}}^{x}$.

Input: Initial point $x$, velocity $v$, step size $h$

// Step 1: Update $v$ (Implicit)

Find $v_{\frac{1}{2}}$ such that $v_{\frac{1}{2}}\leftarrow
v-\frac{h}{2}\frac{\partial H(x,v)}{\partial x}$.

// Step 2: Update $x$ (Implicit)

Find $x_{1}$ such that  $\displaystyle x_{1}$
$\displaystyle=x+\frac{h}{2}\left(\frac{\mathcal{\partial}H}{\mathcal{\partial}v}\left(x,v_{\frac{1}{2}}\right)+\frac{\mathcal{\partial}H}{\mathcal{\partial}v}\left(x_{1},v_{\frac{1}{2}}\right)\right).$

// Step 3: Update $v$ (Explicit)

Set $v_{1}\leftarrow v_{\frac{1}{2}}-\frac{h}{2}\frac{\partial H}{\partial
x}\left(x_{1},v_{\frac{1}{2}}\right)$.

Output: $x_{1,}v_{1}$

Algorithm 3 Generalized Leapfrog Method

####  7.2.1 Second-order

######  Lemma 35.

For $x\in\mathbb{\mathcal{M}}_{\rho}$ and $v\in V_{\text{good}}^{x}$, let
$g=g(x)$ and $h$ step size of LM with $h^{2}\left(n+\sqrt{M_{1}}\right)\leq
10^{-10}$. Let $(\overline{x},\overline{v})$ be the point obtained from RHMC
discretized by LM with the step size $h$ and initial condition $(x,v)$.

  1. 1.

$\left\|x-\overline{x}\right\|_{g}=O\left(h\sqrt{n}+h^{2}\left(n+\sqrt{M_{1}}\right)\right)$
and $\left\|v-\overline{v}\right\|_{g^{-1}}=O\left(h\left(\left\|\nabla
f(\bar{x})\right\|_{g^{-1}}+n+\sqrt{M_{1}}\right)\right)$.

  2. 2.

$C_{x}(x,v)=O\left(n+\sqrt{M_{1}}\right)$.

  3. 3.

$C_{v}(x,v)=O\left(\sqrt{n+\sqrt{M_{1}}}\left(n+\sqrt{M_{1}}+M_{2}^{*}\right)\right)$.

##### Proof of 1.

Let $\bar{x}=x_{1}=\overline{T}_{x}(v)$ and $v_{1}(=\bar{v}),v_{\frac{1}{2}}$
be the velocity obtained from LM with the initial condition $(x,v)$ and the
step size $h$. Let $g_{1}=g(x_{1})$. As $v_{\frac{1}{2}}\to v$ as $h\to 0$, we
can take $h_{0}>0$ such that
$h\left(\left\|v_{\frac{1}{2}}\right\|_{g^{-1}}+\left\|g_{1}^{-1}v_{\frac{1}{2}}\right\|_{g}\right)\leq\frac{2}{1000}$
for $h\leq h_{0}$ with the equality held at $h=h_{0}$. Thus for $h\leq h_{0}$
we have $h\left\|v_{\frac{1}{2}}\right\|_{g^{-1}}\leq\frac{1}{500}$.

From the first step of Algorithm 3 and Lemma 27, for step size $h\leq h_{0}$
it follows from $x\in\mathbb{\mathcal{M}}_{\rho}$ that

$\left\|v_{\frac{1}{2}}\right\|_{g^{-1}}\leq\left\|v\right\|_{g^{-1}}+\frac{h}{2}\left(\sqrt{M_{1}}+n+\left\|v_{\frac{1}{2}}\right\|_{g^{-1}}^{2}\right).$

Multiplying $h$ to both sides, and using
$h\left\|v_{\frac{1}{2}}\right\|_{g^{-1}}\leq\frac{1}{500}$ and $v\in
V_{\text{good}}^{x}$,

$\displaystyle h\left\|v_{\frac{1}{2}}\right\|_{g^{-1}}$ $\displaystyle\leq
150h\sqrt{n}+\frac{h^{2}}{2}\left(\sqrt{M_{1}}+n\right)+\frac{h^{2}}{2}\left\|v_{\frac{1}{2}}\right\|_{g^{-1}}^{2}$
$\displaystyle\leq
150h\sqrt{n}+\frac{h^{2}}{2}\left(\sqrt{M_{1}}+n\right)+\frac{1}{1000}h\left\|v_{\frac{1}{2}}\right\|_{g^{-1}},$

and thus

$\left\|v_{\frac{1}{2}}\right\|_{g^{-1}}\leq
200\sqrt{n}+h^{2}\left(n+\sqrt{M_{1}}\right).$ (7.10)

From the second step of Algorithm 3 and Lemma 27, for step size $h\leq h_{0}$

$\displaystyle\left\|x_{1}-x\right\|_{g}$
$\displaystyle\leq\frac{h}{2}\left(\left\|\frac{\mathcal{\partial}H}{\mathcal{\partial}v}\left(x,v_{\frac{1}{2}}\right)\right\|_{g}+\left\|\frac{\mathcal{\partial}H}{\mathcal{\partial}v}\left(x_{1},v_{\frac{1}{2}}\right)\right\|_{g}\right)\leq\frac{h}{2}\left(\left\|v_{\frac{1}{2}}\right\|_{g^{-1}}+\left\|g_{1}^{-1}v_{\frac{1}{2}}\right\|_{g}\right),$

and thus its is obvious that $\left\|x-x_{1}\right\|_{g}\leq\frac{1}{500}$. We
now lower bound $h_{0}$ as follows:

$\displaystyle\frac{1}{500}$
$\displaystyle=h_{0}\left(\left\|v_{\frac{1}{2}}\right\|_{g^{-1}}+\left\|g_{1}^{-1}v_{\frac{1}{2}}\right\|_{g}\right)\leq
h_{0}\left(\left\|v_{\frac{1}{2}}\right\|_{g^{-1}}+1.1\left\|v_{\frac{1}{2}}\right\|_{g^{-1}}\right)$
$\displaystyle\leq 3h_{0}\left\|v_{\frac{1}{2}}\right\|_{g^{-1}}$
$\displaystyle\leq 600h_{0}\sqrt{n}+3h_{0}^{2}\left(n+\sqrt{M_{1}}\right),$

where in the first inequality we switched the local norm at from $x_{1}$ to
$x$ due to $\left\|x_{1}-x\right\|_{g}\leq\frac{1}{500}$ and used (7.10) in
the last inequality. Therefore,
$h_{0}\geq\frac{1}{10^{4}\sqrt{n+\sqrt{M_{1}}}}$ and for step size
$h\leq\frac{1}{10^{5}\sqrt{n+\sqrt{M_{1}}}}$ we have

$\displaystyle\left\|x-x_{1}\right\|_{g}$ $\displaystyle\leq
600h\sqrt{n}+31h^{2}\left(n+\sqrt{M_{1}}\right),$
$\displaystyle\left\|v-v_{\frac{1}{2}}\right\|_{g^{-1}}$ $\displaystyle\leq
h\left(20000n+\sqrt{M_{1}}\right).$

Similarly, we can bound $\left\|v_{1}-v_{\frac{1}{2}}\right\|_{g^{-1}}$ by
Lemma 27:

$\displaystyle\left\|v_{1}-v_{\frac{1}{2}}\right\|_{g^{-1}}$
$\displaystyle\leq\frac{h}{2}\left\|\frac{\mathcal{\partial}H}{\mathcal{\partial}x}\left(x_{1},v_{\frac{1}{2}}\right)\right\|_{g^{-1}}\leq\frac{h}{2}\left(\left\|\nabla
f(x_{1})\right\|_{g^{-1}}+n+\left\|v_{\frac{1}{2}}\right\|_{g_{\mathrm{mid}}^{-1}}^{2}\right)$
$\displaystyle\leq 40000h\left(\left\|\nabla
f(x_{1})\right\|_{g^{-1}}+n+\sqrt{M_{1}}\right),$

and thus by adding it to the inequality for
$\left\|v-v_{\frac{1}{2}}\right\|_{g^{-1}}$ we have

$\left\|v-v_{1}\right\|_{g^{-1}}\leq 40000\left(\left\|\nabla
f(x_{1})\right\|_{g^{-1}}+n+\sqrt{M_{1}}\right).$

##### Proof of 2.

For $t\in[0,h]$, let $(x_{t},v_{t})$ be the Hamiltonian curve of the ideal
RHMC at time $t$ starting from $(x,v)$. Recall that

$\displaystyle T_{x}(v)$ $\displaystyle=x+\int_{0}^{h}\frac{\partial
H}{\partial v}(x_{t},v_{t})dt=x+\int_{0}^{h}g_{t}^{-1}v_{t}dt,$
$\displaystyle\overline{T}_{x}(v)$
$\displaystyle=\overline{x}=x+\frac{h}{2}\left(g^{-1}+g_{1}^{-1}\right)v_{\frac{1}{2}}.$

Thus,

$\displaystyle\left\|T_{x}(v)-\overline{T}_{x}(v)\right\|_{g}$
$\displaystyle=\left\|\left(x+\int_{0}^{h}g_{t}^{-1}v_{t}dt\right)-\left(x+\frac{h}{2}\left(g^{-1}+g_{1}^{-1}\right)v_{\frac{1}{2}}\right)\right\|_{g}$
$\displaystyle=\left\|\int_{0}^{h}\left(g_{t}^{-1}v_{t}-\frac{1}{2}\left(g^{-1}+g_{1}^{-1}\right)v_{\frac{1}{2}}\right)dt\right\|_{g}$
$\displaystyle\leq
h\max_{t\in[0,h]}\left\|g_{t}^{-1}v_{t}-\frac{1}{2}\left(g^{-1}+g_{1}^{-1}\right)v_{\frac{1}{2}}\right\|_{g}.$

By Lemma 26-11,

$\displaystyle\left\|g_{t}^{-1}v_{t}-\frac{1}{2}\left(g^{-1}+g_{1}^{-1}\right)v_{\frac{1}{2}}\right\|_{g}$
$\displaystyle\leq\frac{1}{2}\left\|g_{t}^{-1}v_{t}-g^{-1}v_{\frac{1}{2}}\right\|_{g}+\frac{1}{2}\left\|g_{t}^{-1}v_{t}-g_{1}^{-1}v_{\frac{1}{2}}\right\|_{g}$
$\displaystyle\lesssim\left\|v_{t}-v_{\frac{1}{2}}\right\|_{g^{-1}}+\left\|x_{t}-x\right\|_{g}\left\|v_{t}\right\|_{g^{-1}}$
$\displaystyle\quad+\left\|v_{t}-v_{\frac{1}{2}}\right\|_{g^{-1}}+\left\|x_{t}-x_{1}\right\|_{g}\left\|v_{t}\right\|_{g^{-1}}$
$\displaystyle\lesssim\left(\left\|v_{t}-v\right\|_{g^{-1}}+\left\|v-v_{\frac{1}{2}}\right\|_{g^{-1}}\right)$
$\displaystyle\quad+\left(\left\|x_{t}-x\right\|_{g}+\left\|x-x_{1}\right\|_{g}\right)\left(\left\|v_{t}-v\right\|_{g^{-1}}+\left\|v\right\|_{g^{-1}}\right).$

Using our bounds on $\left\|x_{1}-x\right\|_{g},\left\|x_{t}-x\right\|_{g}$
and
$\left\|v\right\|_{g^{-1}},\left\|v_{t}-v\right\|_{g^{-1}},\left\|v-v_{\frac{1}{2}}\right\|_{g^{-1}}$,
we conclude that
$\max_{t\in[0,h]}\left\|g_{t}^{-1}v_{t}-g_{\mathrm{mid}}^{-1}v_{\mathrm{mid}}\right\|_{g_{\mathrm{mid}}}\leq
10^{4}h\left(n+\sqrt{M_{1}}\right)$ and thus

$\left\|T_{x}(v)-\overline{T}_{x}(v)\right\|_{g}\leq
10^{4}h^{2}\left(n+\sqrt{M_{1}}\right).$

##### Proof of 3.

From the algorithm,

$\displaystyle v_{h}$ $\displaystyle=v-\int_{0}^{h}\frac{\partial H}{\partial
x}(x_{t},v_{t})dt,$ $\displaystyle\bar{v}$
$\displaystyle=v_{\frac{1}{2}}-\frac{h}{2}\frac{\mathcal{\partial}H}{\mathcal{\partial}x}\left(x_{1},v_{\frac{1}{2}}\right)$
$\displaystyle=v-\frac{h}{2}\frac{\mathcal{\partial}H}{\mathcal{\partial}x}\left(x_{1},v_{\frac{1}{2}}\right)-\frac{h}{2}\frac{\partial
H(x,v)}{\partial x}.$

Thus,

$\displaystyle\left\|v_{h}-\bar{v}\right\|_{g^{-1}}$
$\displaystyle=\left\|\int_{0}^{h}\frac{1}{2}\left(\frac{\mathcal{\partial}H}{\mathcal{\partial}x}(x_{t},v_{t})-\frac{\mathcal{\partial}H}{\mathcal{\partial}x}\left(x_{1},v_{\frac{1}{2}}\right)\right)dt+\int_{0}^{h}\frac{1}{2}\left(\frac{\mathcal{\partial}H}{\mathcal{\partial}x}(x_{t},v_{t})-\frac{\mathcal{\partial}H}{\mathcal{\partial}x}\left(x,v\right)\right)dt\right\|_{g^{-1}}$
$\displaystyle\leq\frac{h}{2}\underbrace{\max_{t\in[0,h]}\left\|\frac{\mathcal{\partial}H}{\mathcal{\partial}x}(x_{t},v_{t})-\frac{\mathcal{\partial}H}{\mathcal{\partial}x}\left(x_{1},v_{\frac{1}{2}}\right)\right\|_{g^{-1}}}_{F}++\frac{h}{2}\underbrace{\max_{t\in[0,h]}\left\|\frac{\mathcal{\partial}H}{\mathcal{\partial}x}(x_{t},v_{t})-\frac{\mathcal{\partial}H}{\mathcal{\partial}x}\left(x,v\right)\right\|_{g^{-1}}}_{S}.$

For $\delta_{v}=\left\|v_{t}-v_{\frac{1}{2}}\right\|_{g^{-1}}$ and
$\delta_{x}=\left\|x_{t}-x_{1}\right\|_{g}$, we use Proposition 26-12 to show
that

$\displaystyle F$
$\displaystyle\lesssim\max_{t\in[0,h]}\left(\delta_{v}+\delta_{x}\left\|v_{\frac{1}{2}}\right\|_{g^{-1}}\right)\left(\left\|v_{\frac{1}{2}}\right\|_{g^{-1}}+\left\|v_{t}\right\|_{g^{-1}}\right)+M_{2}^{*}\delta_{x}.$

In a similar way, $S$ can be bounded as follows:

$S\lesssim h\left(n+\sqrt{M_{1}}\right)^{3/2}+M_{2}^{*}\delta_{x}.$

Using our bounds on $\left\|x_{1}-x\right\|_{g},\left\|x_{t}-x\right\|_{g}$
and
$\left\|v\right\|_{g^{-1}},\left\|v_{t}-v\right\|_{g^{-1}},\left\|v-v_{\frac{1}{2}}\right\|_{g^{-1}}$,

$F+S\lesssim
h\left(n+\sqrt{M_{1}}\right)^{3/2}+h\sqrt{n+\sqrt{M_{1}}}M_{2}^{*}.$

Substituting the bounds on $F$ and $S$, we can conclude that

$\left\|v_{h}-\bar{v}\right\|_{g^{-1}}\leq
10^{10}\left(\left(n+\sqrt{M_{1}}\right)^{3/2}+\sqrt{n+\sqrt{M_{1}}}M_{2}^{*}\right)h^{2}.$

####  7.2.2 Sensitivity

We show that for some step size $h$ the generalized Leapfrog integrator is
sensitive at $(x,v)$ for $x\in\mathbb{\mathcal{M}}_{\rho}$ and $v\in
V_{\text{good}}^{x}$.

######  Lemma 36.

For $x\in\mathbb{\mathcal{M}}_{\rho}$ and $v\in V_{\text{good}}^{x}$, LM is
sensitive at $(x,v)$ if step size $h$ satisfies
$h^{2}\leq\min\left(\frac{10^{-10}}{n^{2}\left(n+\sqrt{M_{1}}\right)},\frac{10^{-5}}{\sqrt{n}\bar{R}_{1}}\right)$
and the step-size conditions in Proposition 28.

###### Proof.

Let $\overline{T}_{x}(v)=\bar{x}=x_{1}$. We lower bound
$\left|D\overline{T}_{x}(v)\right|$. It suffices to look into the determinant
of composition of first two steps in Algorithm 3, since the third step only
changes $v$. The first two steps are

$\displaystyle v_{\frac{1}{2}}$
$\displaystyle=v-\frac{h}{2}\frac{\mathcal{\partial}H}{\mathcal{\partial}x}(x,v_{\frac{1}{2}}),$
$\displaystyle x_{1}$
$\displaystyle=x+\frac{h}{2}\left(\frac{\mathcal{\partial}H}{\mathcal{\partial}v}(x,v_{\frac{1}{2}})+\frac{\mathcal{\partial}H}{\mathcal{\partial}v}(x_{1},v_{\frac{1}{2}})\right).$

Differentiating the first equation with respect to $v$, we have

$\frac{\mathcal{\partial}v_{\frac{1}{2}}}{\mathcal{\partial}v}=I_{n}-\frac{h}{2}\frac{\mathcal{\partial}^{2}H}{\mathcal{\partial}x\mathcal{\partial}v}\left(x,v_{\frac{1}{2}}\right)\frac{\mathcal{\partial}v_{\frac{1}{2}}}{\mathcal{\partial}v},$

and so

$\left(I+\frac{h}{2}\frac{\mathcal{\partial}^{2}H}{\mathcal{\partial}x\mathcal{\partial}v}\left(x,v_{\frac{1}{2}}\right)\right)\frac{\mathcal{\partial}v_{\frac{1}{2}}}{\mathcal{\partial}v}=I_{n}.$
(7.11)

Differentiating the second equation with respect to $v$, we obtain

$\frac{\mathcal{\partial}x_{1}}{\mathcal{\partial}v}=\frac{h}{2}\left(\frac{\mathcal{\partial}^{2}H}{\mathcal{\partial}v^{2}}\left(x,v_{\frac{1}{2}}\right)\frac{\mathcal{\partial}v_{\frac{1}{2}}}{\mathcal{\partial}v}+\frac{\mathcal{\partial}^{2}H}{\mathcal{\partial}x\mathcal{\partial}v}\left(x_{1},v_{\frac{1}{2}}\right)\frac{\mathcal{\partial}x_{1}}{\mathcal{\partial}v}+\frac{\mathcal{\partial}^{2}H}{\mathcal{\partial}v^{2}}\left(x_{1},v_{\frac{1}{2}}\right)\frac{\mathcal{\partial}v_{\frac{1}{2}}}{\mathcal{\partial}v}\right).$

Collecting all $\mathcal{\partial}x_{1}/\mathcal{\partial}v$ terms from this
equation, for $g=g(x)$ and $g_{1}=g(x_{1})$

$\displaystyle\left(I_{n}-\frac{h}{2}\frac{\mathcal{\partial}^{2}H}{\mathcal{\partial}x\mathcal{\partial}v}\left(x_{1},v_{\frac{1}{2}}\right)\right)\frac{\mathcal{\partial}x_{1}}{\mathcal{\partial}v}$
$\displaystyle=\frac{h}{2}\left(\frac{\mathcal{\partial}^{2}H}{\mathcal{\partial}v^{2}}\left(x,v_{\frac{1}{2}}\right)+\frac{\mathcal{\partial}^{2}H}{\mathcal{\partial}v^{2}}\left(x_{1},v_{\frac{1}{2}}\right)\right)\frac{\mathcal{\partial}v_{\frac{1}{2}}}{\mathcal{\partial}v}$
$\displaystyle=\frac{h}{2}\left(g^{-1}+g_{1}^{-1}\right)\left(I+\frac{h}{2}\frac{\mathcal{\partial}^{2}H}{\mathcal{\partial}x\mathcal{\partial}v}\left(x,v_{\frac{1}{2}}\right)\right)^{-1},$

where we used (7.11). Hence,

$\displaystyle\frac{\mathcal{\partial}x_{1}}{\mathcal{\partial}v}$
$\displaystyle=h\left(I_{n}-\frac{h}{2}\frac{\mathcal{\partial}^{2}H}{\mathcal{\partial}x\mathcal{\partial}v}\left(x_{1},v_{\frac{1}{2}}\right)\right)^{-1}\left(\frac{g^{-1}+g_{1}^{-1}}{2}\right)\left(I_{n}+\frac{h}{2}\frac{\mathcal{\partial}^{2}H}{\mathcal{\partial}x\mathcal{\partial}v}\left(x,v_{\frac{1}{2}}\right)\right)^{-1}$
$\displaystyle=h\left(I_{n}-\frac{h}{2}g_{1}^{-1}Dg_{1}\left[g_{1}^{-1}v_{\frac{1}{2}}\right]\right)^{-1}\left(\frac{g^{-1}+g_{1}^{-1}}{2}\right)\left(I_{n}+\frac{h}{2}g^{-1}Dg\left[g^{-1}v_{\frac{1}{2}}\right]\right)^{-1}$
$\displaystyle=hg_{1}^{\frac{1}{2}}\left(I_{n}-\frac{h}{2}g_{1}^{-\frac{1}{2}}Dg_{1}\left[g_{1}^{-1}v_{\frac{1}{2}}\right]g_{1}^{-\frac{1}{2}}\right)^{-1}g_{1}^{-\frac{1}{2}}\left(\frac{g^{-1}+g_{1}^{-1}}{2}\right)g^{\frac{1}{2}}\left(I_{n}+\frac{h}{2}g^{-\frac{1}{2}}Dg\left[g^{-1}v_{\frac{1}{2}}\right]g^{-\frac{1}{2}}\right)^{-1}g^{-\frac{1}{2}}.$

Due to the concavity of log-determinant in the set of positive definite
matrices, we have

$\displaystyle\log\left|\frac{g^{-1}+g_{1}^{-1}}{2}\right|$
$\displaystyle\geq\frac{1}{2}\left(\log\left|g^{-1}\right|+\log\left|g_{1}^{-1}\right|\right)=\log\frac{1}{\sqrt{\left|g\right|\left|g_{1}\right|}},$

and thus

$\displaystyle\left|D\overline{T}_{x}(v)\right|$
$\displaystyle=\left|\frac{\mathcal{\partial}x_{1}}{\mathcal{\partial}v}\right|$
$\displaystyle\geq\frac{h^{n}}{\sqrt{\left|g\right|\left|g_{1}\right|}}\left|I_{n}-\frac{h}{2}g_{1}^{-\frac{1}{2}}Dg_{1}\left[g_{1}^{-1}v_{\frac{1}{2}}\right]g_{1}^{-\frac{1}{2}}\right|^{-1}\left|I_{n}+\frac{h}{2}g^{-\frac{1}{2}}Dg\left[g^{-1}v_{\frac{1}{2}}\right]g^{-\frac{1}{2}}\right|^{-1}.$

For $E\stackrel{{\scriptstyle\mathrm{{\scriptscriptstyle
def}}}}{{=}}\frac{h}{2}g^{-\frac{1}{2}}Dg\left[g^{-1}v_{\frac{1}{2}}\right]g^{-\frac{1}{2}}$,
as bounded in (7.8), we have that

$\displaystyle\left\|E\right\|_{2}$
$\displaystyle\leq\frac{h}{2}\left\|v_{\frac{1}{2}}\right\|_{g^{-1}},$
$\displaystyle\mathrm{Tr}(E)$ $\displaystyle\lesssim
hn\left\|v_{\frac{1}{2}}\right\|_{g^{-1}},$
$\displaystyle\left\|E\right\|_{F}$
$\displaystyle\leq\frac{h\sqrt{n}}{2}\left\|v_{\frac{1}{2}}\right\|_{g^{-1}}.$

Due to $h^{2}\leq\frac{10^{-10}}{n^{2}\left(n+\sqrt{M_{1}}\right)}$, it
follows from (7.10) that $\left\|v_{\frac{1}{2}}\right\|_{g^{-1}}\leq
O\left(\sqrt{n+\sqrt{M_{1}}}\right)$. This condition also allows us to make
all these three quantities smaller than $10^{-5}$. By Lemma 60,

$\left|I_{n}-\frac{h}{2}g_{1}^{-\frac{1}{2}}Dg_{1}\left[g_{1}^{-1}v_{\frac{1}{2}}\right]g_{1}^{-\frac{1}{2}}\right|^{-1}\geq
1-10^{-8}.$

Similarly, we obtain

$\left|I_{n}+\frac{h}{2}g^{-\frac{1}{2}}Dg\left[g^{-1}v_{\frac{1}{2}}\right]g^{-\frac{1}{2}}\right|^{-1}\geq
1-10^{-8},$

and thus
$\left|D\overline{T}_{x}(v)\right|\geq\frac{(1-10^{-6})h^{n}}{\sqrt{\left|g(x)\right|\left|g(\bar{x})\right|}}$.
Using the step-size conditions in Proposition 28, we use Proposition 30 to
show that $\overline{T}_{x}$ is sensitive at $(x,v)$. ∎

##  8 Convergence rate of RHMC in polytopes

In this section, we present the mixing times of the ideal and discretized RHMC
for an exponential density in a polytope. We set $f(x)=\alpha^{\top}x$ for
$\alpha\in\mathbb{R}^{n}$. For a full-rank matrix $A\in\mathbb{R}^{m\times n}$
and $b\in\mathbb{R}^{m}$, the polytope is represented by
$\left\\{x\in\mathbb{R}^{n}:Ax\geq b\right\\}$, equipped with the logarithmic
barrier $\phi(x)=-\sum_{i=1}^{m}\log(a_{i}^{\top}x-b_{i}),$ where $a_{i}$ is
the $i^{th}$ row of $A$ and $b_{i}$ is the $i^{th}$ entry of $b$. We can check
by direct computation that the logarithmic barriers are highly self-
concordant. We view this polytope as the Hessian manifold
$\mathbb{\mathcal{M}}$ induced by the local norm $g(x)=\nabla^{2}\phi(x)$. We
denote a slack vector by $s_{x}=Ax-b\in\mathbb{R}^{m}$ and its diagonalization
by $S_{x}=\mathrm{Diag}(s_{x})\in\mathbb{R}^{m\times m}$. We also define
$A_{x}=S_{x}^{-1}A$ and $s_{v}=A_{x}v$ for $v\in T_{x}\mathbb{\mathcal{M}}$,
where $T_{x}\mathbb{\mathcal{M}}$ is endowed with the local metric $g$. One
can check by direct computation that $\nabla^{2}\phi(x)=A_{x}^{\top}A_{x}$.

In this setting, we can compute all the parameters we have defined, obtaining
the mixing time of RHMC discretized by a sensitive numerical integrator.

###  8.1 Isoperimetry of convex set

An isoperimetry inequality is one of the two main ingredients for bounding the
mixing rate. We use the Riemannian version of some isoperimetry inequality. To
state it, we need another distance called _Hilbert distance_ in addition to
Riemannian distance $d_{\phi}$.

######  Definition 37.

For a convex body $\mathcal{K}$, the cross-ratio distance
$d_{\mathcal{K}}(x,y)$ between $x$ and $y$ is

$d_{\mathcal{K}}(x,y)=\frac{|x-y||p-q|}{|p-x||y-q|},$

where $p$ and $q$ are on the boundary of $\mathcal{K}$ such that $p,x,y,q$ are
on the straight line $\overline{xy}$ and are in order. The Hilbert distance
$d_{H}$ between $x,y\in\mathcal{K}$ is

$d_{H}(x,y)=\log(1+d_{\mathcal{K}}(x,y))=\log\left(1+\frac{|x-y||p-q|}{|p-x||y-q|}\right).$

For sets $X$ and $Y$, we define $d_{*}(X,Y)=\inf_{x\in X,y\in Y}d_{*}(x,y)$
for $*\in\\{\mathcal{K},H,\phi\\}$.

######  Lemma 38 ([Vem05], Theorem 4.4).

Let $\pi$ be a log-concave distribution supported on a convex body
$\mathcal{K}$. Let $S_{1},S_{2},S_{3}$ be a partition of $\mathcal{K}$. Then,

$\pi(S_{3})\geq d_{\mathcal{K}}(S_{1},S_{2})\pi(S_{1})\pi(S_{2}).$

The following lemma is a generalization of Theorem 26 in [LV17] to a subset
$\mathcal{K}^{\prime}$.

######  Lemma 39.

Let $\pi$ be a log-concave distribution supported on a convex body
$\mathcal{K}$, and $\phi$ a self-concordant barrier of $\mathcal{K}$. Let
$\mathcal{K}^{\prime}$ be a convex subset of $\mathcal{K}$, and
$S_{1},S_{2},S_{3}$ a partition of $\mathcal{K}^{\prime}$. Then

$\pi(S_{3})\pi(\mathcal{K}^{\prime})\geq\frac{d_{\phi}(S_{1},S_{2})}{G}\pi(S_{1})\pi(S_{2}),$

where $G=\sup_{x,y\in\mathcal{K}}\frac{d_{\phi}(x,y)}{d_{H}(x,y)}$.

###### Proof.

Applying Lemma 38 to the distribution $\pi_{\mathcal{K}^{\prime}}$ defined by
$\pi$ restricted to $\mathcal{K}^{\prime}$, we have

$\pi(S_{3})\pi(\mathcal{K}^{\prime})\geq
d_{\mathcal{K}^{\prime}}(S_{1},S_{2})\pi(S_{1})\pi(S_{2}).$

Due to $\mathcal{K}^{\prime}\subseteq\mathcal{K}$, one can check
$d_{\mathcal{K}^{\prime}}(S_{1},S_{2})\geq d_{\mathcal{K}}(S_{1},S_{2})$ by
simple algebra. As $d_{\mathcal{K}}(x,y)\geq d_{H}(x,y)$, it follows that

$\pi(S_{3})\pi(\mathcal{K}^{\prime})\geq\frac{d_{\phi}(S_{1},S_{2})}{\frac{d_{\phi}(S_{1},S_{2})}{d_{H}(S_{1},S_{2})}}\pi(S_{1})\pi(S_{2})\geq\frac{d_{\phi}(S_{1},S_{2})}{G}\pi(S_{1})\pi(S_{2}).$

∎

We now define the symmetric self-concordance parameter of the barrier $\phi$.

######  Definition 40 ([LLV20]).

For a convex body $\mathcal{K}\subseteq\mathbb{R}^{n}$, the symmetric self-
concordance parameter $\bar{\nu}_{\phi}$ of $\mathcal{K}$ is the smallest
number such that for any $x\in\mathcal{K}$

$D(x)\subseteq\mathcal{K}\cap(2x-\mathcal{K})\subseteq\sqrt{\bar{\nu}_{\phi}}D(x),$

where
$D(x)=\left\\{y\in\mathbb{R}^{n}:\left\|y-x\right\|_{\nabla^{2}\phi(x)}\leq
1\right\\}$ is the Dikin ellipsoid at $x$.

In general, it is known that $\bar{\nu}_{\phi}=O(\nu_{\phi}^{2})$ for the
self-concordance parameter $\nu_{\phi}$ (see Definition 54), but a tighter
bound of $\bar{\nu}_{\phi}=O(\nu_{\phi})$ holds for important barriers such as
the logarithmic barrier and Lee-Sidford barrier [LS14].

######  Lemma 41 ([LLV20], Lemma 2.3).

$d_{\phi}(x,y)\lesssim\sqrt{\bar{\nu}_{\phi}}d_{H}(x,y)$ for any
$x,y\in\mathcal{K}$.

Using Lemma 39 and 41 together, we have

$\pi(S_{3})\pi(\mathcal{K}^{\prime})\geq\frac{d_{\phi}(S_{1},S_{2})}{\sqrt{\bar{\nu}_{\phi}}}\pi(S_{1})\pi(S_{2}),$

and it implies that the isoperimetry of $\mathcal{K}^{\prime}$ is at least
$1/\sqrt{\bar{\nu}_{\phi}}$. As $\nu_{\phi}=O(m)$ for the logarithmic barrier,
$\psi_{\mathcal{K}^{\prime}}\geq 1/\sqrt{m}$ for a convex subset
$\mathcal{K}^{\prime}$.

###  8.2 Good region $\mathbb{\mathcal{M}}_{\rho}$

Taking a proper good region $\mathbb{\mathcal{M}}_{\rho}$ plays an important
role in establishing a condition-number independent mixing rate of RHMC for an
exponential density in a polytope. To this end, we set our good region to

$\displaystyle\mathbb{\mathcal{M}}_{\rho}$
$\displaystyle\stackrel{{\scriptstyle\mathrm{{\scriptscriptstyle
def}}}}{{=}}\left\\{x\in\mathbb{\mathcal{M}}:\left\|\alpha\right\|_{g(x)^{-1}}^{2}\leq
10n^{2}\log^{2}\left(\frac{1}{\rho}\right)\right\\}.$

To establish the isoperimetry of $\mathbb{\mathcal{M}}_{\rho}$ following
Section 8.1, we check its convexity in the following lemma. Note that the
assumption in the lemma is satisfied by the logarithmic barriers.

######  Lemma 42.

If the fourth directional derivative of $\phi$ is positive (i.e.,
$D^{4}\phi[a,a,b,b]\geq 0$), then $\mathbb{\mathcal{M}}_{\rho}$ is convex.

###### Proof.

Let
$\Upsilon(x):=\alpha^{\top}g(x)^{-1}\alpha=\alpha^{\top}(\nabla^{2}\phi(x))^{-1}\alpha$.
It suffices to show that $\Upsilon(x)$ is convex. Note that

$\frac{\partial\Upsilon(x)}{\partial
x_{i}}=\alpha^{\top}g(x)^{-1}\frac{\partial g(x)}{\partial
x_{i}}g(x)^{-1}\alpha,$

and thus its directional derivative in $h=(h_{1},...,h_{n})$ is

$\nabla\Upsilon(x)\cdot h=\sum_{i}h_{i}\left(s(x)^{\top}\frac{\partial
g(x)}{\partial x_{i}}s(x)\right),$

where $s(x):=g(x)^{-1}\alpha$. Note that

$\displaystyle\frac{\partial}{\partial x_{j}}\left(s(x)^{\top}\frac{\partial
g(x)}{\partial x_{i}}s(x)\right)$
$\displaystyle=s(x)^{\top}\frac{\partial^{2}g(x)}{\partial x_{i}\partial
x_{j}}s(x)+2s(x)^{\top}\frac{\partial g(x)}{\partial
x_{i}}\left(\frac{\partial s(x)}{\partial x_{j}}\right)$
$\displaystyle=s(x)^{\top}\frac{\partial^{2}g(x)}{\partial x_{i}\partial
x_{j}}s(x)+2s(x)^{\top}\frac{\partial g(x)}{\partial
x_{i}}g(x)^{-1}\frac{\partial g(x)}{\partial x_{j}}s(x).$

Therefore,

$\displaystyle D^{2}\Upsilon(x)[h,h]$
$\displaystyle=\sum_{i,j}h_{i}h_{j}s(x)^{\top}\frac{\partial^{2}g(x)}{\partial
x_{i}\partial x_{j}}s(x)+2\sum_{i,j}\left(\frac{\partial g(x)}{\partial
x_{i}}s(x)h_{i}\right)^{\top}g(x)^{-1}\left(\frac{\partial g(x)}{\partial
x_{j}}s(x)h_{j}\right)$
$\displaystyle=D^{4}\phi[h,h,s(x),s(x)]+2\sum_{i,j}\left(\frac{\partial
g(x)}{\partial x_{i}}s(x)h_{i}\right)^{\top}g(x)^{-1}\left(\frac{\partial
g(x)}{\partial x_{j}}s(x)h_{j}\right).$

The first term is non-negative due to the assumption, and the second term is
also non-negative since $g(x)^{-1}$ is also positive semi-definite. ∎

Next, we show that $\mathbb{\mathcal{M}}_{\rho}$ takes up probability of at
least $1-\rho$ over the stationary distribution $\pi$, where
$\frac{d\pi(x)}{dx}\propto\exp(-\alpha^{\top}x)$.

######  Lemma 43.

$\pi(\mathbb{\mathcal{M}}_{\rho})\geq 1-\rho$.

###### Proof.

Let $g=g(x)$. For $\left\|\alpha\right\|_{g^{-1}}$, note that

$\displaystyle\left\|\alpha\right\|_{g(x)^{-1}}=$
$\displaystyle\max_{\left\|u\right\|_{g(x)}=1}\alpha^{\top}u$ $\displaystyle=$
$\displaystyle\alpha^{\top}x-\min_{\left\|y-x\right\|_{g}=1}\alpha^{\top}y\leq\alpha^{\top}x-\min_{y\in\mathbb{\mathcal{M}}}\alpha^{\top}y,$

where the first equality is due to duality of norms and the last inequality
follows from the well-known fact that the Dikin ellipsoid at $x$ is inside
$\mathbb{\mathcal{M}}$.

By Lemma 61 with $c=\alpha/\left\|\alpha\right\|_{2}$ and
$T=1/\left\|\alpha\right\|_{2}$, we have

$\mathbb{E}_{x\sim\pi^{*}}[\alpha^{\top}x]\leq
n+\min_{y\in\mathbb{\mathcal{M}}}\alpha^{\top}y.$

Then, by Lemma 62 we have
$\mathbb{E}[(\alpha^{\top}x-\min_{y\in\mathbb{\mathcal{M}}}\alpha^{\top}y)^{2}]-\mathbb{E}[\alpha^{\top}x-\min_{y\in\mathbb{\mathcal{M}}}\alpha^{\top}y]^{2}\leq
n$ so that
$\mathbb{E}[(\alpha^{\top}x-\min_{y\in\mathbb{\mathcal{M}}}\alpha^{\top}y)^{2}]\leq
n+n^{2}.$ By Lemma 63, we have

$\Pr_{x\sim\pi}\left[\alpha^{\top}x-\min_{y\in\mathbb{\mathcal{M}}}\alpha^{\top}y>2\left(\log\frac{1}{\rho}+1\right)n\right]\leq\rho.$

∎

###  8.3 Auxiliary function $\ell$ and smoothness parameters $R$

In this region $\mathbb{\mathcal{M}}_{\rho}$ and step size $h$, the parameters
$M_{1},M_{2}$ and $M_{1}^{*},M_{2}^{*}$ (see Definition 23) are computed by

$\displaystyle M_{1}$
$\displaystyle=\max\left(n,\left\|\alpha\right\|_{g(x)^{-1}}^{2}\right)\leq
10n^{2}\log^{2}\left(\frac{1}{\rho}\right),$ $\displaystyle M_{1}^{*}$
$\displaystyle\leq\left\|\alpha\right\|_{g(x)^{-1}}^{2}\leq
10n^{2}\log^{2}\left(\frac{1}{\rho}\right),$ $\displaystyle M_{2}$
$\displaystyle,M_{2}^{*}=0.$

We use the following auxiliary function $\ell$ proposed in [LV18] and
symmetric auxiliary function $\bar{\ell}$:

$\displaystyle\ell(\gamma)=\max_{t\in[0,h]}\bigg{(}\frac{\left\|s_{\gamma^{\prime}(t)}\right\|_{2}}{\sqrt{n}+2M_{1}^{1/4}}+\frac{\left\|s_{\gamma^{\prime}(t)}\right\|_{4}}{2M_{1}^{1/4}}+$
$\displaystyle\frac{\left\|s_{\gamma^{\prime}(t)}\right\|_{\infty}}{\sqrt{\log
n}+2h\sqrt{M_{1}}}\bigg{)}+\frac{\left\|s_{\gamma^{\prime}(0)}\right\|_{2}}{\sqrt{n}}+\frac{\left\|s_{\gamma^{\prime}(0)}\right\|_{4}}{n^{1/4}}+\frac{\left\|s_{\gamma^{\prime}(0)}\right\|_{\infty}}{\sqrt{\log
n}},$
$\displaystyle\bar{\ell}(\gamma)=\max_{t\in[0,h]}\bigg{(}\frac{\left\|s_{\gamma^{\prime}(t)}\right\|_{2}}{\sqrt{n}+2M_{1}^{1/4}}+\frac{\left\|s_{\gamma^{\prime}(t)}\right\|_{4}}{2M_{1}^{1/4}}+$
$\displaystyle\frac{\left\|s_{\gamma^{\prime}(t)}\right\|_{\infty}}{\sqrt{\log
n}+2h\sqrt{M_{1}}}\bigg{)}.$

This measures how fast a Hamiltonian trajectory approaches the facets of a
polytope in the local norm.

We make simple observations based on the self-concordance of $g$.

######  Proposition 44.

Let
$\overline{\mathbb{\mathcal{M}}_{\rho}}\stackrel{{\scriptstyle\mathrm{{\scriptscriptstyle
def}}}}{{=}}\left\\{x\in\mathbb{\mathcal{M}}:\left\|\alpha\right\|_{g(x)^{-1}}^{2}\leq
20n^{2}\log^{2}\left(\frac{1}{\rho}\right)\right\\}$ and $\gamma$ be any
Hamiltonian curve $\gamma$ starting at $x\in\mathbb{\mathcal{M}}_{\rho}$ with
$v\in V_{\text{good}}^{x}$. If step size $h$ satisfies $h^{2}\leq
10^{-11}\min\left(\frac{1}{n\log\frac{1}{\rho}},\frac{1}{C_{x}(x,v)}\right)$,
then $x_{h}$ and $\bar{x}_{h}$ are contained in
$\overline{\mathbb{\mathcal{M}}_{\rho}}$.

###### Proof.

Due to the assumption on the step size, we can use Proposition 28-1, obtaining
$\left\|x-\gamma(t)\right\|_{x}\leq
O\left(t\sqrt{n+\sqrt{M_{1}}}\right)=O\left(t\sqrt{n\log\frac{1}{\rho}}\right)<\frac{1}{4}$.
Also,
$\left\|x-\bar{x}_{h}\right\|_{g}\leq\left\|x-\gamma(h)\right\|_{g}+\left\|\gamma(h)-\bar{x}_{h}\right\|_{g}\leq\frac{1}{4}+h^{2}C_{x}(x,v)\leq\frac{1}{3}$.
The claim follows from the self-concordance of $g(x)$, due to
$\left\|\alpha\right\|_{g(\gamma(h))^{-1}}^{2}\leq\left(1+\left\|x-\gamma(h)\right\|_{x}\right)\left\|\alpha\right\|_{g(x)^{-1}}^{2}\leq
20n^{2}\log^{2}\frac{1}{\rho}$ and
$\left\|\alpha\right\|_{g(\bar{x}_{h})^{-1}}^{2}\leq\frac{4}{3}\left\|\alpha\right\|_{g(x)^{-1}}^{2}\leq
20n^{2}\log^{2}\frac{1}{\rho}$. ∎

As in [LV18], we can represent the parameters $\ell_{0},\ell_{1}$ and the
smoothness parameters $R_{1},R_{2},R_{3}$ in terms of $M_{1}$. The original
proof in [LV18] relies on the fact that $\left\|\nabla
f(\gamma(t))\right\|_{g(\gamma(t))^{-1}}^{2}\leq M_{1}$ for any time
$t\in[0,h]$ and any regular Hamiltonian curves. In our setting, $\left\|\nabla
f(\gamma(t))\right\|_{g(\gamma(t))^{-1}}^{2}\leq 2M_{1}$ for any time
$t\in[0,h]$ if $h^{2}\leq\frac{10^{-11}}{n\log\frac{1}{\rho}}$, we can simply
reproduce Lemma 54~59 by replacing $M_{1}$ by $2M_{1}$.

######  Lemma 45.

Consider a Hamiltonian trajectory $\gamma$ starting at
$x\in\mathbb{\mathcal{M}}_{\rho}$ with an initial (normalized) velocity
randomly chosen from $\mathcal{N}(0,g(x)^{-1})$, with step size $h$ satisfying
$h^{2}n\log\frac{1}{\rho}\leq 10^{-11}$. For $n$ large enough, if $s$
satisfies $sh=O\left(n\right)$, then

$\mathbf{P}_{\gamma}\left(\ell(\gamma)\geq
128\right)\leq\frac{1}{100}\min\left(1,\frac{\ell_{0}}{sh}\right).$

As we shortly see in Lemma 49, we have
$\ell_{1}h=O\left(h^{2}M_{1}^{1/4}\right)=O\left(\frac{1}{\sqrt{n\log\frac{1}{\rho}}}\right)$,
and thus $\ell_{1}$ can be used in place of $s$ in this lemma.

######  Lemma 46.

Let $\gamma$ be a Hamiltonian curve starting at
$x\in\mathbb{\mathcal{M}}_{\rho}$ with $\ell(\gamma)\leq\ell_{0}\leq 256$ and
step size $h$ satisfying $h^{2}n\log\frac{1}{\rho}\leq 10^{-11}$. Then

$\sup_{t\in[0,h]}\left\|\Phi(\gamma,t)\right\|_{F,\gamma(t)}\leq R_{1}$

with $R_{1}=O\left(\sqrt{M_{1}}\right)$.

######  Lemma 47.

Let $\gamma$ be a Hamiltonian curve starting at
$x\in\mathbb{\mathcal{M}}_{\rho}$ with $\ell(\gamma)\leq\ell_{0}\leq 256$ and
step size $h$ satisfying $h^{2}n\log\frac{1}{\rho}\leq 10^{-11}$. For any
$t\in[0,h]$, any curve $c(s)$ starting from $\gamma(t)$ and any velocity field
$v(c(s))$ on $c(s)$ with $v(c(0))=v(\gamma(t))=\gamma^{\prime}(t)$, we have
that

$\left|\frac{d}{ds}\mathrm{Tr}\Phi(v(c(s))\bigg{|}_{s=0}\right|\leq
R_{2}\left(\left\|\frac{dc}{ds}\bigg{|}_{s=0}\right\|_{\gamma(t)}+h\left\|D_{s}v\big{|}_{s=0}\right\|_{\gamma(t)}\right)$

with
$R_{2}=O\left(\sqrt{nM_{1}}+\sqrt{n}M_{1}h^{2}+\frac{M_{1}^{1/4}}{h}+\frac{\sqrt{n\log
n}}{h}\right)$.

######  Lemma 48.

Let $\gamma$ be a Hamiltonian curve starting at
$x\in\mathbb{\mathcal{M}}_{\rho}$ with $\ell(\gamma)\leq\ell_{0}\leq 256$ and
step size $h$ satisfying $h^{2}n\log\frac{1}{\rho}\leq 10^{-11}$. Let
$\zeta(t)$ be the parallel transport of the vector $\gamma^{\prime}(0)$ to
$\gamma(t)$. Then

$\sup_{t\in[0,h]}\left\|\Phi(\gamma,t)\zeta(t)\right\|_{\gamma(t)}\leq R_{3}$

with $R_{3}=O\left(\sqrt{M_{1}\log n}+M_{1}^{3/4}n^{1/4}h\right)$.

######  Lemma 49.

Let $\gamma_{s}$ be a Hamiltonian variation starting at
$x\in\mathbb{\mathcal{M}}_{\rho}$ with $\ell(\gamma_{s})\leq\ell_{0}\leq 256$
and step size $h$ satisfying $h^{2}n\log\frac{1}{\rho}\leq 10^{-11}$. Then

$\left|\frac{d}{ds}\ell(\gamma_{s})\right|\leq
O\left(M_{1}^{1/4}h+\frac{1}{h\sqrt{\log
n}}\right)\left(\left\|\frac{d}{ds}\gamma_{s}(0)\right\|_{\gamma_{s}(0)}+h\left\|D_{s}\gamma_{s}^{\prime}(0)\right\|_{\gamma_{s}(0)}\right),$

and thus $\ell_{1}=O\left(M_{1}^{1/4}h+\frac{1}{h\sqrt{\log n}}\right)$.

For $\bar{\ell}_{0},\bar{\ell}_{1},\bar{R}_{1}$, we can repeat the arguments
so far for regular Hamiltonian curves starting from
$\overline{\mathbb{\mathcal{M}}_{\rho}}$, in which
$\left\|\alpha\right\|_{g(\gamma(t))^{-1}}^{2}$ is within a constant factor of
$M_{1}$. Therefore, these three parameters also have the same bounds in Lemma
45, 46 and 49 up to a multiplicative constant factor.

###  8.4 Convergence rate of RHMC with numerical integrators

Now that we estimated all the parameters, we can put them together and state
the mixing time of RHMC discretized by a sensitive numerical integrator.

See 2

###### Proof.

We first note that
$V_{\text{good}}^{x}=\left\\{v\in\mathbb{R}^{n}:\bar{\ell}\left(\mathrm{Ham}_{x,t}\left(g(x)^{-1}v\right)\right)\leq
128\right\\}$, the measure of which is at least $0.99$ by the definition of
$\bar{\ell}_{0}$. We check the conditions on the step size in Theorem 24. Let
$\rho=\frac{\varepsilon}{2\Lambda}$. We first bound $M_{1},M_{1}^{*}$ by
$20n^{2}\log^{2}\frac{1}{\rho}$ and set $M_{2}$ to $0$. Substituting these to
Lemma 49, 46, 47 and 48, we have

$\displaystyle\ell_{1}$ $\displaystyle\lesssim
h\sqrt{n\log\frac{1}{\rho}}+\frac{1}{h},$ $\displaystyle R_{1}$
$\displaystyle\lesssim n\log\frac{1}{\rho},$ $\displaystyle R_{2}$
$\displaystyle\lesssim
n^{3/2}\log\frac{1}{\rho}+h^{2}n^{5/2}\log^{2}\frac{1}{\rho}+\frac{\sqrt{n\log\frac{1}{\rho}}}{h}+\frac{\sqrt{n\log
n}}{h},$ $\displaystyle R_{3}$ $\displaystyle\lesssim n\sqrt{\log
n}\log\frac{1}{\rho}+hn^{7/4}\log^{3/2}\frac{1}{\rho}.$

Due to $h\leq\frac{10^{-20}}{n^{7/12}\log^{1/2}\frac{1}{\rho}}$, direct
computation leads to
$h^{2}\max\left(R_{1},\bar{R}_{1}\right),h^{5}R_{1}^{2}\ell_{1}/\ell_{0},h^{3}R_{2}+h^{2}R_{3}\lesssim
1$ and $h\lesssim\min\left(1,\frac{\ell_{0}}{\ell_{1}}\right)$. The rest of
conditions on the step size,
$hC_{x}(x,v)\leq\frac{10^{-20}}{\sqrt{n}},\,h^{2}C_{x}(x,v)\leq\frac{10^{-10}}{n\log\frac{1}{\rho}}\text{
and }h^{2}C_{v}(x,v)\leq\frac{10^{-10}}{\sqrt{n\log\frac{1}{\rho}}}$,
guarantee that

$hC_{x}(x,v)\leq\frac{10^{-20}}{\sqrt{n}},\,h^{2}C_{x}(x,v)\leq
10^{-10}\min\left(1,\frac{\bar{\ell}_{0}}{\bar{\ell}_{1}},\frac{1}{n+\sqrt{M_{1}}+\sqrt{M_{1}^{*}}}\right),\
h^{2}C_{v}(x,v)\leq\frac{10^{-10}}{\sqrt{n+\sqrt{M_{1}}}}.$

As the isoperimetry is lower bounded by $\frac{1}{\sqrt{m}}$, Theorem 24
results in the mixing time of
$T=O\left(mh^{-2}\log\frac{\Lambda}{\varepsilon}\right)$ that ensures
$d_{\text{TV}}(\pi_{T},\pi)\leq\varepsilon$. ∎

By setting $C_{x},C_{v}$ to $0$, we can obtain the mixing time of the ideal
RHMC for exponential densities in polytopes.

See 3

####  8.4.1 Implicit midpoint method

In the polytope setting, we can explicitly compute $C_{x}(x,v)$ and
$C_{v}(x,v)$ of IMM in terms of $n$ and $\rho$.

######  Lemma 50.

For $x\in\mathbb{\mathcal{M}}_{\rho}$ and $v\in V_{\text{good}}^{x}$, let $h$
be step size of IMM with $h^{2}n\log\frac{1}{\rho}\leq 10^{-11}$. Then

$C_{x}(x,v)=O\left(n\log\frac{1}{\rho}\right),\
C_{v}(x,v)=O\left(n^{3/2}\log^{3/2}\frac{1}{\rho}\right).$

###### Proof.

By Lemma 33-2, it follows that

$\displaystyle C_{x}(x,v)$ $\displaystyle=O\left(n+\sqrt{M_{1}}\right)\lesssim
n+n\log\frac{1}{\rho}=O\left(n\log\frac{1}{\rho}\right).$

For $C_{v}(x,v)$, we first note that $M_{2}^{*}=0$ due to $\nabla^{2}f(x)=0$.
Thus by Lemma 33-3, we have

$\displaystyle C_{v}(x,v)$
$\displaystyle\lesssim\left(n+\sqrt{M_{1}}\right)^{3/2}=O\left(n^{3/2}\log^{3/2}\frac{1}{\rho}\right).$

∎

We can also specify a sufficient condition on the step size for the
sensitivity of IMM in the polytope setting.

######  Lemma 51.

For $x\in\mathbb{\mathcal{M}}_{\rho}$, $v\in V_{\text{good}}^{x}$ and step
size $h$ with $h^{2}n^{2}\log\frac{1}{\rho}\leq 10^{-10}$, IMM is sensitive at
$(x,v)$.

###### Proof.

Note that $\log\det g(x)$ is convex in $\mathbb{\mathcal{M}}$, since the
volumetric barrier defined by $\log\det\nabla^{2}\phi(x)$ is convex in $x$
(Lemma 1~3 in [Vai96]). Thus, the claim follows from Lemma 34. ∎

Substituting the estimates of $C_{x}(x,v)$ and $C_{v}(x,v)$ as well as the
sufficient condition for the sensitivity to Theorem 2, we prove that the
mixing rate of RHMC discretized by IMM for an exponential density in a
polytope is independent of the condition number and
$\left\|\alpha\right\|_{2}$.

See 4

###### Proof.

We can check that the step size
$h=O\left(\frac{1}{n^{3/2}\log\frac{\Lambda}{\varepsilon}}\right)$ satisfies
all the conditions in Theorem 2. Hence, it suffices to choose
$T=O\left(mn^{3}\log^{3}\frac{\Lambda}{\varepsilon}\right)$ to obtain
$d_{\text{TV}}(\pi_{T},\pi)\leq\varepsilon$. ∎

####  8.4.2 Generalized Leapfrog method

We now compute the mixing rate of RHMC discretized by LM. For LM, we have the
same results on $C_{x}(x,v)$ and $C_{v}(x,v)$ as IMM.

######  Lemma 52.

For $x\in\mathbb{\mathcal{M}}_{\rho}$ and $v\in V_{\text{good}}^{x}$, let $h$
be step size of LM with $h^{2}n\log\frac{1}{\rho}\leq 10^{-10}$. Then

$C_{x}(x,v)=O\left(n\log\frac{1}{\rho}\right),\
C_{v}(x,v)=O\left(n^{3/2}\log^{3/2}\frac{1}{\rho}\right).$

For the sensitivity, LM requires a slightly stronger condition on step size
compared to IMM, which follows from Lemma 36.

######  Lemma 53.

For $x\in\mathbb{\mathcal{M}}_{\rho}$, $v\in V_{\text{good}}^{x}$ and step
size $h$ with $h^{2}n^{3}\log\frac{1}{\rho}\leq 10^{-20}$, LM is sensitive at
$(x,v)$.

We prove that the mixing rate of RHMC discretized by LM for an exponential
density in a polytope with $m$ constraints is also independent of the
condition number.

See 5

###### Proof.

For step size
$h=O\left(\frac{1}{n^{3/2}\log\frac{\Lambda}{\varepsilon}}\right)$, LM is
sensitive in $\mathbb{\mathcal{M}}_{\rho}\times V_{1}^{c}$ by Lemma 36, and we
can use the estimates of $C_{x}$ and $C_{v}$ proven in Lemma 52. Thus, this
step size satisfies all the conditions in Theorem 2. Hence, it suffices to
choose $T=O\left(mn^{3}\log^{3}\frac{\Lambda}{\varepsilon}\right)$ to obtain
$d_{\text{TV}}(\pi_{T},\pi)\leq\varepsilon$. ∎

##### Acknowledgement.

This work was supported in part by NSF awards CCF-1909756, CCF-2007443 and
CCF-2134105.

## References

  * [AC21] Kwangjun Ahn and Sinho Chewi.  Efficient constrained sampling via the mirror-Langevin algorithm.  Advances in Neural Information Processing Systems (NeurIPS), 34:28405–28418, 2021.
  * [BCJ+19] Eli Bingham, Jonathan P. Chen, Martin Jankowiak, Fritz Obermeyer, Neeraj Pradhan, Theofanis Karaletsos, Rohit Singh, Paul A. Szerlip, Paul Horsfall, and Noah D. Goodman.  Pyro: Deep Universal Probabilistic Programming.  The Journal of Machine Learning Research (JMLR), 20:28:1–28:6, 2019\.
  * [CBMR19] Adam D Cobb, Atılım Güneş Baydin, Andrew Markham, and Stephen J Roberts.  Introducing an explicit symplectic integration scheme for Riemannian manifold Hamiltonian Monte Carlo.  arXiv preprint arXiv:1910.06243, 2019.
  * [CCBJ18] Xiang Cheng, Niladri S Chatterji, Peter L Bartlett, and Michael I Jordan.  Underdamped Langevin MCMC: a non-asymptotic analysis.  In Conference on Learning Theory (COLT), pages 300–323. PMLR, 2018\.
  * [CDWY20] Yuansi Chen, Raaz Dwivedi, Martin J Wainwright, and Bin Yu.  Fast mixing of Metropolized Hamiltonian Monte Carlo: Benefits of multi-step gradients.  The Journal of Machine Learning Research (JMLR), 21:92–1, 2020\.
  * [CEL+22] Sinho Chewi, Murat A Erdogdu, Mufan Li, Ruoqi Shen, and Shunshi Zhang.  Analysis of Langevin Monte Carlo from Poincaré to Log-Sobolev.  In Conference on Learning Theory (COLT), pages 1–2. PMLR, 2022\.
  * [CLA+21] Sinho Chewi, Chen Lu, Kwangjun Ahn, Xiang Cheng, Thibaut Le Gouic, and Philippe Rigollet.  Optimal dimension dependence of the Metropolis-adjusted Langevin algorithm.  In Conference on Learning Theory (COLT), pages 1260–1300. PMLR, 2021.
  * [CV16] Ben Cousins and Santosh Vempala.  A practical volume algorithm.  Mathematical Programming Computation, 8(2):133–160, 2016.
  * [CV22] Zongchen Chen and Santosh S Vempala.  Optimal Convergence Rate of Hamiltonian Monte Carlo for Strongly Logconcave Distributions.  Theory of Computing, 18(1):1–18, 2022.
  * [Dal17] Arnak Dalalyan.  Further and stronger analogy between sampling and optimization: Langevin Monte Carlo and gradient descent.  In Conference on Learning Theory (COLT), pages 678–689. PMLR, 2017\.
  * [DCWY18] Raaz Dwivedi, Yuansi Chen, Martin J Wainwright, and Bin Yu.  Log-concave sampling: Metropolis-Hastings algorithms are fast!  In Conference on Learning Theory (COLT), pages 793–797. PMLR, 2018\.
  * [DMM19] Alain Durmus, Szymon Majewski, and Błażej Miasojedow.  Analysis of Langevin Monte Carlo via convex optimization.  The Journal of Machine Learning Research (JMLR), 20(1):2666–2711, 2019.
  * [GC11] Mark Girolami and Ben Calderhead.  Riemann manifold Langevin and Hamiltonian Monte Carlo methods.  Journal of the Royal Statistical Society: Series B (Statistical Methodology), 73(2):123–214, 2011.
  * [GV22] Khashayar Gatmiry and Santosh S Vempala.  Convergence of the Riemannian Langevin Algorithm.  arXiv preprint arXiv:2204.10818, 2022.
  * [HHIL06] Ernst Hairer, Marlis Hochbruck, Arieh Iserles, and Christian Lubich.  Geometric numerical integration.  Oberwolfach Reports, 3(1):805–882, 2006.
  * [Jia21] Qijia Jiang.  Mirror Langevin Monte Carlo: the case under isoperimetry.  Advances in Neural Information Processing Systems, 34:715–725, 2021\.
  * [JLLV21] He Jia, Aditi Laddha, Yin Tat Lee, and Santosh Vempala.  Reducing isotropy and volume to KLS: an $O^{*}(n^{3}\psi^{2})$ volume algorithm.  In Proceedings of the 53rd Annual ACM SIGACT Symposium on Theory of Computing (STOC), pages 961–974, 2021.
  * [KLD+16] Zachary A King, Justin Lu, Andreas Dräger, Philip Miller, Stephen Federowicz, Joshua A Lerman, Ali Ebrahim, Bernhard O Palsson, and Nathan E Lewis.  BiGG Models: A platform for integrating, standardizing and sharing genome-scale models.  Nucleic acids research, 44(D1):D515–D522, 2016.
  * [KLS97] Ravi Kannan, László Lovász, and Miklós Simonovits.  Random walks and an $O^{*}(n^{5})$ volume algorithm for convex bodies.  Random Structures & Algorithms, 11(1):1–50, 1997.
  * [KLSV22] Yunbum Kook, YinTat Lee, Ruoqi Shen, and Santosh Vempala.  Sampling with Riemannian Hamiltonian Monte Carlo in a Constrained Space.  In Advances in Neural Information Processing Systems, 2022.
  * [KN12] Ravindran Kannan and Hariharan Narayanan.  Random walks on polytopes and an affine interior point method for linear programming.  Mathematics of Operations Research, 37(1):1–20, 2012.
  * [KV06] Adam Tauman Kalai and Santosh Vempala.  Simulated annealing for convex optimization.  Mathematics of Operations Research, 31(2):253–266, 2006.
  * [LLV20] Aditi Laddha, Yin Tat Lee, and Santosh Vempala.  Strong self-concordance and sampling.  In Proceedings of the 52nd Annual ACM SIGACT Symposium on Theory of Computing (STOC), pages 1212–1222, 2020.
  * [LNP12] Nathan E Lewis, Harish Nagarajan, and Bernhard O Palsson.  Constraining the metabolic genotype–phenotype relationship using a phylogeny of in silico methods.  Nature Reviews Microbiology, 10(4):291–305, 2012.
  * [Lov99] László Lovász.  Hit-and-run mixes fast.  Mathematical programming, 86:443–461, 1999.
  * [LS93] László Lovász and Miklós Simonovits.  Random walks in a convex body and an improved volume algorithm.  Random structures & algorithms, 4(4):359–412, 1993.
  * [LS14] Yin Tat Lee and Aaron Sidford.  Path finding methods for linear programming: Solving linear programs in $O(\sqrt{\text{rank}})$ iterations and faster algorithms for maximum flow.  In 2014 IEEE 55th Annual Symposium on Foundations of Computer Science (FOCS), pages 424–433. IEEE, 2014.
  * [LST20] Yin Tat Lee, Ruoqi Shen, and Kevin Tian.  Logsmooth gradient concentration and tighter runtimes for metropolized Hamiltonian Monte Carlo.  In Conference on Learning Theory (COLT), pages 2565–2597. PMLR, 2020.
  * [LTVW22] Ruilin Li, Molei Tao, Santosh S Vempala, and Andre Wibisono.  The mirror Langevin algorithm converges with vanishing bias.  In International Conference on Algorithmic Learning Theory (ALT), pages 718–742. PMLR, 2022.
  * [LV07] László Lovász and Santosh Vempala.  The geometry of logconcave functions and sampling algorithms.  Random Structures & Algorithms, 30(3):307–358, 2007.
  * [LV17] Yin Tat Lee and Santosh S Vempala.  Geodesic walks in polytopes.  In Proceedings of the 49th Annual ACM SIGACT Symposium on theory of Computing (STOC), pages 927–940, 2017.
  * [LV18] Yin Tat Lee and Santosh S Vempala.  Convergence rate of Riemannian Hamiltonian Monte Carlo and faster polytope volume computation.  In Proceedings of the 50th Annual ACM SIGACT Symposium on Theory of Computing (STOC), pages 1115–1121, 2018.
  * [LV22] Yin Tat Lee and Santosh S Vempala.  The manifold joys of sampling.  In 49th International Colloquium on Automata, Languages, and Programming (ICALP). Schloss Dagstuhl-Leibniz-Zentrum für Informatik, 2022\.
  * [Nes03] Yurii Nesterov.  Introductory lectures on convex optimization: A basic course, volume 87.  Springer Science & Business Media, 2003.
  * [Ngu13] Van Hoang Nguyen.  Dimensional variance inequalities of Brascamp-Lieb type and a local approach to dimensional Prékopa’s theorem.  arXiv preprint arXiv:1302.4589, 2013.
  * [NN94] Yurii Nesterov and Arkadii Nemirovskii.  Interior-point polynomial algorithms in convex programming.  SIAM, 1994.
  * [NT+02] Yurii E Nesterov, Michael J Todd, et al.  On the Riemannian geometry defined by self-concordant barriers and interior-point methods.  Foundations of Computational Mathematics, 2(4):333–361, 2002.
  * [SL19] Ruoqi Shen and Yin Tat Lee.  The randomized midpoint method for log-concave sampling.  Advances in Neural Information Processing Systems (NeurIPS), 32, 2019.
  * [Smi84] Robert L Smith.  Efficient Monte Carlo procedures for generating points uniformly distributed over bounded regions.  Operations Research, 32(6):1296–1308, 1984.
  * [Sta20] Stan Development Team.  RStan: the R interface to Stan, 2020.  R package version 2.21.2.
  * [SWF16] John Salvatier, Thomas V Wiecki, and Christopher Fonnesbeck.  Probabilistic programming in Python using PyMC3.  PeerJ Computer Science, 2:e55, 2016.
  * [TSF+13] Ines Thiele, Neil Swainston, Ronan MT Fleming, Andreas Hoppe, Swagatika Sahoo, Maike K Aurich, Hulda Haraldsdottir, Monica L Mo, Ottar Rolfsson, Miranda D Stobbe, et al.  A community-driven global reconstruction of human metabolism.  Nature biotechnology, 31(5):419–425, 2013.
  * [Vai96] Pravin M Vaidya.  A new algorithm for minimizing convex functions over convex sets.  Mathematical programming, 73(3):291–341, 1996.
  * [Vem05] Santosh Vempala.  Geometric random walks: a survey.  Combinatorial and computational geometry, 52(573-612):2, 2005.
  * [VW19] Santosh Vempala and Andre Wibisono.  Rapid convergence of the unadjusted Langevin algorithm: Isoperimetry suffices.  Advances in neural information processing systems (NeurIPS), 32, 2019.
  * [ZPFP20] Kelvin Shuangjian Zhang, Gabriel Peyré, Jalal Fadili, and Marcelo Pereyra.  Wasserstein control of mirror Langevin Monte Carlo.  In Conference on Learning Theory (COLT), pages 3814–3841. PMLR, 2020.

##  Appendix A Definitions

######  Definition 54 (Self-concordant barrier).

A self-concordant barrier $\phi:K\subset\mathbb{R}^{n}\to\mathbb{R}$ is a
function such that $\phi(x)\to\infty$ as $x\to\mathcal{\partial}K$ and that
$\left|Df^{3}(x)[h,h,h]\right|\leq 2\left(D^{2}f(x)[h,h]\right)^{3/2}$ for all
$x\in K$ and $h\in\mathbb{R}^{n}$. If $\left|Df^{4}(x)[h,h,h,h]\right|\leq
6\left(D^{2}f(x)[h,h]\right)^{2}$ is also satisfied for all $h$, then $\phi$
is called a highly self-concordant barrier.

######  Definition 55 (Self-concordance parameter).

For a self-concordant function $\phi$, the self-concordance parameter of
$\phi$ is the smallest non-negative real number $\nu_{\phi}$ such that

$|D\phi(x)[h]|^{2}\leq\nu_{\phi}D^{2}\phi(x)[h,h],$

where $Df(x)[h]$ is the directional derivative of $f$ along direction $h$ and
$D^{2}f(x)[h_{1},h_{2}]$ is the second-order directional derivative of $f$
along directions $h_{1}$ and $h_{2}$.

######  Definition 56 (Riemannian length and distance).

Let $\phi:\mathbb{R}^{n}\rightarrow\mathbb{R}$ be a self-concordant function.
For all $x\in\mathbb{R}^{d}$, we define the local norm induced by
$\nabla^{2}\phi(x)$ by

$\left\|h\right\|_{\nabla^{2}\phi(x)}=\sqrt{h^{\top}\nabla^{2}\phi(x)h}.$

For any smooth curve $c:[0,1]\rightarrow\mathbb{R}^{n}$, we define the length
of the curve as

$L_{\phi}(c)=\int_{0}^{1}\left\|\frac{d}{dt}c(t)\right\|_{\nabla^{2}\phi(c(t))}dt.$

For any $x,y\in\mathbb{R}^{d}$, we define the distance $d_{\phi}(x,y)$ to be
the infimum of the lengths of all piecewise smooth curves with $c(0)=x$ and
$c(1)=y$.

######  Definition 57 (Total variation distance).

For probability distributions $P$ and $Q$ supported on $K$, the total
variation distance (TV distance) is defined by

$d_{\text{TV}}(P,Q)=\sup_{A\subset K}\left(P(A)-Q(A)\right).$

##  Appendix B Lemmas

######  Lemma 58.

For $n\in\mathbb{N}$ and matrix $X\in\mathbb{R}^{2n\times 2n}$ of the form

$X=\left[\begin{array}[]{cc}C&I_{n}\\\ -C^{2}+R&-C\end{array}\right]$

with a symmetric matrix $C\in\mathbb{R}^{n\times n}$ and matrix
$R\in\mathbb{R}^{n\times n}$, we have

$\displaystyle X^{2n}$ $\displaystyle=\left[\begin{array}[]{cc}R^{n}&0\\\
R^{n}C-CR^{n}&R^{n}\end{array}\right],$ $\displaystyle X^{2n+1}$
$\displaystyle=\left[\begin{array}[]{cc}R^{n}C&R^{n}\\\
R^{n+1}-CR^{n}C&-CR^{n}\end{array}\right].$

The claim immediately follows from induction.

######  Lemma 59 ([LV18], Lemma 7).

In the Euclidean coordinate, the Hamiltonian equations in (2.1) can be
represented via the second-order ODE as follows:

$\displaystyle D_{t}\frac{dx}{dt}$ $\displaystyle=\mu(x),$
$\displaystyle\frac{dx}{dt}(0)$ $\displaystyle\sim\mathcal{N}(0,g(x)^{-1}),$

where $D_{t}$ is the covariant derivative along the Hamiltonian trajectory
$x(t)$ and $\mu(x)\stackrel{{\scriptstyle\mathrm{{\scriptscriptstyle
def}}}}{{=}}-g(x)^{-1}\nabla
f(x)-\frac{1}{2}g(x)^{-1}\mathrm{Tr}\left[g(x)^{-1}Dg(x)\right]$ .

######  Lemma 60 ([LV18], Lemma 64).

For matrix $E\in\mathbb{R}^{n\times n}$ with
$\left\|E\right\|_{2}<\frac{1}{4}$, we have

$\left|\log\det\left(I+E\right)-\mathrm{Tr}E\right|\leq\left\|E\right\|_{F}^{2}.$

######  Lemma 61 ([KV06], Lemma 4.1).

For a unit vector $c\in\mathbb{R}^{n}$, constant $T$ and convex set
$K\subset\mathbb{R}^{n}$, we have

$\mathbb{E}_{x\sim\pi}\left[c^{\top}x\right]\leq nT+\min_{x\in K}c^{\top}x,$

where $\pi$ is a probability density proportional to
$e^{-\frac{c^{\top}x}{T}}$.

######  Lemma 62 ([Ngu13], Corollary 6).

Let $\pi$ be a log-concave density proportional to $\exp(-V)$ on
$\mathbb{R}^{n}$. Then,

$\text{Var}_{x\sim\pi}\left(V(x)\right)\leq n.$

######  Lemma 63 ([LV07], Lemma 5.17).

Let $X\in\mathbb{R}^{n}$ be randomly chosen from a log-concave distribution.
Then for any $R>1$,

$\mathbf{P}\left(\left|X\right|>R\sqrt{\mathbb{E}X^{2}}\right)<e^{-R+1}.$

######  Lemma 64 ([NT+02], Lemma 3.1).

Suppose $\phi:\mathbb{R}^{n}\rightarrow\mathbb{R}$ is self-concordant and
$\mathcal{K}\subset\mathbb{R}^{n}$ is convex. For any $x,y\in\mathcal{K}$, i

  * •

If $d_{\phi}(x,y)\leq\delta-\delta^{2}<1$ for some $0<\delta<1$, then
$\left\|y-x\right\|_{\nabla^{2}\phi(x)}\leq\delta$.

  * •

If $\delta=\left\|x-y\right\|_{\nabla^{2}\phi(x)}<1$, then
$\delta-\frac{1}{2}\delta^{2}\leq d_{\phi}(x,y)\leq-\log(1-\delta)$.


* [43] C. Alexandrou, S. Bacchio, M. Constantinou, J. Finkenrath, K. Hadjiyiannakou, K. Jansen, G. Koutsou, and A. Vaquero Aviles-Casco, Nucleon axial, tensor, and scalar charges and $\sigma$-terms in lattice QCD, Phys. Rev. D 102 (2020), no. 5 054517, [arXiv:1909.00485].
  * [44] J. S. Lee, Revisiting Supernova 1987A Limits on Axion-Like-Particles, arXiv:1808.10136.
  * [45] BNL-E949 Collaboration, A. V. Artamonov et al., Study of the decay $K^{+}\to\pi^{+}\nu\bar{\nu}$ in the momentum region $140<P_{\pi}<199$ MeV/c, Phys. Rev. D 79 (2009) 092004, [arXiv:0903.0030].
  * [46] KOTO Collaboration, J. K. Ahn et al., Search for the $K_{L}\\!\to\\!\pi^{0}\nu\overline{\nu}$ and $K_{L}\\!\to\\!\pi^{0}X^{0}$ decays at the J-PARC KOTO experiment, Phys. Rev. Lett. 122 (2019), no. 2 021802, [arXiv:1810.09655].
  * [47] M. J. Dolan, F. J. Hiskens, and R. R. Volkas, Advancing Globular Cluster Constraints on the Axion-Photon Coupling, arXiv:2207.03102.
  * [48] K. Abazajian et al., CMB-S4 Science Case, Reference Design, and Project Plan, arXiv:1907.04473.
  * [49] E. Masso, F. Rota, and G. Zsembinszki, On axion thermalization in the early universe, Phys. Rev. D 66 (2002) 023004, [hep-ph/0203221].
  * [50] Planck Collaboration, N. Aghanim et al., Planck 2018 results. VI. Cosmological parameters, Astron. Astrophys. 641 (2020) A6, [arXiv:1807.06209]. [Erratum: Astron.Astrophys. 652, C4 (2021)].
  * [51] Simons Observatory Collaboration, P. Ade et al., The Simons Observatory: Science goals and forecasts, JCAP 02 (2019) 056, [arXiv:1808.07445].
  * [52] C. Dvorkin et al., The Physics of Light Relics, in 2022 Snowmass Summer Study, 3, 2022.  arXiv:2203.07943.
  * [53] W. L. Xu, J. B. Muñoz, and C. Dvorkin, Cosmological constraints on light but massive relics, Phys. Rev. D 105 (2022), no. 9 095029, [arXiv:2107.09664].
  * [54] L. Wolfenstein, Parametrization of the Kobayashi-Maskawa Matrix, Phys. Rev. Lett. 51 (1983) 1945.
  * [55] L. Husdal, On Effective Degrees of Freedom in the Early Universe, Galaxies 4 (2016), no. 4 78, [arXiv:1609.04979].
  * [56] K. G. Chetyrkin, J. H. Kuhn, and M. Steinhauser, RunDec: A Mathematica package for running and decoupling of the strong coupling and quark masses, Comput. Phys. Commun. 133 (2000) 43–65, [hep-ph/0004189].