% --- Section Content : Unsoftened Collision Benchmarks ($\varepsilon = 0$) Against State of the Art --- \subsection{Unsoftened Collision Benchmarks ($\varepsilon = 0$) Against State of the Art} \label{sec:unsoftened_benchmarks} We evaluate our unsoftened normal-cone symplectic integrator (built on the cotangent lift of the Fulton--MacPherson deformation to the normal cone) against standard Cartesian integrators and classical symplectic splitting methods (Cartesian Verlet, Strang, Yoshida~4, and direct Blanes--Moan compositions) without any force cutoff ($d_0 = 0$). To establish the empirical validity, convergence properties, and structural boundaries of our unsoftened normal-cone symplectic integrator, we subject the discrete evolution operators constructed in Section~\ref{ch:ch05_higher_order_symplectic} and implemented in Section~\ref{ch:ch06_algorithmic_architecture} to a rigorous battery of numerical experiments. Following the methodological principles articulated by the Paris Observatory celestial mechanics school \cite{laskar2001large,blanes2013optimized,farres2013new_families}, all convergence benchmarks are formulated under the invariant computational cost metric $\tau/s = T_{\mathrm{tot}}/N_{\mathrm{eval}}$, ensuring an unbiassed, Pareto-optimal comparison across integrators possessing differing stage complexities. The benchmark suite is executed simultaneously across two twin arithmetic engines: \begin{enumerate}[label=(\roman*)] \item Standard IEEE~754 double precision (\texttt{Float64}, 53-bit significand, machine precision $\varepsilon_{\mathrm{mach}} \approx 2.22 \times 10^{-16}$), capturing the operational domain of real-time astrodynamical simulations and the onset of Brouwer roundoff diffusion \cite{brouwer1937measures}; \item High-precision software emulation (\texttt{Float128}, 35 decimal digits / 116-bit significand via \texttt{mpmath}, $\varepsilon_{\mathrm{mach}} \approx 1.00 \times 10^{-35}$), providing an unpolluted, machine-floor-free ground truth to confirm formal symplectic convergence orders down to $10^{-30}\text{--}10^{-31}$, central-collision and pericenter invariants down to $10^{-33}\text{--}10^{-38}$, and exact Lie--Bargmann Noether conservation down to $10^{-34}\text{--}10^{-40}$ without numerical saturation. \end{enumerate} % ----------------------------------------------------------------------------- % SUBSECTION 7.1.1 : DUAL-PRECISION CONVERGENCE & INVERTED FORCE COST % ----------------------------------------------------------------------------- \subsubsection{Dual-Precision Symplectic Convergence and the Farr\'es--Laskar Inverted-Cost Metric} \label{sub:laskar_dual_convergence} In standard numerical studies, error is commonly displayed against the macro time-step $h$. However, as demonstrated by Farr\'es, Laskar, Blanes, Casas, Makazaga, and Murua \cite{farres2013new_families, blanes2013optimized}, comparing splitting schemes of differing stages $s$ at identical $h$ obscures their true computational overhead. In chained first-same-as-last (FSAL) execution, the number of force evaluations over a total physical duration $T_{\mathrm{tot}}$ satisfies $N_{\mathrm{eval}} = s \cdot (T_{\mathrm{tot}} / h)$. We therefore adopt the inverted computational cost per evaluation: \begin{equation} \label{eq:laskar_inverted_cost} \frac{\tau}{s} = \frac{T_{\mathrm{tot}}}{N_{\mathrm{eval}}} = \frac{h}{s}, \end{equation} where $s = 1$ for the 2nd-order Strang propagator $\mathcal{S}_2$ (Section~\ref{sec:strang_propagator}), $s = 3$ for the 4th-order Yoshida composition $\mathcal{Y}_4$ (Section~\ref{sec:suzuki_barrier_mcgehee}), $s = 6$ for the 7-stage optimal 4th-order scheme $\mathrm{BM}_4^{(7)}$ and its processed counterpart $\mathrm{PBM}_4$ (Section~\ref{sec:bm4_seven_stages}), and $s = 10$ for the 11-stage optimal 6th-order scheme $\mathrm{BM}_6^{(11)}$ (Section~\ref{sec:bm6_eleven_stages}); here $s = 6$ and $s = 10$ count the effective non-trivial inter-cluster force evaluations per macro-step under the direct ABA/RKN stage-merging formulation of Blanes and Moan \cite{blanes2002practical} (whereas a generic composition of $m = 7$ or $m = 11$ Strang sub-steps $\mathcal{S}_2(w_j\tau)$ with boundary FSAL merging evaluates $m$ inter-cluster forces per macro-step, as detailed in Definition~\ref{def:bm4_seven_stages}). To disentangle asymptotic truncation order from floating-point roundoff and reference-trajectory saturation, our protocol combines two complementary windows on the Figure-8 three-body choreography of Chenciner and Montgomery \cite{chenciner2000remarkable} under the canonical zero-hysteresis Voronoi atlas ($\gamma = 1.0, \gamma_{\mathrm{build}} = 1.0$, with radial-velocity tie-breaking $\dot{r}_{12}(0) = -1.0414 < 0$ on the exact initial equidistance boundary $r_{02}(0) = r_{12}(0)$): \begin{enumerate}[label=(\roman*)] \item a high-resolution dual-precision asymptotic sweep over $T_{\mathrm{conv}} = 0.16 < T_{\mathrm{orb}}/6 \approx 1.054$ across 10 dyadic octaves $h \in [0.16, 3.125 \times 10^{-4}]$ (extended to 13 dyadic octaves down to $h = 3.90625 \times 10^{-5}$ for $\mathrm{BM}_6^{(11)}$ in \texttt{Float128}) against an unoptimized 35-digit \texttt{Float128} $\mathrm{BM}_6^{(11)}$ Richardson-extrapolated $\mathcal{O}(h^8)$ reference trajectory ($h_{\mathrm{ref}} \in \{7.8125 \times 10^{-5}, 3.90625 \times 10^{-5}\}$), and \item a full-period equal-cost benchmark over $T = 6.0$ across all 6 democratic Voronoi chart transitions per period at the common inverted cost $h_{\mathrm{eff}} = h / s = 0.005$ ($N_{\mathrm{eval}} = 1{,}200$ force evaluations, using $h = 0.005$ for $\mathcal{S}_2$, $h = 0.015$ for $\mathcal{Y}_4$, $h = 0.030$ for $\mathrm{BM}_4^{(7)}$ and $\mathrm{PBM}_4$, and $h = 0.050$ for $\mathrm{BM}_6^{(11)}$), benchmarked against an unoptimized 35-digit \texttt{Float128} $\mathrm{BM}_6^{(11)}$ reference integration at $h_{\mathrm{ref}} = 0.01$ ($N_{\mathrm{eval}} = 6{,}000$). \end{enumerate} Crucially, all \texttt{Float64} evaluations execute through the accelerated micro-kernels of Section~\ref{sec:computational_optimizations}, so that comparison against the unoptimized \texttt{Float128} reference simultaneously certifies both symplectic convergence and exact numerical consistency between the optimized and reference solver paths. Furthermore, throughout both precision backends, $\mathrm{BM}_4^{(7)}$ and $\mathrm{PBM}_4$ employ the unified 7-stage symmetric Blanes--Moan kernel of Definition~\ref{def:bm4_seven_stages} coupled to the matching homological corrector $\mathcal{Z}_4^{\mathrm{BM4}}(h) \in \mathcal{P}_1$ with residual moment $c_4^{\mathrm{BM4}} = +1.760005010041013 \times 10^{-4}$. Figure~\ref{fig:laskar_dual_convergence} and Table~\ref{tab:convergence_orders_table} establish four quantitative results: \begin{enumerate} \item \textbf{Exact Asymptotic Orders in \texttt{Float128} down to $10^{-31}$} : In 35-digit quadruple precision, linear regressions on $\log_{10}(\max_t |\Delta H(t)|/|H_0|)$ versus $\log_{10} h_{\mathrm{eff}}$ (Figure~\ref{fig:laskar_dual_convergence}a,b) recover the exact theoretical orders without any high-precision floor saturation: \begin{equation} \begin{aligned} p_{\mathrm{Strang}} &= 2.000, \qquad p_{\mathrm{Yoshida4}} = 3.999, \qquad p_{\mathrm{BM4}} = 4.000, \\ p_{\mathrm{PBM4}} &= 4.000, \qquad p_{\mathrm{BM6}} = 5.999. \end{aligned} \end{equation} Both the Hamiltonian error $\max_t |\Delta H(t)|/|H_0|$ (Figure~\ref{fig:laskar_dual_convergence}a) and the Euclidean trajectory error $\|\mathbf{q}(T) - \mathbf{q}_{\mathrm{ref}}^{\mathrm{f128}}(T)\|_2$ (Figure~\ref{fig:laskar_dual_convergence}b) descend monotonically across thirty decades, reaching $3.74 \times 10^{-31}$ in energy and $1.89 \times 10^{-31}$ in position for $\mathrm{BM}_6^{(11)}$. \item \textbf{Equal-Cost Efficiency ($h_{\mathrm{eff}} = 0.005$) versus Equal-Macro-Step Gain ($h = 0.01$)} : As emphasized by Farr\'es et al.\ \cite{farres2013new_families, blanes2013optimized}, comparing integrators at an identical macro-step $h = 0.01$ yields a $180.9\times$ energy-error reduction ($2.87 \times 10^{-11}$ versus $5.19 \times 10^{-9}$) for $\mathrm{BM}_4^{(7)}$ over $\mathcal{Y}_4$, but requires twice as many force evaluations ($s = 6$ versus $s = 3$). When normalized to the strictly identical evaluation budget $h_{\mathrm{eff}} = h / s = 0.005$ ($N_{\mathrm{eval}} = 1{,}200$ over $T = 6.0$ across 6 Voronoi Jacobi chart transitions, Table~\ref{tab:convergence_orders_table}), $\mathrm{BM}_4^{(7)}$ ($h = 0.030$) outperforms $\mathcal{Y}_4$ ($h = 0.015$) by a net factor of $11.25\times$ in energy ($3.33 \times 10^{-9}$ versus $3.75 \times 10^{-8}$) and $14.5\times$ in trajectory accuracy ($6.77 \times 10^{-8}$ versus $9.78 \times 10^{-7}$), while $\mathrm{PBM}_4$ ($h = 0.030$) achieves a $6.85\times$ energy gain ($5.47 \times 10^{-9}$) over $\mathcal{Y}_4$ across the 6 non-hierarchical chart transitions ($\varepsilon \sim 1$), and $\mathrm{BM}_6^{(11)}$ ($h = 0.050$) achieves a $6{,}354\times$ gain in energy ($5.90 \times 10^{-12}$) and a $1.60 \times 10^5\times$ gain in position ($6.11 \times 10^{-12}$) over $\mathcal{Y}_4$. \item \textbf{Non-Hierarchical Choreography ($\varepsilon \sim \mathcal{O}(1)$) versus 3D Hierarchical Perturbations ($\varepsilon \ll 1$)} : On the equal-mass Figure-8 choreography (Figure~\ref{fig:laskar_dual_convergence}a,b), the three bodies exchange roles symmetrically with no small perturbation parameter ($\varepsilon \sim \mathcal{O}(1)$). Consequently, the non-perturbative Lie-bracket terms $\mathcal{O}(\varepsilon^2 h^4)$ in the shadow Hamiltonian have unit weight, so $\mathrm{PBM}_4$ rightly exhibits asymptotic order $p_{\mathrm{eff}} = 4.000$ in the non-hierarchical regime. Conversely, when the exact same five canonical integrators are evaluated on a 3D eccentric hierarchical triple ($m_1 = m_2 = 0.5, a_{\mathrm{in}} = 1.0, e_{\mathrm{in}} = 0.30$ perturbed by an inclined outer companion $m_3 = 1.0, a_{\mathrm{out}} = 12.0, i = 28^\circ$, where $\eta_S \approx a_{\mathrm{in}} / a_{\mathrm{out}} = 1/12 \ll 1$), Panels~\ref{fig:laskar_dual_convergence}(c,d) demonstrate that $\mathrm{PBM}_4$ ($s = 6$) eliminates the leading $\mathcal{O}(\varepsilon h^4)$ error term and detaches from $\mathrm{BM}_4^{(7)}$ along an effective 6th-order slope ($p_{\mathrm{eff}} \approx 6.25$) parallel to $\mathrm{BM}_6^{(11)}$ ($s = 10$), achieving up to a $2{,}329\times$ energy-error reduction ($1.53 \times 10^{-14}$ versus $3.55 \times 10^{-11}$ in \texttt{Float128}, $1{,}762\times$ in \texttt{Float64}, Panel~c) and a $1{,}640\times$ trajectory-position error reduction ($8.29 \times 10^{-14}$ versus $1.36 \times 10^{-10}$, Panel~d) over bare $\mathrm{BM}_4^{(7)}$ at $h_{\mathrm{eff}} = 8.18 \times 10^{-3}$. \item \textbf{Brouwer Roundoff Bounce in \texttt{Float64}} : In IEEE~754 double precision, the deterministic $\mathcal{O}(h^p)$ symplectic truncation error dominates for $h_{\mathrm{eff}} \gtrsim 10^{-3}$ and reaches the machine floor $\sim 3.45 \times 10^{-16}$ at $h_{\mathrm{eff}} \approx 10^{-3}$; for finer steps $h_{\mathrm{eff}} < 10^{-3}$ ($N_{\mathrm{steps}} = T/h \gg 10^3$), truncation falls beneath IEEE~754 precision and the measured error rebounds as $h_{\mathrm{eff}}^{-1/2} \propto \sqrt{N_{\mathrm{steps}}}$. In accordance with Brouwer's law \cite{brouwer1937measures}, accumulated floating-point roundoff across the $s_{\mathrm{eval}}\,N_{\mathrm{steps}}$ force evaluations executes an unbiased random walk $\sigma \propto \varepsilon_{\mathrm{mach}} \sqrt{s_{\mathrm{eval}}\,N_{\mathrm{steps}}}$, confirming that 35-digit \texttt{Float128} arithmetic is essential to resolve the unpolluted asymptotic regime of 6th-order symplectic compositions. \end{enumerate} \pythonrefright{py_bench_gen01}{{benchmarks/run, figures/plot}_generator_01_*.py}{Float64 / Float128} \begin{figure}[!htbp] \centering \includegraphics[width=\textwidth]{chapters/ch07_benchmarks_and_lean4/figures/fig_ch07_laskar_dual_precision_convergence.png} \caption{\textbf{Dual-precision symplectic convergence and Farr\'es--Laskar inverted-cost frontiers ($\mathcal{S}_2, \mathcal{Y}_4, \mathrm{BM}_4^{(7)}, \mathrm{PBM}_4, \mathrm{BM}_6^{(11)}$; $\varepsilon = 0, d_0 = 0, \gamma = 1.0$).} \textbf{(a)}~Figure-8 energy error $\max_t |\Delta H(t)| / |H_0|$ vs.\ $h_{\mathrm{eff}} = h / s_{\mathrm{eval}}$ (\texttt{Float64} solid, \texttt{Float128} dashed). \textbf{(b)}~Figure-8 trajectory error $\|\mathbf{q}(T) - \mathbf{q}_{\mathrm{ref}}^{\mathrm{f128}}(T)\|_2$ vs.\ $h_{\mathrm{eff}}$. \textbf{(c)}~3D hierarchical triple ($\varepsilon \ll 1, e_{\mathrm{in}} = 0.30, i = 28^\circ, a_{\mathrm{out}} = 12.0$): energy error showing $\mathrm{PBM}_4$ order-$4 \to 6$ elevation ($2{,}329\times$ gain over $\mathrm{BM}_4^{(7)}$). \textbf{(d)}~3D hierarchical triple trajectory error ($1{,}640\times$ gain for $\mathrm{PBM}_4$).} \label{fig:laskar_dual_convergence} \end{figure} \begin{table}[!htbp] \centering \caption{\textbf{Effective convergence orders, equal-cost accuracy ($h_{\mathrm{eff}} = h/s = 0.005, N_{\mathrm{eval}} = 1{,}200$ over $T = 6.0$, Voronoi $\gamma = 1.0$), and dual-precision error floors on the Figure-8 choreography (\S\ref{sub:laskar_dual_convergence}).} Here $s$ ($m$) denotes effective FSAL force evaluations per macro-step under direct ABA/RKN merging (generic $\mathcal{S}_2$ composition), and $\|\Delta\mathbf{q}(T)\|_2 = \|\mathbf{q}(T) - \mathbf{q}_{\mathrm{ref}}^{\mathrm{f128}}(T)\|_2$.} \label{tab:convergence_orders_table} \vspace{0.3em} \renewcommand{\arraystretch}{1.2} \resizebox{\textwidth}{!}{% \begin{tabular}{lccccccc} \toprule \textbf{Integrator} & \makecell{\textbf{Stages}\\$s$ ($m$)} & \makecell{\textbf{Order}\\\textbf{(Th.)}} & \makecell{$p_{\mathrm{eff}}$\\\textbf{(\texttt{Float128})}} & \makecell{$\max_t |\Delta H|/|H_0|$\\\textbf{at } $h_{\mathrm{eff}} = 0.005$} & \makecell{$\|\Delta\mathbf{q}(T)\|_2$\\\textbf{at } $h_{\mathrm{eff}} = 0.005$} & \makecell{\textbf{Min Floor}\\\textbf{(\texttt{Float64})}} & \makecell{\textbf{Min Floor}\\\textbf{(\texttt{Float128})}} \\ \midrule $\mathcal{S}_2$ (Strang) & 1 (1) & 2 & $2.000$ & $2.81 \times 10^{-6}$ & $5.26 \times 10^{-5}$ & $5.61 \times 10^{-9}$ & $5.61 \times 10^{-9}$ \\ $\mathcal{Y}_4$ (Yoshida) & 3 (3) & 4 & $3.999$ & $3.75 \times 10^{-8}$ & $9.78 \times 10^{-7}$ & $1.38 \times 10^{-14}$ & $4.95 \times 10^{-15}$ \\ $\mathrm{BM}_4^{(7)}$ & 6 (7) & 4 & $4.000$ & $3.33 \times 10^{-9}$ ($11.25\times$) & $6.77 \times 10^{-8}$ ($14.5\times$) & $6.90 \times 10^{-16}$ & $2.74 \times 10^{-17}$ \\ $\mathrm{PBM}_4$ & 6 (7) & 4 (eff.\ 6) & $4.000$ ($6.25$ in 3D) & $5.47 \times 10^{-9}$ ($6.85\times$) & $2.09 \times 10^{-7}$ ($4.67\times$) & $8.63 \times 10^{-16}$ & $3.51 \times 10^{-17}$ \\ $\mathrm{BM}_6^{(11)}$ & 10 (11) & 6 & $5.999$ & $5.90 \times 10^{-12}$ ($6{,}354\times$) & $6.11 \times 10^{-12}$ ($1.60 \times 10^5\times$) & $3.45 \times 10^{-16}$ & $3.74 \times 10^{-31}$ \\ \midrule \textbf{Ref.\ $\mathrm{BM}_6^{(11)}$ (\texttt{Float128}, unopt.)} & 10 (11) & 6 & $5.999$ & $6.60 \times 10^{-17}$ ($h_{\mathrm{eff}}=10^{-3}$) & $0$ (Reference) & --- & $1.89 \times 10^{-31}$ ($\|\Delta\mathbf{q}\|_2$) \\ \bottomrule \end{tabular}% } \end{table} % ----------------------------------------------------------------------------- % SUBSECTION 7.1.2 : TWO-CHANNEL ERROR LAW, ORDER 6->4 CROSSOVER & QUADRATIC FLOOR % ----------------------------------------------------------------------------- \subsubsection{Experimental Proof of the Two-Channel Error Law, Order $6 \to 4$ Crossover, and the Hall-Cohomology Quadratic Floor (Section~\ref{sub:universal_order6_obstruction} and Section~\ref{sec:symplectic_correctors})} \label{sub:pbm4_quadratic_floor_proof} Section~\ref{sub:laskar_dual_convergence} revealed a fundamental dichotomy in the behavior of the processed integrator $\mathrm{PBM}_4$: on the non-hierarchical Figure-8 choreography ($\varepsilon \sim \mathcal{O}(1)$, Figure~\ref{fig:laskar_dual_convergence}a,b), $\mathrm{PBM}_4$ scales with asymptotic order $p_{\mathrm{eff}} = 3.994 \approx 4$, whereas on the 3D eccentric hierarchical triple ($\varepsilon \ll 1$, Figure~\ref{fig:laskar_dual_convergence}c,d), it detaches from $\mathrm{BM}_4^{(7)}$ along an effective 6th-order trajectory ($p_{\mathrm{eff}} = 6.25$). The algebraic theory developed in Section~\ref{sec:symplectic_correctors} (Theorems~\ref{thm:universal_order6_obstruction}, \ref{thm:c4_corrector_solvability}, and \ref{thm:asymptotic_conormal_scaling}) explains and predicts this transition across the entire $(\varepsilon, h)$ plane. Because the 2-stage palindromic Lie conjugator $\mathcal{Z}_4^{\mathrm{BM4}}(h)$ generated by the cotangent-lifted corrector $K_4 = c_4^{\mathrm{BM4}} h^4 D_A^3 B \in \mathcal{P}_1$ identically annihilates the leading grade-$(5,1)$ linear tidal commutator $h^4 H_5^{\mathrm{linear}} = c_4^{\mathrm{BM4}} h^4 D_A^4 B = \mathcal{O}(\varepsilon h^4)$, both the relative Hamiltonian error and the phase-space trajectory error obey the two-channel composite error law: \begin{equation} \label{eq:composite_error_law_exp} \frac{|\Delta H(T)|}{|H_0|} \sim C_6\,\varepsilon\,h^6 + C_4\,\varepsilon^2\,h^4, \qquad \|\mathbf{q}(T) - \mathbf{q}_{\mathrm{ref}}(T)\|_2 \sim \tilde{C}_6\,\varepsilon\,h^6 + \tilde{C}_4\,\varepsilon^2\,h^4, \end{equation} where $\varepsilon = (m_{\mathrm{pert}}/M_{\mathrm{bin}})(a_{\mathrm{in}}/R_{\mathrm{out}})^3 = \eta_S^3$ parametrizes the dimensionless cubic tidal perturbation strength, the first term $\mathcal{O}(\varepsilon h^6)$ stems from the uncompensated grade-$(7,1)$ seventh-order commutators, and the irreducible quadratic floor $\mathcal{O}(\varepsilon^2 h^4)$ originates jointly from the non-vanishing grade-$(5,2)$ Hall cohomology class $[D_B D_A^3 B] \in H^1(D_A, \mathcal{L}_{5,2}) \cong \mathbb{R}$ (Theorem~\ref{thm:universal_order6_obstruction}, Tier~2) and the cross-term $\{K_4, W_{\mathrm{apex}}^{\mathrm{screen}}\} \in \mathcal{P}_0$. Figure~\ref{fig:pbm4_quadratic_floor} subjects Equation~\eqref{eq:composite_error_law_exp} to a four-part quantitative verification on the 3D eccentric hierarchical triple ($m_0 = m_1 = 0.5, a_{\mathrm{in}} = 1.0, e_{\mathrm{in}} = 0.30$ perturbed by an outer companion $m_2 = 1.0$ on an inclined orbit $i = 28^\circ$), comparing exclusively $\mathrm{PBM}_4$ against its unprocessed base scheme $\mathrm{BM}_4^{(7)}$ at identical evaluation cost ($s_{\mathrm{eval}} = 6$): \begin{enumerate}[label=(\alph*)] \item \textbf{Order $6 \to 4$ Energy Error Crossover and $100\times$/Decade Quadratic Floor in \texttt{Float128} (Figure~\ref{fig:pbm4_quadratic_floor}a)} : Sweeping macro-steps $h \in [0.01, 0.32]$ ($h_{\mathrm{eff}} = h/6 \in [1.67 \times 10^{-3}, 5.33 \times 10^{-2}]$) over $T_{\mathrm{cross}} = 1.60$ in 35-digit \texttt{Float128} arithmetic across three tidal hierarchy decades $\varepsilon \in \{10^{-2}, 10^{-3}, 10^{-4}\}$ ($R_{\mathrm{out}} \in \{2.92, 6.30, 13.57\}$), unprocessed $\mathrm{BM}_4^{(7)}$ exhibits a uniform 4th-order decay $\mathcal{O}(\varepsilon h^4)$ whose curves are separated vertically by $10\times$ per decade in $\varepsilon$. In contrast, above the crossover scale $h_{\mathrm{cross}}(\varepsilon) = \sqrt{(C_4/C_6)\,\varepsilon}$, $\mathrm{PBM}_4$ follows the 6th-order channel $\mathcal{O}(\varepsilon h^6)$, achieving a maximum energy accuracy gain over $\mathrm{BM}_4^{(7)}$ of $19.1\times$ at $\varepsilon = 10^{-2}$, \textbf{$262.5\times$} at $\varepsilon = 10^{-3}$ ($5.94 \times 10^{-12}$ versus $1.56 \times 10^{-9}$ at $h = 0.08$), and \textbf{$2{,}537\times$} at $\varepsilon = 10^{-4}$ ($4.62 \times 10^{-15}$ versus $1.17 \times 10^{-11}$ at $h = 0.04$). Below $h_{\mathrm{cross}}(\varepsilon)$, the grade-$(5,2)$ Hall-cohomology floor $\mathcal{O}(\varepsilon^2 h^4)$ takes over and bends the $\mathrm{PBM}_4$ curves back onto slope $4.00$. At the finest step $h = 0.01$ ($h_{\mathrm{eff}} = 1.67 \times 10^{-3}$), the $\mathrm{PBM}_4$ energy error descends from $1.49 \times 10^{-13}$ ($\varepsilon = 10^{-2}$) to $2.30 \times 10^{-15}$ ($\varepsilon = 10^{-3}$) and $2.52 \times 10^{-17}$ ($\varepsilon = 10^{-4}$, a $1{,}815\times$ floor reduction against $4.58 \times 10^{-14}$ for $\mathrm{BM}_4^{(7)}$), demonstrating the exact $100\times$ reduction per decade in $\varepsilon$ ($\propto \varepsilon^2$) predicted by Theorem~\ref{thm:universal_order6_obstruction}. \item \textbf{Unclipped Local Empirical Order $p_{\mathrm{emp}}(h) : 6 \to 4$ and Square-Root Crossover Shift (Figure~\ref{fig:pbm4_quadratic_floor}b)} : Differentiating Equation~\eqref{eq:composite_error_law_exp} with respect to $\ln h$ yields the exact analytical local-order transition law: \begin{equation} \label{eq:local_order_transition_law} p_{\mathrm{th}}(h;\varepsilon) = \frac{\partial \ln\bigl(C_6\,\varepsilon\,h^6 + C_4\,\varepsilon^2\,h^4\bigr)}{\partial \ln h} = 4 + \frac{2}{1 + \bigl(h_{\mathrm{cross}}(\varepsilon)/h\bigr)^2}, \qquad h_{\mathrm{cross}}(\varepsilon) = \sqrt{\frac{C_4}{C_6}\,\varepsilon}. \end{equation} Evaluating the unclipped finite-difference trajectory slope directly from the raw \texttt{Float128} data: \[ p_{\mathrm{emp}}\bigl(\sqrt{h_k h_{k+1}}\bigr) = \frac{\ln\bigl(E_q(h_k)/E_q(h_{k+1})\bigr)}{\ln(h_k/h_{k+1})}, \qquad E_q(h) = \bigl\|\mathbf{q}(T; h) - \mathbf{q}_{\mathrm{ref}}^{\mathrm{f128}}(T)\bigr\|_2 \] (where $E_q(h)$ denotes the Euclidean phase-space position error of Panel~(d)) closely tracks the universal sigmoid~\eqref{eq:local_order_transition_law}. As $\varepsilon$ decreases by a factor of $100$ from $10^{-2}$ to $10^{-4}$, the midpoint inflection $p_{\mathrm{th}}(h_{\mathrm{cross}}) = 5.0$ shifts leftward by a factor of $\sqrt{100} = 10$ in step size: for $\varepsilon = 10^{-2}$, $p_{\mathrm{emp}}$ relaxes rapidly from $5.17$ down to $4.28 \to 4.05 \to 4.01 \to 4.00$; for $\varepsilon = 10^{-3}$, it transitions smoothly along $6.12 \to 5.60 \to 4.76 \to 4.20 \to 4.05$; and for $\varepsilon = 10^{-4}$, it remains anchored in the 6th-order regime ($6.23 \to 6.03 \to 5.84 \to 5.32$) across nearly the entire sweep before descending to $4.53$ at the smallest step pair, while $\mathrm{BM}_4^{(7)}$ stays flat at $p_{\mathrm{emp}} = 3.99 \pm 0.02$. \item \textbf{Two-Regime Tidal Perturbation Scaling versus $\varepsilon = \eta_S^3$ at Fixed $h$ (Figure~\ref{fig:pbm4_quadratic_floor}c)} : Holding the macro-step fixed at $h = 0.08$ ($h_{\mathrm{eff}} = 0.0133$, executed in \texttt{Float64} with the accelerated micro-kernels of Section~\ref{sec:computational_optimizations}) and sweeping the tidal parameter across 3.5 decades $\varepsilon \in [10^{-5}, 3 \times 10^{-2}]$ exposes the dual role of $\varepsilon$ in Equation~\eqref{eq:composite_error_law_exp}, separated by the crossover threshold $\varepsilon_{\mathrm{cross}}(h) = (C_6/C_4)\,h^2 \approx 6.4 \times 10^{-3}$. Whereas $\mathrm{BM}_4^{(7)}$ scales linearly as $\mathcal{O}(\varepsilon h^4) \propto \varepsilon^1$ across the entire domain ($2.09 \times 10^{-8}$ at $\varepsilon = 3 \times 10^{-2}$ down to $2.23 \times 10^{-11}$ at $\varepsilon = 10^{-5}$), $\mathrm{PBM}_4$ exhibits two distinct asymptotic regimes: (i)~in the moderate-hierarchy regime $\varepsilon \in [10^{-3}, 3 \times 10^{-2}]$ ($\varepsilon \gtrsim \varepsilon_{\mathrm{cross}}$), the second channel dominates and $\mathrm{PBM}_4$ follows the quadratic floor scaling $\mathcal{O}(\varepsilon^2 h^4) \propto \varepsilon^2$ (empirical slope $1.84$), expanding its accuracy advantage over $\mathrm{BM}_4^{(7)}$ from $5.3\times$ at $\varepsilon = 3 \times 10^{-2}$ to \textbf{$236.1\times$} at $\varepsilon = 10^{-3}$ ($7.32 \times 10^{-12}$ versus $1.73 \times 10^{-9}$); (ii)~in the deep-hierarchy regime $\varepsilon \in [10^{-5}, 3 \times 10^{-4}]$ ($\varepsilon \ll \varepsilon_{\mathrm{cross}}$), the quadratic floor drops below the 6th-order truncation term and $\mathrm{PBM}_4$ transitions onto the pre-asymptotic channel $\mathcal{O}(\varepsilon h^6) \propto \varepsilon^1$ (empirical slope $0.93$), maintaining a constant $\mathcal{O}(h^{-2})$ accuracy plateau of $215\times$--$236\times$ over $\mathrm{BM}_4^{(7)}$ down to $\varepsilon = 10^{-5}$ ($1.04 \times 10^{-13}$ versus $2.23 \times 10^{-11}$). \item \textbf{Phase-Space Trajectory Position Error Crossover $\|\mathbf{q}(T) - \mathbf{q}_{\mathrm{ref}}^{\mathrm{f128}}(T)\|_2$ (Figure~\ref{fig:pbm4_quadratic_floor}d)} : Benchmarking the Euclidean phase-space configuration error at $T_{\mathrm{cross}} = 1.60$ against an independent, unoptimized 35-digit \texttt{Float128} $\mathrm{BM}_6^{(11)}$ reference trajectory (Richardson-extrapolated at $h_{\mathrm{ref}} \in \{0.01, 0.005\}$) confirms that the boundary readout map $\mathbf{z}^{\mathrm{phys}} = \mathcal{Z}_4^{\mathrm{BM4}}(h)(\tilde{\mathbf{z}})$ transfers the exact same two-channel law from the shadow Hamiltonian to physical orbital trajectories. Above $h_{\mathrm{cross}}(\varepsilon)$, $\mathrm{PBM}_4$ follows the 6th-order trajectory slope and delivers up to a \textbf{$1{,}335\times$} position-accuracy gain over $\mathrm{BM}_4^{(7)}$ at $\varepsilon = 10^{-4}$ ($2.49 \times 10^{-17}$ versus $3.33 \times 10^{-14}$ at $h = 0.01$, and $134.1\times$ at $\varepsilon = 10^{-3}$). Below $h_{\mathrm{cross}}(\varepsilon)$, the trajectory error bends onto the 4th-order floor $\mathcal{O}(\varepsilon^2 h^4)$ with the exact same $100\times$-per-decade vertical separation in $\varepsilon$: at $h = 0.01$, $\|\Delta\mathbf{q}(T)\|_2$ drops from $2.49 \times 10^{-13}$ ($\varepsilon = 10^{-2}$) to $2.30 \times 10^{-15}$ ($\varepsilon = 10^{-3}$) and $2.49 \times 10^{-17}$ ($\varepsilon = 10^{-4}$). \end{enumerate} \pythonrefright{py_bench_gen02}{{benchmarks/run, figures/plot}_generator_02_*.py}{Float64 / Float128} \begin{figure}[!htbp] \centering \includegraphics[width=\textwidth]{chapters/ch07_benchmarks_and_lean4/figures/fig_ch07_pbm4_quadratic_floor_crossover.png} \caption{\textbf{Two-channel error law, order $6 \to 4$ crossover, and Hall-cohomology quadratic floor $\mathcal{O}(\varepsilon^2 h^4)$ for $\mathrm{PBM}_4$ vs.\ $\mathrm{BM}_4^{(7)}$ on the 3D eccentric hierarchical triple ($e_{\mathrm{in}} = 0.30, i = 28^\circ$).} \textbf{(a)}~\texttt{Float128} energy error $|\Delta H(T)|/|H_0|$ vs.\ $h_{\mathrm{eff}} = h/6$ across $\varepsilon \in \{10^{-2}, 10^{-3}, 10^{-4}\}$. \textbf{(b)}~Unclipped local empirical order $p_{\mathrm{emp}}(h)$ vs.\ $\sqrt{h_k h_{k+1}}$ overlaid with $p_{\mathrm{th}}(h;\varepsilon)$~\eqref{eq:local_order_transition_law}. \textbf{(c)}~Tidal scaling $\max_t |\Delta H(t)|/|H_0|$ vs.\ $\varepsilon = \eta_S^3$ at $h = 0.08$ (\texttt{Float64}). \textbf{(d)}~\texttt{Float128} trajectory error $\|\mathbf{q}(T) - \mathbf{q}_{\mathrm{ref}}^{\mathrm{f128}}(T)\|_2$ vs.\ $h_{\mathrm{eff}}$.} \label{fig:pbm4_quadratic_floor} \end{figure} % ----------------------------------------------------------------------------- % SUBSECTION 7.1.3 : UNCONSTRAINED SINGULARITIES & LIMIT CASES % ----------------------------------------------------------------------------- \subsubsection{Dynamics Across Singularities: Unsoftened Limit Cases ($\varepsilon = 0, d_0 = 0$)} \label{sub:unsoftened_limit_cases} A defining requirement of this monograph is the deterministic traversal of physical singularities without the introduction of ad-hoc smoothing lengths ($d_0 = 0$) or artificial force cutoffs. Rather than grouping disparate dynamical systems onto a single composite axis, we structure our evaluation around four dedicated benchmark plates (Figures~\ref{fig:l1_pericenter_conditioning}--\ref{fig:l4_adaptive_multiscale_burrau_kozai}), each isolating one of the four singularity-resolution mechanisms established in Sections~\ref{ch:ch04_fulton_macpherson_splitting}--\ref{ch:ch06_algorithmic_architecture}. \medskip\noindent\textbf{Limit Case L1 : Analytic Kepler--Sundman--Stumpff Pericenter Conditioning Across $1 - e \in [10^{-12}, 10^{-1}]$ and $e = 1.0$ ($N=2$, Figure~\ref{fig:l1_pericenter_conditioning}).} In unregularized Cartesian integrators (such as velocity Verlet), the pericenter acceleration diverges as $r_{\min}^{-2} = [a(1-e)]^{-2}$, causing a catastrophic energy error explosion scaling as $(1-e)^{-5/2}$ ($6.61 \times 10^2$ at $1-e = 10^{-3}$ and $6.55 \times 10^2$ at $1-e = 10^{-12}$, Figure~\ref{fig:l1_pericenter_conditioning}b). Even within a universal Sundman--Stumpff propagator, naive evaluation of $\rho_{\mathrm{raw}} = \rho_0 c_0(z) + \sigma_0 \tau c_1(z) + GM\tau^2 c_2(z)$ at pericenter suffers subtractive cancellation of order $(1-e)^{-2}\varepsilon_{\mathrm{mach}}$ in total energy ($\sim 1.39 \times 10^{-6}$ at $1-e = 10^{-5}$ and $\mathcal{O}(1)$ for $1-e \le 10^{-8}$ in \texttt{Float64}). By enforcing the exact angular-momentum/energy quadratic identity coupled with incremental Wronskian renormalization (Section~\ref{sub:kepler_sundman_analytic_flow}, Proposition~\ref{prop:battin_unimodular_projection}, and Theorem~\ref{thm:leaf_apex_decoupling} in Section~\ref{sub:leaf_apex_decoupling}): \[ \rho_{\mathrm{final}} = \frac{\sigma_{\mathrm{end}}^2 + \|\mathbf{L}\|^2/\mu^2}{2GM - \alpha\rho_{\mathrm{raw}}}, \qquad \delta\dot{g} = \frac{\dot{f}g - \delta f}{f}, \] and integrating over a complete orbital period $1.0\,T_{\mathrm{orb}} = \sqrt{2}\pi \approx 4.4429$ using an odd step count $N_{\mathrm{steps}} = 129$ ($\Delta t = T_{\mathrm{orb}}/129 \approx 0.03444 \le 0.05\,\tau_{\mathrm{dyn}}(a_C)$, $k_{\mathrm{sym}}(t) \equiv 0$) so that central Step~64 straddles $t = 0.5\,T_{\mathrm{orb}}$ inside the regularized Sundman chart down to $r_{\min} = 10^{-12}\text{ AU}$ (Figure~\ref{fig:l1_pericenter_conditioning}a), the conditioned $N=2$ leaf propagator eliminates pericenter loss across all 12 decades, maintaining $\max_{t\in[0,T_{\mathrm{orb}}]} |\Delta H|/|H_0| \le 1.39 \times 10^{-13}$ ($8.17 \times 10^{-14}$ at $1-e = 10^{-12}$ and $9.95 \times 10^{-14}$ at the rectilinear collision limit $e = 1.0$) in \texttt{Float64}, and $\le 1.53 \times 10^{-33}$ ($2.71 \times 10^{-34}$ at $1-e = 10^{-12}$ and $5.50 \times 10^{-38}$ at $e = 1.0$) in \texttt{Float128}. \begin{figure}[!htbp] \centering \includegraphics[width=\textwidth]{chapters/ch07_benchmarks_and_lean4/figures/fig_ch07_l1_pericenter_conditioning.png} \caption{\textbf{Limit Case L1: Analytic Kepler--Sundman--Stumpff pericenter conditioning across $1 - e \in [10^{-12}, 10^{-1}]$ and $e = 1.0$ ($N=2$).} \textbf{(a)}~Radial separation $r(t)$ vs.\ pericenter-centered normalized time $\Delta\tilde{t} = (t - 0.5\,T_{\mathrm{orb}})/T_{\mathrm{orb}}$ (\texttt{symlog} scale) down to $r_{\min} = 10^{-12}\text{ AU}$ within Step~64 ($N_{\mathrm{steps}} = 129, k_{\mathrm{sym}}(t) \equiv 0$). \textbf{(b)}~Full-orbit relative energy error $\max_t |\Delta H(t)|/|H_0|$ across 12 decades in $1-e$ (\texttt{Float64} and \texttt{Float128}).} \label{fig:l1_pericenter_conditioning} \end{figure} \medskip\noindent\textbf{Limit Case L2 : Asymptotic Free-Fall Regularity of the $\mathrm{PBM}_4$ ($\mathcal{Z}_4^{\mathrm{BM4}}$) Cotangent Corrector Across 12 Decades $r \in [10^{-12}, 10^0]\text{ AU}$ (Figure~\ref{fig:l2_pbm4_free_fall_regularity}).} Why does the cotangent boundary corrector $\mathcal{Z}_4^{\mathrm{BM4}}(h_{\mathrm{sub}})$ remain completely nonsingular when a binary pair plunges toward collision ($r \to 0^+$) despite its momentum dependence on the diverging free-fall impulse $\|\mathbf{w}(r)\| = \mu\sqrt{2GM/r - \alpha} \propto r^{-1/2}$? By Theorems~\ref{thm:c4_corrector_solvability} and~\ref{thm:asymptotic_conormal_scaling} (Equations~\eqref{eq:k4_conormal_generator_explicit}--\eqref{eq:msym_dynamic_quadrupole_tensor}) and Theorem~\ref{thm:shadow_alignment} (Equations~\eqref{eq:shadow_entry_operator} and~\eqref{eq:corrector_verlet_step}), the corrector is governed by the scalar generating potential and dynamic quadrupole tensor: \[ \begin{aligned} K_4(\mathbf{r},\mathbf{w}) &= \operatorname{Tr}\bigl(\mathbf{T}_{\mathrm{apex}}(\mathbf{R})\,\mathbf{M}_{\mathrm{sym}}(\mathbf{r},\mathbf{w})\bigr) \in \mathcal{P}_1, \\ \mathbf{M}_{\mathrm{sym}}(\mathbf{r},\mathbf{w}) &= c_4^{\mathrm{BM4}} h_{\mathrm{sub}}^4 \frac{GM_{ab}}{r^3}\biggl[2(\mathbf{r}\otimes\mathbf{w} + \mathbf{w}\otimes\mathbf{r}) - \frac{3(\mathbf{r}\cdot\mathbf{w})}{r^2}\mathbf{r}\otimes\mathbf{r}\biggr], \end{aligned} \] generating the corrective position shift $\Delta\mathbf{r}_{\mathrm{corr}} = +\nabla_{\mathbf{w}}K_4$ and momentum shift $\Delta\mathbf{w}_{\mathrm{corr}} = -\nabla_{\mathbf{r}}K_4$. Because the local dyadic sub-step size scales with the Keplerian dynamical timescale $h_{\mathrm{sub}}(r) = h_0 (r/r_{\mathrm{ref}})^{3/2} \propto \tau_{\mathrm{dyn}}(r) = \sqrt{r^3/(GM_{ab})} \propto r^{3/2}$ ($h_{\mathrm{sub}}^4 \propto r^6$), the prefactor $h_{\mathrm{sub}}(r)^4 \frac{GM_{ab}}{r^3} \propto r^{+3}$ overwhelms both the $r^{-3}$ Keplerian Hessian and the $w(r) \propto r^{-1/2}$ free-fall velocity divergence. Sweeping $r \in [10^{-12}, 10^0]\text{ AU}$ across 12 orders of magnitude ($N=3, R_2 = \sqrt{29}\text{ AU}, h_0 = 0.01$) confirms the exact analytical extinction laws: \begin{equation} \label{eq:pbm4_free_fall_power_laws} |K_4(\mathbf{r},\mathbf{w})| \propto r^{+3.5000} = \mathcal{O}(r^{7/2}) \to 0, \qquad \|\Delta\mathbf{r}_{\mathrm{corr}}\| = \|\nabla_{\mathbf{w}}K_4\| \propto r^{+4.0000} = \mathcal{O}(r^4) \to 0, \end{equation} \begin{equation} \label{eq:pbm4_free_fall_momentum_laws} \|\Delta\mathbf{w}_{\mathrm{corr}}\| = \|\nabla_{\mathbf{r}}K_4\| \propto r^{+2.5000} = \mathcal{O}(r^{5/2}) \to 0, \qquad \frac{\|\Delta\mathbf{w}_{\mathrm{corr}}\|}{\|\mathbf{w}\|} \propto r^{+3.0000} = \mathcal{O}(r^3) \to 0, \end{equation} with the generating potential $|K_4|$ decaying from $1.21 \times 10^{-14}$ down to $1.40 \times 10^{-56}$ and the corrective position shift $\|\Delta\mathbf{r}_{\mathrm{corr}}\|$ decaying from $6.04 \times 10^{-14}\text{ AU}$ down to $6.04 \times 10^{-62}\text{ AU}$ at $r = 10^{-12}\text{ AU}$ (Figure~\ref{fig:l2_pbm4_free_fall_regularity}a), while the absolute momentum shift $\|\Delta\mathbf{w}_{\mathrm{corr}}\|$ decays to $4.06 \times 10^{-44}$ and the relative momentum shift $\|\Delta\mathbf{w}_{\mathrm{corr}}\|/\|\mathbf{w}\|$ decays to $4.06 \times 10^{-50}$ (Figure~\ref{fig:l2_pbm4_free_fall_regularity}b), proving that $\mathcal{Z}_4^{\mathrm{BM4}}(h_{\mathrm{sub}})$ smoothly reduces to the exact identity in both position and momentum as $r \to 0^+$. \pythonrefright{py_bench_gen03}{benchmarks/{run_L1_eccentricity, run_L2_mcgehee}.py}{Float64 / Float128} \begin{figure}[!htbp] \centering \includegraphics[width=\textwidth]{chapters/ch07_benchmarks_and_lean4/figures/fig_ch07_l2_pbm4_free_fall_regularity.png} \caption{\textbf{Limit Case L2: Asymptotic free-fall regularity of the $\mathrm{PBM}_4$ ($\mathcal{Z}_4^{\mathrm{BM4}}$) cotangent corrector across $r \in [10^{-12}, 10^0]\text{ AU}$.} \textbf{(a)}~Power-law extinction of generating potential $|K_4(\mathbf{r},\mathbf{w})| \propto r^{+3.5000}$ and position shift $\|\Delta\mathbf{r}_{\mathrm{corr}}\| \propto r^{+4.0000}$. \textbf{(b)}~Extinction of absolute momentum shift $\|\Delta\mathbf{w}_{\mathrm{corr}}\| \propto r^{+2.5000}$ and relative momentum shift $\|\Delta\mathbf{w}_{\mathrm{corr}}\|/\|\mathbf{w}\| \propto r^{+3.0000}$.} \label{fig:l2_pbm4_free_fall_regularity} \end{figure} \medskip\noindent\textbf{Limit Case L3 : Exact McGehee Exceptional Divisor Traversals ($\rho = 0$), Conormal Bounce Involution $\mathcal{B}_C^2 = \mathrm{Id}$, and Round-Trip Reversibility for Central Collisions $k \in \{2, 3, 4\}$ (Figure~\ref{fig:l3_mcgehee_central_collisions}).} Under the normal-cone radial blow-up $\rho_C = \sqrt{I_C}$ and McGehee regularized radial velocity $v_{\mathrm{McG}} = \dot{\rho}_C\sqrt{\rho_C}$ (Section~\ref{sec:suzuki_barrier_mcgehee}, Lemma~\ref{thm:mcgehee_bounce_unitary}, and Sections~\ref{sec:synchronous_supernodes} and~\ref{sub:synchronous_supernodes_flow}), simultaneous $k$-body central collisions ($d_0 = 0$) are resolved on the blown-up collision manifold by the unitary conormal reflection involution $\mathcal{B}_C(\mathbf{s}_0, -v_{\mathrm{coll}}\mathbf{s}_0) = (\mathbf{s}_0, +v_{\mathrm{coll}}\mathbf{s}_0)$ ($\mathcal{B}_C^2 = \mathrm{Id}$). Figure~\ref{fig:l3_mcgehee_central_collisions} separates the three canonical central-collision multiplicities into three independent columns over a complete infall--bounce--return cycle $[0, 2\,t_{\mathrm{coll}}]$ at $\Delta t = 0.025$, where the discrete circle (\texttt{Float64}) and triangle (\texttt{Float128}) markers in the bottom panels (d--f) evaluate the exact round-trip time-reversibility error $|\Delta H_{\mathrm{FB}}(t_{\mathrm{turn}})|/|H_0|$ after integrating forward $[0, t_{\mathrm{turn}}]$ and backward $[t_{\mathrm{turn}}, 0]$ back to $t = 0$ at four turnaround checkpoints straddling the collision time $t_{\mathrm{coll}}$: \begin{enumerate}[label=(\roman*)] \item \emph{Column~1 --- Binary head-on collision} ($k=2$, $r_0 = 2.0\text{ AU}$, $m_i = 1.0$, $t_{\mathrm{coll}}^{(k=2)} = \pi/\sqrt{2} = 2.2214$, Figure~\ref{fig:l3_mcgehee_central_collisions}a,d): Traverses $r = 0$ via the Sundman--Levi-Civita chart ($\mathcal{B}_C: -2.000 \mapsto +2.000$) with $\max |\Delta H|/|H_0| = 1.05 \times 10^{-13}$ in \texttt{Float64} and $6.60 \times 10^{-38}$ in \texttt{Float128}, and 100-step ($t_{\mathrm{turn}} = 2.50 > t_{\mathrm{coll}}^{(k=2)}$) round-trip reversibility ($\Phi_{-\Delta t}^{100} \circ \Phi_{\Delta t}^{100}$ crossing $r=0$ twice across $N_{\mathrm{eval}} = 1{,}424$ leaf evaluations) of $|\Delta H_{\mathrm{FB}}|/|H_0| = 2.48 \times 10^{-13} \le 2\,N_{\mathrm{eval}}\,\varepsilon_{\mathrm{mach}} = 3.16 \times 10^{-13}$ in \texttt{Float64} and $1.10 \times 10^{-33}$ in \texttt{Float128}; \item \emph{Column~2 --- Equilateral Lagrange homothetic triple collision} ($k=3$, $R_0 = 1.20\text{ AU}$, $m_i = 1.0$, $GM_{\mathrm{eff}} = 1/\sqrt{3}$, $t_{\mathrm{coll}}^{(k=3)} = 1.9216$, Figure~\ref{fig:l3_mcgehee_central_collisions}b,e): Crosses the 3-body collision divisor $\rho = 0$ ($\mathcal{B}_C: -1.075 \mapsto +1.075$) with $\max |\Delta H|/|H_0| = 1.04 \times 10^{-13}$ in \texttt{Float64} and $8.91 \times 10^{-34}$ in \texttt{Float128}, and 100-step ($t_{\mathrm{turn}} = 2.50 > t_{\mathrm{coll}}^{(k=3)}$) round-trip reversibility $|\Delta H_{\mathrm{FB}}|/|H_0| = 4.55 \times 10^{-14}$ in \texttt{Float64} and $5.65 \times 10^{-34}$ in \texttt{Float128} (maximum across all four checkpoints $1.47 \times 10^{-13}$ and $1.43 \times 10^{-33}$ at $N = 80$); \item \emph{Column~3 --- Symmetric square $D_4$ quadruple collision} ($k=4$, $R_0 = 1.00\text{ AU}$, $m_i = 1.0$, $GM_{\mathrm{eff}} = \frac{1}{4}+\frac{1}{\sqrt{2}} = 0.9571$, $t_{\mathrm{coll}}^{(k=4)} = 1.1353$, Figure~\ref{fig:l3_mcgehee_central_collisions}c,f): Integrated up to $T_{\max} = 2.50 > 2\,t_{\mathrm{coll}}^{(k=4)}$, genuinely crossing the simultaneous 4-body origin $\rho = 0$ ($\mathcal{B}_C: -1.384 \mapsto +1.384$) and rebounding with $\max |\Delta H|/|H_0| = 7.39 \times 10^{-14}$ in \texttt{Float64} and $1.46 \times 10^{-33}$ in \texttt{Float128}, 100-step ($t_{\mathrm{turn}} = 2.50$) round-trip reversibility across the collision of $|\Delta H_{\mathrm{FB}}|/|H_0| = 1.76 \times 10^{-13}$ (\texttt{Float64}) and $1.93 \times 10^{-33}$ (\texttt{Float128}, position residual $\|\Delta\mathbf{q}_{\mathrm{FB}}\|_2 = 3.52 \times 10^{-13}$ in \texttt{Float64} and $3.86 \times 10^{-33}$ in \texttt{Float128}), and dihedral $D_4$ point-group symmetry defect $\Delta_{D_4} = 7.00 \times 10^{-16}$ in \texttt{Float64} and $0.00$ in \texttt{Float128}. \end{enumerate} \begin{figure}[!htbp] \centering \includegraphics[width=\textwidth]{chapters/ch07_benchmarks_and_lean4/figures/fig_ch07_l3_mcgehee_central_collisions.png} \caption{\textbf{Limit Case L3: Exact McGehee exceptional divisor traversals ($\rho = 0$) and conormal bounce involution $\mathcal{B}_C^2 = \mathrm{Id}$ for central collisions $k \in \{2, 3, 4\}$.} \textbf{(a--c)}~McGehee blown-up phase portraits $(\rho, \dot{\rho}\sqrt{\rho})$ with infall ($\leftarrow$), vertical reflection $\mathcal{B}_C$ at $\rho = 0$ ($\uparrow$), and rebound ($\rightarrow$) for $k = 2$ Sundman, $k = 3$ Lagrange, and $k = 4$ square $D_4$. \textbf{(d--f)}~Dual-precision continuous Hamiltonian error $|\Delta H(t)|/|H_0|$ and round-trip reversibility markers $|\Delta H_{\mathrm{FB}}(t_{\mathrm{turn}})|/|H_0|$ across $t = t_{\mathrm{coll}}$.} \label{fig:l3_mcgehee_central_collisions} \end{figure} \medskip\noindent\textbf{Limit Case L4 : Adaptive Multiscale Resolution Without Exact Collision --- Chaotic Burrau 3-4-5 Scattering (L4a) and 3D Kozai--Lidov Secular Resonance (L4b) (Figure~\ref{fig:l4_adaptive_multiscale_burrau_kozai}).} In general non-homothetic $N$-body dynamics without a massive central host, extreme close encounters and secular resonances occur with small but nonzero impact parameters ($r_{\min} > 0$) and are resolved on the outer Fulton--MacPherson / Jacobi chart atlas by time-symmetric dyadic sub-stepping and shadow gauge connections (Sections~\ref{sec:dyadic_multirate_scheduling}, \ref{sec:dyadic_scheduler}, and~\ref{sub:synchronous_supernode_anti_zeno}): \begin{enumerate}[label=(\roman*)] \item \emph{L4a --- Pythagorean 3-4-5 chaotic three-body scattering of Burrau \cite{burrau1913berechnung}} ($m = \{3,4,5\}$, $t \in [0, 12.0]$, $\Delta t = 0.04$, Figure~\ref{fig:l4_adaptive_multiscale_burrau_kozai}a--b): Across successive chaotic triple flybys reaching $r_{\min} = 1.15 \times 10^{-2}\text{ AU}$, the dyadic scheduler dynamically ascends and descends $k_{\mathrm{sym}}(t) \in \{0, \dots, 7\}$ (zero saturation under $k_{\max} = 8$) in exact synchrony with $\min_{i