% ==============================================================================
% chapters/ch05_higher_order_symplectic/sec06_dyadic_multirate_scheduling/content.tex
% Section 5.6: Dyadic Multirate Symplectic Scheduling and Skeel--Gear Synchronization
% ==============================================================================

\subsection{Dyadic Multirate Symplectic Scheduling and Skeel--Gear Synchronization}
\label{sec:dyadic_multirate_scheduling}

In gravitational $N$-body systems characterized by hierarchical clustering, high eccentricity, and large mass ratios, physical interaction timescales span multiple orders of magnitude, ranging from days for hardened inner binaries to millennia for outer spectators~\cite{aarseth2003gravitational, heggie1975binary}. Advancing the entire phase space with a global uniform time step governed by the instantaneous minimum timescale $\Delta t = \eta \min_{i<j} \tau_{ij}$ incurs a prohibitive computational penalty, scaling as $\mathcal{O}(N_{\mathrm{steps}} \cdot N^2)$. Conversely, naive asynchronous multi-rate schemes, in which particles or Jacobi nodes advance under independent uncoordinated step sizes, break the symplectic 2-form of phase space and trigger the unremovable $\mathcal{O}(\delta^{-5/2})$ Lie commutator catastrophe established in Proposition~\ref{prop:excluded_middle_dilemma} (Section~\ref{sec:synchronous_supernodes}).

To resolve this dichotomy, this section establishes the mathematical theory and algorithmic architecture of the \emph{dyadic multirate symplectic scheduler}. By restricting local sub-steps to a discrete dyadic lattice $h_k = \Delta t_{\mathrm{macro}} / 2^k$, enforcing macro-boundary synchronization in the sense of Skeel \& Gear (1992)~\cite{skeel1992does}, decoupling desingularized leaf pairs from singular apex refinement (the Leaf/Apex Decoupling Theorem), and freezing the Jacobi tree chart throughout each macroblock, we construct a multirate integrator that guarantees:
\begin{enumerate}
    \item Strict symplecticity on the extended phase space $T^*\mathbf{Q}[N] \times \mathbb{R}_t \times \mathbb{R}_E$ to machine precision IEEE 754 float64;
    \item Exact microscopic time-reversibility $\mathcal{S}(-h) = \mathcal{S}(h)^{-1}$;
    \item Absence of artificial dyadic saturation ($k = 17$) on tight binaries;
    \item Continuity and long-term boundedness of the modified shadow Hamiltonian $\widetilde{H}$.
\end{enumerate}

% ------------------------------------------------------------------------------
\subsubsection{The Discrete Dyadic Lattice and Skeel--Gear Macro-Boundary Synchronization}
\label{sub:dyadic_lattice_skeel_gear}

Let $\Delta t_{\mathrm{macro}} > 0$ denote the base macroscopic time step of the global integration scheme. To combine multi-rate flexibility with the rigid algebraic symmetries required by symplectic composition, local sub-steps are constrained to a discrete dyadic lattice:

\begin{defbox}[Dyadic Multirate Lattice]{def:conormal_dyadic_lattice}
The \emph{dyadic multirate lattice of depth $k_{\max} \in \mathbb{N}$} is the set of discrete step sizes:
\begin{equation}
    \mathbb{D}(k_{\max}) \coloneqq \left\{ h_k = \frac{\Delta t_{\mathrm{macro}}}{2^k} \;\middle|\; k \in \mathbb{N}, \; 0 \le k \le k_{\max} \right\}.
    \label{eq:dyadic_lattice_def}
\end{equation}
A macroblock over the physical interval $[t_n, t_n + \Delta t_{\mathrm{macro}}]$ assigned dyadic level $k$ is subdivided into exactly $2^k$ identical substeps:
\begin{equation}
    t_n^{(m)} = t_n + m \, h_k, \qquad m \in \{0, 1, \dots, 2^k\}.
    \label{eq:dyadic_substep_timestamps}
\end{equation}
By construction, the dyadic partition satisfies the exact completeness identity:
\begin{equation}
    \sum_{m=1}^{2^k} h_k = 2^k \cdot \left( \frac{\Delta t_{\mathrm{macro}}}{2^k} \right) = \Delta t_{\mathrm{macro}},
    \label{eq:dyadic_completeness}
\end{equation}
which evaluates bit-exactly to machine precision without cumulative floating-point summation drift by computing intermediate timestamps via $t_n + m \cdot (\Delta t_{\mathrm{macro}} / 2^k)$ rather than recursive addition $t \leftarrow t + h_k$.
\end{defbox}

The fundamental challenge of multirate symplectic integration lies in the preservation of geometric conservation laws across step size transitions. Skeel \& Gear (1992)~\cite{skeel1992does} established that variable time-stepping in Hamiltonian systems generically destroys symplectic conservation unless the step adaptation occurs synchronously at pre-assigned macro-boundaries:

\begin{propbox}[Intra-Block Invariance and Macroblock Symplectic Propagation]{prop:skeel_gear_boundary_sync}\label{thm:skeel_gear_boundary_sync}
Following the multirate Hamiltonian framework of Skeel and Gear~\cite{skeel1992does}, let $\mathcal{S}_h : T^*\mathbf{Q}[N] \to T^*\mathbf{Q}[N]$ be a symmetric palindromic symplectic integrator of order $2p$ ($\mathcal{S}_{-h} = \mathcal{S}_h^{-1}$) advancing across macroblocks $I_n = [t_n, t_{n+1}]$ of duration $\Delta t_{\mathrm{macro}}$.
If the dyadic refinement level $k_n \in \{0, \dots, k_{\max}\}$ and the Jacobi chart $\mathcal{T}_n \in \mathrm{Trees}(N)$ remain frozen throughout macroblock $I_n$, then:
\begin{enumerate}
    \item \emph{Exact Microscopic Reversibility}: The macroblock propagator $\Psi_{\Delta t_{\mathrm{macro}}}^{(n)} \coloneqq (\mathcal{S}_{h_{k_n}})^{2^{k_n}}$ is an exact symplectomorphism on $T^*\mathbf{Q}[N]$ preserving the time-reversibility identity $(\Psi_{\Delta t_{\mathrm{macro}}}^{(n)})^{-1} = \Psi_{-\Delta t_{\mathrm{macro}}}^{(n)}$.
    \item \emph{Suppression of Internal Commutator Spikes}: Freezing chart $\mathcal{T}_n$ inside $I_n$ eliminates mid-step Lie commutator jumps $\Delta H_{\mathrm{flip}} \sim h_{k_n} \{W_{\mathrm{apex}}^{\mathrm{screen},(A)}, W_{\mathrm{apex}}^{\mathrm{screen},(B)}\} \ne 0$. Chart transitions are executed exclusively at macro-boundaries $t_n$ via the covariant gauge operator $T_{A \to B}^{(\mathrm{PBM4})}$ (Theorem~\ref{thm:pbm4_gauge_operator_exact}).
\end{enumerate}
\end{propbox}

\begin{proofbox}[Proof of Proposition~\ref{prop:skeel_gear_boundary_sync}]
Symplecticity follows from the group property of canonical transformations: the $2^{k_n}$-fold composition of symplectic maps $\mathcal{S}_{h_{k_n}}$ is strictly symplectic. Microscopic reversibility follows from $(\mathcal{S}_{h_{k_n}})^{-2^{k_n}} = (\mathcal{S}_{-h_{k_n}})^{2^{k_n}} = \Psi_{-\Delta t_{\mathrm{macro}}}^{(n)}$ by palindromic symmetry~\cite{skeel1992does}. If a chart change occurs at an interior substep $m \in \{1, \dots, 2^{k_n}-1\}$, the broken symmetry introduces an uncompensated Baker--Campbell--Hausdorff cross-bracket $\frac{1}{2} h_{k_n}^2 [X_{V_A}, X_{V_B}] \ne 0$, destroying shadow energy preservation. Restricting transitions to macro-boundaries restores global reversibility.
\end{proofbox}

% ------------------------------------------------------------------------------
\subsubsection{Scale-Invariant Dynamical Timescales and Predictive Encounter Anticipation}
\label{sub:conormal_timescales_anticipation}

Having established the stability of macro-boundary synchronization, we now formalize the physical observables governing the selection of the dyadic exponent $k \in \{0, \dots, k_{\max}\}$. Rather than relying on arbitrary metric distance tables (such as empirical thresholds $r < 0.015$), the dyadic scheduler dynamically evaluates three scale-invariant, Galilean-invariant characteristic timescales for every interacting pair $(i, j)$:

\begin{defbox}[Scale-Invariant Interaction Timescales]{def:conormal_timescales}
Let $\mathbf{q}_i, \mathbf{q}_j \in \mathbb{R}^3$ and $\mathbf{v}_i, \mathbf{v}_j \in \mathbb{R}^3$ be the Cartesian positions and velocities of bodies $i$ and $j$, with relative separation vector $\mathbf{r}_{ij} = \mathbf{q}_i - \mathbf{q}_j$, distance $r_{ij} = \|\mathbf{r}_{ij}\|$, relative velocity $\mathbf{v}_{\mathrm{rel}} = \mathbf{v}_i - \mathbf{v}_j$, and radial velocity:
\begin{equation}
    v_{\mathrm{rad}} = \frac{\mathbf{r}_{ij} \cdot \mathbf{v}_{\mathrm{rel}}}{r_{ij}}.
    \label{eq:radial_velocity_def}
\end{equation}
We define:
\begin{enumerate}
    \item \emph{Kinematic Fly-by Timescale}:
    \begin{equation}
        \tau_{\mathrm{kin}}(i, j) \coloneqq \frac{r_{ij}}{\max(\|\mathbf{v}_{\mathrm{rel}}\|, \varepsilon_v)},
        \label{eq:tau_kin_def}
    \end{equation}
    measuring the transit time across the local interaction sphere at constant relative velocity.
    \item \emph{Dynamical Free-Fall Timescale}:
    \begin{equation}
        \tau_{\mathrm{dyn}}(i, j) \coloneqq \sqrt{\frac{r_{ij}^3}{G (m_i + m_j)}},
        \label{eq:tau_dyn_def}
    \end{equation}
    governing the local gravitational orbital frequency and free-fall collapse time.
    \item \emph{Predictive Bilateral Pericenter Encounter Timescale}: For pairs whose osculating encounter envelope satisfies the symmetric macro-window condition $|t_{\mathrm{peri}}| \le 2 |\Delta t_{\mathrm{macro}}|$ (where $t_{\mathrm{peri}} \coloneqq - (\mathbf{r}_{ij} \cdot \mathbf{v}_{\mathrm{rel}}) / \|\mathbf{v}_{\mathrm{rel}}\|^2$ in the rectilinear flyby approximation, or $|M| \le n \Delta t_{\mathrm{macro}}$ for eccentric orbits), the two-body Keplerian energy and specific angular momentum:
    \begin{equation}
        E_{2b} = \frac{1}{2} \|\mathbf{v}_{\mathrm{rel}}\|^2 - \frac{G(m_i + m_j)}{r_{ij}}, \qquad \mathbf{L} = \mathbf{r}_{ij} \times \mathbf{v}_{\mathrm{rel}},
        \label{eq:kepler_2b_invariants}
    \end{equation}
    predict the osculating pericenter distance $r_{\mathrm{peri}} = a(1-e)$ and pericenter velocity $v_{\mathrm{peri}} = \sqrt{\|\mathbf{v}_{\mathrm{rel}}\|^2 + 2 G(m_i + m_j)(r_{\mathrm{peri}}^{-1} - r_{ij}^{-1})}$. The predictive encounter timescale is:
    \begin{equation}
        \tau_{\mathrm{peri}}(i, j) \coloneqq \min\left( \frac{r_{\mathrm{peri}}}{v_{\mathrm{peri}}}, \sqrt{\frac{r_{\mathrm{peri}}^3}{G(m_i + m_j)}} \right).
        \label{eq:tau_peri_def}
    \end{equation}
    Because $|t_{\mathrm{peri}}|$ is an even function under phase-space time reversal $\mathcal{R}(\mathbf{q}, \mathbf{p}) = (\mathbf{q}, -\mathbf{p})$, the bilateral activation window ensures identical refinement sequences on approach and departure, strictly preserving microscopic time-reversibility.
\end{enumerate}
\end{defbox}

\begin{propbox}[Scale and Galilean Invariance of the Dyadic Scheduler]{prop:dyadic_scale_invariance}
Under the gravitational homothety group of Sundman~\cite{sundman1912memoire} and McGehee~\cite{mcgehee1974triple}:
\begin{equation}
    \mathbf{q} \longmapsto \lambda \mathbf{q}, \qquad t \longmapsto \lambda^{3/2} t, \qquad \mathbf{v} \longmapsto \lambda^{-1/2} \mathbf{v} \quad (\lambda > 0),
    \label{eq:gravitational_homothety_action}
\end{equation}
all three timescales scale with identical homogeneity:
\begin{equation}
    \tau_{\mathrm{kin}} \longmapsto \lambda^{3/2} \tau_{\mathrm{kin}}, \qquad \tau_{\mathrm{dyn}} \longmapsto \lambda^{3/2} \tau_{\mathrm{dyn}}, \qquad \tau_{\mathrm{peri}} \longmapsto \lambda^{3/2} \tau_{\mathrm{peri}}.
    \label{eq:timescale_scaling}
\end{equation}
Consequently, for any safety factor $\eta_{\mathrm{safety}} \in (0, 1)$, the dimensionless ratio:
\begin{equation}
    \mathcal{R} \coloneqq \frac{\Delta t_{\mathrm{macro}}}{\eta_{\mathrm{safety}} \min_{i<j} \tau(i, j)}
    \label{eq:dimensionless_refinement_ratio}
\end{equation}
is strictly invariant under homotheties ($\mathcal{R}(\lambda \mathbf{q}) = \mathcal{R}(\mathbf{q})$). The resulting continuous dyadic level:
\begin{equation}
    k = \min\left( k_{\max}, \max\left( k_{\mathrm{default}}, \left\lceil \log_2 \mathcal{R} \right\rceil \right) \right)
    \label{eq:continuous_dyadic_level_formula}
\end{equation}
is an absolute conformal invariant, preserving numerical resolution identically from compact stellar triples to galactic clusters.
\end{propbox}

\begin{proofbox}[Proof of Proposition~\ref{prop:dyadic_scale_invariance}]
By direct substitution of the homothety transformations:
$\tau_{\mathrm{kin}} = (\lambda r) / (\lambda^{-1/2} v) = \lambda^{3/2} \tau_{\mathrm{kin}}$, and $\tau_{\mathrm{dyn}} = \sqrt{(\lambda r)^3 / GM} = \lambda^{3/2} \tau_{\mathrm{dyn}}$. In the predictive pericenter formula, $E_{2b} \mapsto \lambda^{-1} E_{2b}$ and $L \mapsto \lambda^{1/2} L$, yielding eccentricity $e = \sqrt{1 + 2 E L^2 / (GM)^2} \mapsto e$ (strictly invariant) and semi-major axis $a = -GM / (2E) \mapsto \lambda a$. Thus $r_{\mathrm{peri}} \mapsto \lambda r_{\mathrm{peri}}$ and $v_{\mathrm{peri}} \mapsto \lambda^{-1/2} v_{\mathrm{peri}}$, confirming $\tau_{\mathrm{peri}} \mapsto \lambda^{3/2} \tau_{\mathrm{peri}}$. Since $\Delta t_{\mathrm{macro}} \mapsto \lambda^{3/2} \Delta t_{\mathrm{macro}}$, the ratio $\mathcal{R}$ is the quotient of two quantities scaling as $\lambda^{3/2}$, hence $\mathcal{R} \mapsto \mathcal{R}$.
\end{proofbox}

% ------------------------------------------------------------------------------
\subsubsection{The Leaf/Apex Decoupling Theorem and Singular Pericenter Isolation}
\label{sub:leaf_apex_decoupling}

A naive application of predictive pericenter refinement ($\tau_{\mathrm{peri}} \to 0$ as $r_{\mathrm{peri}} \to 0$) to \emph{all} particle pairs in an $N$-body system introduces a catastrophic computational pathology: for any tightly bound, eccentric Keplerian binary (e.g., $e \ge 0.999$, $r_{\mathrm{peri}} \sim 10^{-5}$), the scheduler would saturate at the maximum depth $k_{\max} = 17$ ($2^{17} = 131\,072$ substeps per macro-step), even though the binary is already integrated analytically to machine precision!

This obstruction is resolved by the \emph{Leaf/Apex Decoupling Theorem}:

\begin{thmbox}[Leaf/Apex Decoupling Theorem]{thm:leaf_apex_decoupling}
Under the Fulton--MacPherson desingularization of terminal pairs (Fulton and MacPherson~\cite{fulton1994compactification}), let $\mathcal{T} \in \mathrm{Trees}(N)$ be an adapted binary Jacobi tree, and let $\operatorname{LeafPairs}(\mathcal{T}) \subset \binom{\{1, \dots, N\}}{2}$ be the set of mutual two-body binary pairs occupying the terminal leaf nodes of $\mathcal{T}$.
\begin{enumerate}
    \item \emph{Exact Analytical Desingularization and Invariant Timescale}: For every leaf pair $C = (i, j) \in \operatorname{LeafPairs}(\mathcal{T})$ with negative binding energy $E_C = \frac{1}{2\mu_C} \|\boldsymbol{\pi}_C\|^2 - \frac{G m_i m_j}{\|\mathbf{r}_C\|} < 0$, the internal Keplerian Hamiltonian:
    \begin{equation}
        H_C(\mathbf{r}_C, \boldsymbol{\pi}_C) = \frac{\|\boldsymbol{\pi}_C\|^2}{2\mu_C} - \frac{G m_i m_j}{\|\mathbf{r}_C\|}
        \label{eq:leaf_pair_hamiltonian}
    \end{equation}
    is integrated analytically by the universal Stumpff--Sundman diffeomorphism (Section~\ref{sec:suzuki_barrier_mcgehee}). Consequently, the local flow $\exp(h X_{H_C})$ generates \emph{zero operator splitting error} and zero truncation error across arbitrarily close pericenters ($r_{\mathrm{peri}} \to 0$, $e \to 1$). The dynamical timescale of the bound leaf pair is governed by its invariant semi-major axis $a_C = - G m_i m_j / (2 E_C)$:
    \begin{equation}
        \tau_{\mathrm{dyn}}(C) = \sqrt{\frac{a_C^3}{G(m_i + m_j)}} = \frac{P_{\mathrm{orb}}(C)}{2\pi},
    \end{equation}
    which remains strictly bounded and non-vanishing throughout eccentric pericenter passages.
    \item \emph{Hill Isolation Guard and Two-Regime Quadrupolar Tidal Decoupling}: Let $M_C = m_i + m_j$ and $\mathbf{R}_{\mathrm{cm},C} = (m_i \mathbf{q}_i + m_j \mathbf{q}_j) / M_C$. For $N = 2$, the apex remainder vanishes identically ($W_{\mathrm{apex}}^{\mathrm{screen}} \equiv 0$) and $\tau_{\mathrm{pair}}(i, j) \coloneqq \tau_{\mathrm{dyn}}(C)$. For $N \ge 3$, under the instantaneous Hill sphere isolation condition for all external perturbers $k \notin \{i, j\}$ at distance $R_{\mathrm{ext},k} = \|\mathbf{q}_k - \mathbf{R}_{\mathrm{cm},C}\|$:
    \begin{equation}
        r_{ij} \le \frac{1}{3} R_{\mathrm{Hill}}(C, k) \coloneqq \frac{1}{3} R_{\mathrm{ext},k} \left( \frac{M_C}{3(M_C + m_k)} \right)^{1/3} \qquad \text{and} \qquad R_{\mathrm{ext},k} \ge 3 r_{ij}, \qquad \forall k \notin \{i, j\},
    \end{equation}
    the external tidal perturbation $W_{\mathrm{apex}}^{\mathrm{screen}}$ on pair $C$ is strictly quadrupolar ($\mathcal{O}(r_C / R_{\mathrm{ext},k}^3)$) and uniformly bounded across pericenter. Defining the Hill-clamped effective semi-major axis $a_{\mathrm{eff}}(C) \coloneqq \min\!\bigl(a_C,\, \min_{k \notin \{i, j\}} \tfrac{1}{3} R_{\mathrm{Hill}}(C, k)\bigr)$, the characteristic orbital size $\ell_C \coloneqq \max\bigl(r_{ij}, a_{\mathrm{eff}}(C)\bigr)$, the effective orbital timescale $\tau_{\mathrm{dyn,eff}}(C) \coloneqq \sqrt{a_{\mathrm{eff}}(C)^3 / (G M_C)}$, the local encounter timescale $\tau_{\mathrm{loc}}(i, j) \coloneqq \min\bigl(\tau_{\mathrm{kin}}(i, j), \tau_{\mathrm{dyn}}(i, j), \tau_{\mathrm{peri}}(i, j)\bigr)$ at effective distance $r_{\mathrm{loc}}(i, j) \coloneqq \min(r_{ij}, r_{\mathrm{eff}})$, and the dimensionless quadrupolar tidal parameter:
    \begin{equation}
        \varepsilon_{\mathrm{tidal}}(C) \coloneqq \max_{k \notin \{i, j\}} \frac{m_k}{M_C} \left( \frac{\ell_C}{R_{\mathrm{ext},k}} \right)^3,
    \end{equation}
    the pairwise scheduling timescale is governed by the two-regime tidal decoupling law:
    \begin{equation}
        \tau_{\mathrm{pair}}(i, j) \coloneqq \begin{cases} 
            \tau_{\mathrm{dyn}}(C), & \begin{aligned}[t] &\text{if } C = (i, j) \in \operatorname{LeafPairs}(\mathcal{T}),\; E_C < 0, \\[-2pt] &\text{and } \bigl(N = 2 \text{ or } \varepsilon_{\mathrm{tidal}}(C) \le 10^{-30}\bigr), \end{aligned} \\[6pt]
            \min\!\Bigl( \tau_{\mathrm{dyn}}(C),\; \max\bigl( \tau_{\mathrm{reg},1}(i, j),\, \tau_{\mathrm{reg},2}(C) \bigr) \Bigr), & \begin{aligned}[t] &\text{if } C = (i, j) \in \operatorname{LeafPairs}(\mathcal{T}),\; E_C < 0, \\[-2pt] &N \ge 3, \text{ and Hill-isolated}, \end{aligned} \\[6pt]
            \tau_{\mathrm{loc}}(i, j) = \min\left( \tau_{\mathrm{kin}}(i, j), \tau_{\mathrm{dyn}}(i, j), \tau_{\mathrm{peri}}(i, j) \right), & \text{otherwise},
        \end{cases}
        \label{eq:leaf_apex_decoupling_rule}
    \end{equation}
    where the instantaneous quadrupolar scale $\tau_{\mathrm{reg},1}$ and the periastron-invariant tidal floor $\tau_{\mathrm{reg},2}$ are given by:
    \begin{equation}
        \tau_{\mathrm{reg},1}(i, j) \coloneqq \frac{\tau_{\mathrm{loc}}(i, j)}{\Bigl( \varepsilon_{\mathrm{tidal}}(C) \bigl( \min(r_{\mathrm{loc}}, \ell_C) / \ell_C \bigr)^2 \Bigr)^{1/4}}, \qquad
        \tau_{\mathrm{reg},2}(C) \coloneqq \frac{0.1 \, \eta_{\mathrm{safety}}^3}{\varepsilon_{\mathrm{tidal}}(C)^{3/4}} \, \tau_{\mathrm{dyn,eff}}(C).
        \label{eq:two_regime_tidal_modulator}
    \end{equation}
    \item \emph{Protection Against Algorithmic Paralysis}: Because $\tau_{\mathrm{reg},2}(C)$ depends exclusively on $a_{\mathrm{eff}}(C)$ and $\varepsilon_{\mathrm{tidal}}(C)$ and is strictly independent of $r_{ij}$ and $r_{\mathrm{peri}}$, as $r_{ij} \to r_{\mathrm{peri}} \to 0$ ($e \to 1$) the vanishing branch $\tau_{\mathrm{reg},1} = \mathcal{O}(r_{\mathrm{loc}}) \to 0$ is superseded in $\max(\tau_{\mathrm{reg},1}, \tau_{\mathrm{reg},2})$ by the strictly positive floor $\tau_{\mathrm{reg},2}(C) > 0$. Under rule~\eqref{eq:leaf_apex_decoupling_rule}, active leaf pairs require only modest orbital resolution ($k \le 7$, $h_k / P_{\mathrm{orb}} \sim 0.05$), while extreme dyadic refinement ($k \to 17$) is deployed exclusively when an external body plunges into the non-integrable screened apex potential $W_{\mathrm{apex}}^{\mathrm{screen}}$ or breaches the Hill sphere.
\end{enumerate}
\end{thmbox}
\pythonrefright{py_system}{system/dyadic_scheduler.py}{Float64 / Float128}

\begin{proofbox}[Proof of Theorem~\ref{thm:leaf_apex_decoupling}]
By the Fulton--MacPherson Hamiltonian decomposition (Section~\ref{sec:screened_apex}), the global Hamiltonian partitions as:
\begin{equation}
    H = T_0(\boldsymbol{\pi}_0) + \sum_{v \in \operatorname{Comp}} \frac{\|\boldsymbol{\pi}_v\|^2}{2\mu_v} + \sum_{C \in \operatorname{Leaf}} H_C + W_{\mathrm{apex}}^{\mathrm{screen}}.
\end{equation}
The mutual Newtonian potential $-G m_i m_j / r_{ij}$ of any leaf pair $C = (i, j)$ is entirely contained within $H_C$ and is \emph{strictly absent} from the apex remainder $W_{\mathrm{apex}}^{\mathrm{screen}}$ by definition of the screened partition. In the operator splitting $\mathcal{S}_h = \exp(\frac{h}{2} X_V) \exp(h X_T) \exp(\frac{h}{2} X_V)$, the Keplerian sub-flow $\exp(h X_{H_C})$ is executed in closed form via Stumpff functions without spatial discretization. Thus, no truncation error proportional to $r_{\mathrm{peri}}^{-p}$ is generated by $H_C$. The only error involving pair $(i, j)$ arises from tidal coupling with third bodies $k \notin \{i, j\}$ through $W_{\mathrm{apex}}^{\mathrm{screen}}$.
As long as the Hill isolation guard $r_{ij} \le \frac{1}{3} R_{\mathrm{Hill}}(C, k)$ and $R_{\mathrm{ext},k} \ge 3 r_{ij}$ is respected, the tidal force exerted by third bodies scales as $\mathcal{O}(\varepsilon_{\mathrm{tidal}} r_{ij} / \ell_C)$ (Equation~\eqref{eq:quadrupole_leading_order}), which vanishes linearly as $r_{ij} \to 0$. Away from pericenter ($r_{\mathrm{loc}} \sim \ell_C$), matching the fourth-order splitting commutator $\mathcal{O}(\varepsilon_{\mathrm{tidal}} (r_{\mathrm{loc}}/\ell_C)^2 (h/\tau_{\mathrm{loc}})^4)$ yields the first regime $\tau_{\mathrm{reg},1} = \tau_{\mathrm{loc}} (\varepsilon_{\mathrm{tidal}} (r_{\mathrm{loc}}/\ell_C)^2)^{-1/4}$. As $r_{\mathrm{loc}} \to 0$ near periastron, the integrated tidal impulse over the analytical Sundman passage is bounded by the semi-major axis scale $a_{\mathrm{eff}}(C)$, yielding the periastron-invariant second regime $\tau_{\mathrm{reg},2}(C) = 0.1 \eta_{\mathrm{safety}}^3 \varepsilon_{\mathrm{tidal}}^{-3/4} \tau_{\mathrm{dyn,eff}}(C) > 0$. Taking $\min(\tau_{\mathrm{dyn}}(C), \max(\tau_{\mathrm{reg},1}, \tau_{\mathrm{reg},2}))$ therefore resolves external quadrupolar tides near apocenter while rigorously preventing artificial step-size collapse at periastron ($r_{\mathrm{peri}}(C) \to 0$). Conversely, if a third body breaches the Hill sphere ($r_{ij} > \frac{1}{3} R_{\mathrm{Hill}}(C, k)$ or $R_{\mathrm{ext},k} < 3 r_{ij}$), or for non-leaf pairs $(i, k)$ whose interaction is mediated by $W_{\mathrm{apex}}^{\mathrm{screen}}$, strong three-body tidal deformation occurs, requiring the full local and predictive encounter timescale $\tau_{\mathrm{loc}}(i, j)$ until an atlas chart rotation transfers the tightest pair to a leaf.
\end{proofbox}

% ------------------------------------------------------------------------------
\subsubsection{Shadow Hamiltonian Continuity and Phase Preservation}
\label{sub:shadow_hamiltonian_continuity}

We conclude by demonstrating that the dyadic multirate scheduler preserves the asymptotic properties of the modified (shadow) Hamiltonian across variable step sizes:

\begin{propbox}[Interfacial Shadow Energy Bounds and Quasi-Reversibility]{prop:dyadic_shadow_continuity}\label{thm:dyadic_shadow_continuity}
Following the adiabatic invariant and backward error analysis of Benettin and Giorgilli~\cite{benettin1994hamiltonian} and Hairer, Lubich and Wanner~\cite{hairer2006geometric}, let $\mathcal{S}_{h_k}$ be a symmetric symplectic scheme of order $2p$ advancing on the dyadic multirate lattice. On each macroblock $I_n = [t_n, t_{n+1}]$ of duration $\Delta t_{\mathrm{macro}}$ with level $k_n$, the numerical flow interpolates the autonomous modified Hamiltonian:
\begin{equation}
    \widetilde{H}_{k_n}(\mathbf{r}, \boldsymbol{\pi}) = H(\mathbf{r}, \boldsymbol{\pi}) + h_{k_n}^{2p} H_{2p+1}(\mathbf{r}, \boldsymbol{\pi}) + \mathcal{O}(h_{k_n}^{2p+2}).
    \label{eq:modified_hamiltonian_dyadic_block}
\end{equation}
Then:
\begin{enumerate}
    \item \emph{Interfacial Energy Jump Bound}: At the macro-boundary $t_n$ separating block $I_{n-1}$ (level $k_{n-1}$) and block $I_n$ (level $k_n$), the shadow energy discontinuity is bounded by:
    \begin{equation}
        \left| \widetilde{H}_{k_n}(\mathbf{z}(t_n)) - \widetilde{H}_{k_{n-1}}(\mathbf{z}(t_n)) \right| \le \left| 2^{-2p k_n} - 2^{-2p k_{n-1}} \right| \Delta t_{\mathrm{macro}}^{2p} \, \|H_{2p+1}(\mathbf{z}(t_n))\|.
        \label{eq:interfacial_shadow_jump}
    \end{equation}
    \item \emph{Adiabatic Confinement vs.\ Chaotic Random-Walk Diffusion}:
    \begin{itemize}
        \item For symmetric, isolated close encounters, the bilateral time-reversible selection of timescales ensures that inward and outward level progressions match identically ($k_{n_{\mathrm{peri}} + m} = k_{n_{\mathrm{peri}} - m}$). The interfacial jumps cancel telescopically to higher order ($\sum \Delta \widetilde{H} = \mathcal{O}(\Delta t_{\mathrm{macro}}^{2p+1})$), bounding physical energy fluctuations within an envelope $|\Delta E| \le C \Delta t_{\mathrm{macro}}^{2p}$ over exponentially long times.
        \item Under generic asymmetric multi-body encounters or chaotic flybys where chart transitions and level switches do not pair symmetrically, the interfacial energy jumps behave as zero-mean quasi-random perturbations. The cumulative energy error exhibits slow diffusive random-walk growth:
        \begin{equation}
            \Delta E_{\mathrm{rms}}(t) \le C \, \sqrt{\frac{t}{\Delta t_{\mathrm{macro}}}} \cdot \Delta t_{\mathrm{macro}}^{2p},
        \end{equation}
        preserving the superior $\mathcal{O}(\sqrt{t})$ Brownian envelope characteristic of symplectic integrators over the linear secular drift $\mathcal{O}(t \cdot \Delta t)$ of non-geometric schemes.
    \end{itemize}
\end{enumerate}
\end{propbox}

\begin{proofbox}[Proof of Proposition~\ref{prop:dyadic_shadow_continuity}]
Subtracting $\widetilde{H}_{k_{n-1}}$ from $\widetilde{H}_{k_n}$ at the interface point $\mathbf{z}(t_n)$ immediately yields~\eqref{eq:interfacial_shadow_jump} from the dyadic scaling $h_k = \Delta t_{\mathrm{macro}} / 2^k$. Under time-reversible orbit segments, timescales depend only on even functions of time ($r(t)$ and $\|\mathbf{v}(t)\|$), yielding identical dyadic sequences on approach and departure. The interfacial jumps telescope symmetrically, leaving only higher-order non-adiabatic errors bounded by the Benettin--Giorgilli theorem~\cite{benettin1994hamiltonian, hairer2006geometric}.
Under chaotic or asymmetric multi-body encounters, the interface jump signs and amplitudes decouple, acting as uncorrelated perturbations with variance $\sigma_{\Delta E}^2 \sim \mathcal{O}(\Delta t_{\mathrm{macro}}^{4p})$. Summing over $N_{\mathrm{steps}} = t / \Delta t_{\mathrm{macro}}$ independent steps yields root-mean-square diffusion $\Delta E_{\mathrm{rms}} = \sqrt{N_{\mathrm{steps}}} \, \sigma_{\Delta E} \sim \mathcal{O}(\sqrt{t \cdot \Delta t_{\mathrm{macro}}^{4p-1}})$, establishing the Brownian diffusion envelope.
\end{proofbox}

\noindent The complete operational architecture of the dyadic multirate scheduler is summarized in the formal algorithmic specification below:

\begin{algobox}{Dyadic Multirate Symplectic Scheduler (Skeel--Gear \& Leaf-Decoupled)}\label{alg:dyadic_multirate_scheduler}
\begin{enumerate}
    \item \textbf{Input}: Macro-step $\Delta t_{\mathrm{macro}}$, tree $\mathcal{T}$, state $\mathbf{z}_n = (\mathbf{q}_n, \mathbf{p}_n)$, limits $(k_{\mathrm{default}}, k_{\max})$, safety $\eta_{\mathrm{safety}}$.
    \item \textbf{Leaf Identification}: Extract terminal binary pairs $\operatorname{LeafPairs}(\mathcal{T}) = \{(i_c, j_c)\}$.
    \item \textbf{Pairwise Timescale Scanning}: For all $1 \le i < j \le N$:
    \begin{itemize}
        \item Compute kinematic $\tau_{\mathrm{kin}} = r_{ij} / \|\mathbf{v}_{\mathrm{rel}}\|$ and dynamical $\tau_{\mathrm{dyn}} = \sqrt{r_{ij}^3 / GM_{ij}}$, and initialize $\tau_{\mathrm{loc}} = \min(\tau_{\mathrm{kin}}, \tau_{\mathrm{dyn}})$, $r_{\mathrm{loc}} = r_{ij}$.
        \item If not ($C = (i, j) \in \operatorname{LeafPairs}(\mathcal{T})$ and $N = 2$ and $E_{2b} < 0$) and $t_{\mathrm{Kepler}} \le 2 |\Delta t_{\mathrm{macro}}|$:
        \begin{itemize}
            \item Compute effective intra-block distance $r_{\mathrm{eff}}$ and velocity $v_{\mathrm{eff}}$ via Kepler 2-body invariants, and update $\tau_{\mathrm{loc}} = \min(\tau_{\mathrm{loc}}, r_{\mathrm{eff}}/v_{\mathrm{eff}}, \sqrt{r_{\mathrm{eff}}^3/GM_{ij}})$ and $r_{\mathrm{loc}} = \min(r_{ij}, r_{\mathrm{eff}})$.
        \end{itemize}
        \item If $C = (i, j) \in \operatorname{LeafPairs}(\mathcal{T})$ with $E_{2b} < 0$ and Hill-isolated ($r_{ij} \le \frac{1}{3} R_{\mathrm{Hill}}(C, k)$ and $R_{\mathrm{ext},k} \ge 3 r_{ij}$ for all $k \notin \{i, j\}$): set $\tau(i, j) = \tau_{\mathrm{dyn}}(C) = \sqrt{a_C^3/GM_C}$ if $N = 2$ (or $\varepsilon_{\mathrm{tidal}}(C) \le 10^{-30}$), else set $\tau(i, j) = \min\bigl(\tau_{\mathrm{dyn}}(C), \max(\tau_{\mathrm{reg},1}(i, j), \tau_{\mathrm{reg},2}(C))\bigr)$ via~\eqref{eq:two_regime_tidal_modulator}.
        \item Else: set $\tau(i, j) = \tau_{\mathrm{loc}}$.
    \end{itemize}
    \item \textbf{Dyadic Exponent Determination}:
    \begin{equation*}
        \tau_{\min} = \min_{i<j} \tau(i, j), \qquad k = \min\left( k_{\max}, \max\left( k_{\mathrm{default}}, \left\lceil \log_2 \frac{\Delta t_{\mathrm{macro}}}{\eta_{\mathrm{safety}} \tau_{\min}} \right\rceil \right) \right).
    \end{equation*}
    \item \textbf{Substep Propagation}: Set $h_k = \Delta t_{\mathrm{macro}} / 2^k$. Advance state $\mathbf{z}$ across $2^k$ substeps via $\mathcal{S}_{h_k}$, holding chart $\mathcal{T}$ and step $h_k$ strictly frozen.
    \item \textbf{Macro-Boundary Closure}: At $t_{n+1} = t_n + \Delta t_{\mathrm{macro}}$, re-evaluate Jacobi tree transitions and apply gauge connection $T_{A \to B}^{(\mathrm{PBM4})}$ if an atlas rotation is triggered.
\end{enumerate}
\end{algobox}

% ------------------------------------------------------------------------------
\subsection{Conclusion of Section~\ref{ch:ch05_higher_order_symplectic} and Algorithmic Transition to Section~\ref{ch:ch06_algorithmic_architecture}}
\label{sub:ch05_conclusion_transition}

Section~\ref{ch:ch05_higher_order_symplectic} has established the complete mathematical architecture of high-order symplectic integration on the Fulton--MacPherson compactified configuration manifold $\mathbf{Q}[N]$ under zero gravitational softening ($\epsilon = 0$):
\begin{enumerate}
    \item \textbf{Geometric Splitting and Regularized Leaf Flow (\S\ref{sec:strang_propagator}--\S\ref{sec:suzuki_barrier_mcgehee}):} Starting from the palindromic 7-operator Strang splitting $\mathcal{S}_2(h)$ on $\overline{\mathcal{P}}_{\mathrm{ext}}$, we reconciled the Sheng--Suzuki negative-step barrier with unsoftened binary collisions ($r_C \to 0$) by confining retrograde sub-steps $a_i < 0$ to analytical Stumpff--Sundman leaf solvers on the McGehee collision manifold $\mathbf{M}_C$ ($\mathcal{B}_C^2 = \operatorname{Id}_{\mathbf{M}_C}$).
    \item \textbf{High-Order Blanes--Moan Propagators (\S\ref{sec:bm4_seven_stages}--\S\ref{sec:bm6_eleven_stages}):} Replacing Yoshida's 3-stage scheme with the 7-stage fourth-order composition $\mathrm{BM}_4^{(7)}$ (suppressing the order-5 Lie error norm by $>400\times$) and the 11-stage sixth-order composition $\mathrm{BM}_6^{(11)}$ established exponential Benettin--Giorgilli shadow Hamiltonian confinement on $\mathcal{D}_{\mathrm{reg}}$ and exact discrete preservation of all 10 Lie--Bargmann invariants.
    \item \textbf{Universal Obstruction and Cotangent-Lifted Boundary Processing (\S\ref{sec:symplectic_correctors}):} Using the bi-graded Poisson--Lie filtration $\mathcal{L}_{m,k}$ (Proposition~\ref{prop:bigraded_bch_conormal}), we proved the Universal Obstruction Theorem~\ref{thm:universal_order6_obstruction}: no pure-configuration generator $K_4 \in \mathcal{P}_0$ can cancel the even-parity grade-5 error $H_5$ (in particular $D_A^4 W_{\mathrm{apex}}^{\mathrm{screen}}$), and no single-commutator corrector can eliminate the quadratic channel $\mathcal{L}_{5,2}$ globally. Lifting the generator to the cotangent bundle $K_4 = c_4^{\mathrm{BM4}} h^4 D_A^3 B \in \mathcal{P}_1$ (realized via $\mathcal{Z}_4^{\mathrm{BM4}}(h)$) to eliminate the linear tidal error $\mathcal{O}(\varepsilon h^4) \in \mathcal{L}_{5,1}$ yielded the processed scheme $\mathrm{PBM}_4$, achieving pre-asymptotic effective order~6 on $\mathcal{D}_{\mathrm{reg}}$ ($\varepsilon \ll h^2$) above the quadratic apex floor $\mathcal{O}(\varepsilon^2 h^4) \in \mathcal{L}_{5,2}$, equipped with the covariant atlas gauge transition operator $T_{A \to B}^{(\mathrm{PBM4})}$ (Theorem~\ref{thm:pbm4_gauge_operator_exact}).
    \item \textbf{Dyadic Multirate Synchronization and Leaf/Apex Decoupling (\S\ref{sec:dyadic_multirate_scheduling}):} Embedding these propagators into a Skeel--Gear dyadic time lattice $h_k = \Delta t_{\mathrm{macro}} / 2^k$ with intra-block chart freezing and inner/outer Leaf/Apex decoupling (Theorem~\ref{thm:leaf_apex_decoupling}) eliminated the $2^{k_{\max}}$ evaluation bottleneck for bound hierarchical binaries while preserving time-reversible shadow energy stability.
\end{enumerate}

With the analytical foundations, obstruction theorems, gauge transitions, and Lean~4 kernel certificates complete, the treatise shifts in Section~\ref{ch:ch06_algorithmic_architecture} from formal geometric numerical analysis to concrete computational realization. Operating on a hierarchical multirate operator-splitting schedule, the algorithmic architecture of Section~\ref{ch:ch06_algorithmic_architecture} coordinates local Keplerian Stumpff--Sundman leaf solvers, composite center-of-mass drifts, screened apex tidal kicks, and covariant Jacobi atlas transitions in vectorized dual precision, delivering deterministic, reproducible, and verifiable high-precision gravitational simulations.