Imagen de Portada

A Universal Recursive Feedback Framework for Stabilizing Chaotic Dynamical Systems: Applications to Classical and Quantum Physics

Fuente: Zenodo

Guardado en:

Detalles Bibliográficos
Autor principal:	McEvoy, Adam L
Formato:	Recurso digital
Lenguaje:	inglés
Publicado:	Zenodo 2026
Materias:	dynamical systems chaos control feedback stabilization three-body problem quantum mechanics predictive control thedr
Acceso en línea:	Acceder al recurso
Etiquetas:	Agregar Etiqueta Sin Etiquetas, Sea el primero en etiquetar este registro!

_version_	1866901897336586240
author	McEvoy, Adam L
author_facet	McEvoy, Adam L
contents	<p># A Multi-Layer Fractal Correction Engine for Stabilizing Chaotic Dynamical Systems: Theory, Implementation, and Numerical Validation</p> <p>Author: Adam L McEvoy</p> <p>Date: March 2026</p> <p>Keywords: fractal correction engine, chaotic dynamical systems, curvature-based feedback, three-body problem, double pendulum, Lorenz attractor, energy conservation, Lyapunov exponent reduction, multi-layer correction, phase-space prediction</p> <p>---</p> <p>## Abstract</p> <p>I present the Fractal Correction Engine (FCE), a curvature-based recursive feedback framework that stabilizes chaotic dynamical systems by exploiting the fractal self-similarity of phase-space trajectories. The FCE operates on the principle that any orbit, wave, or waveform admits a pi-scaled curvature decomposition whose phase structure repeats with measurable regularity. By binning observed curvature values against a normalized phase coordinate $\varphi(t) = (\theta(t)/\pi) \bmod 2$, the engine constructs a predictive model of expected curvature at each phase and applies a perpendicular correction that steers the trajectory toward its own fractal attractor without directly injecting or removing kinetic energy.</p> <p>I extend the single-layer FCE to a multi-layer architecture in which secondary FCE instances monitor scalar features of the dynamics --- total energy $E(t)$, minimum inter-body distance $d_{\min}(t)$, and Lorenz lobe-transition timing $x(t)$ --- and modulate the primary correction in real time. This layered approach is motivated by the observation that the FCE sometimes needs to be applied a second time to features that affect the primary dynamics, analogous to the dual-correction architecture demonstrated in prior five-body celestial simulations.</p> <p>I validate the framework on three canonical chaotic systems: (i) the gravitational three-body problem initialized as a perturbed Chenciner--Montgomery figure-eight orbit, (ii) the double pendulum at high amplitude, and (iii) the Lorenz attractor at standard parameters. Across all three systems the multi-layer FCE produces measurable improvements: 92.5% reduction in energy drift for the double pendulum, 31.5% reduction in maximum energy excursion for the three-body problem, and 14.6% reduction in the maximal Lyapunov exponent for the Lorenz attractor. All simulations are fully reproducible from the companion code repository.</p> <p>---</p> <p>## 1. Introduction</p> <p>### 1.1 Background and Motivation</p> <p>Chaotic dynamical systems --- from planetary orbits to atmospheric convection --- exhibit sensitive dependence on initial conditions (SDIC), making long-term prediction exponentially difficult [1, 2]. The maximal Lyapunov exponent $\lambda_{\max}$ quantifies this sensitivity: two trajectories separated by an infinitesimal perturbation $\delta_0$ diverge as</p> <p>$$\\|\delta(t)\\| \sim \\|\delta_0\\| \, e^{\lambda_{\max} t}$$</p> <p>For $\lambda_{\max} > 0$ the system is chaotic, and the prediction horizon scales as $T_{\text{pred}} \sim \lambda_{\max}^{-1} \ln(\varepsilon / \\|\delta_0\\|)$ where $\varepsilon$ is an acceptable error tolerance [3].</p> <p>Traditional approaches to chaos control fall into two broad categories. OGY-type methods [4] stabilize unstable periodic orbits (UPOs) embedded in the chaotic attractor by applying small perturbations when the trajectory passes near a UPO. Pyragas delayed-feedback methods [5] use a continuous control signal proportional to $\mathbf{X}(t) - \mathbf{X}(t - \tau)$ where $\tau$ is the period of the target UPO. Both approaches require prior knowledge of the target orbit and are restricted to systems whose UPO structure is known.</p> <p>The three-body gravitational problem, first studied by Newton and proven by Poincare to lack general closed-form solutions [1], exemplifies the need for a more universal approach. While special solutions exist --- the Euler--Lagrange collinear configurations, the Chenciner--Montgomery figure-eight [6] --- generic initial conditions lead to chaotic scattering, close encounters, and eventual escape. No general-purpose stabilization method has been demonstrated for arbitrary three-body configurations.</p> <p>### 1.2 The Fractal Correction Engine</p> <p>The Fractal Correction Engine (FCE) introduced in this work takes a fundamentally different approach. Rather than targeting a known periodic orbit, the FCE learns the fractal structure of whatever trajectory the system is already following and uses that structure to predict and correct deviations in real time.</p> <p>The key insight is that the trajectories of chaotic systems, while unpredictable in detail, exhibit fractal self-similarity in their curvature profiles. When curvature is parameterized by a $\pi$-normalized phase variable, recurring patterns emerge that can be exploited for short-horizon prediction. The FCE applies corrections perpendicular to the velocity, which changes the trajectory's curvature without directly altering kinetic energy --- a critical property for Hamiltonian systems where energy conservation is a fundamental constraint.</p> <p>The FCE operates on any orbit, wave, wavelength, or waveform. It uses $\pi$ and local curvature to extract a fractal path that is identical to the observed path and can be used for forward and backward trajectory prediction and wave/interference mapping.</p> <p>### 1.3 Multi-Layer Architecture</p> <p>A single-layer FCE corrects the trajectory in position/velocity space. However, chaotic dynamics involve coupled features that each exhibit their own waveform structure: total energy oscillates, inter-body distances fluctuate, and regime transitions (such as Lorenz lobe switches) occur at semi-regular intervals. Each of these features can itself be treated as a waveform to which the FCE is applied.</p> <p>The multi-layer FCE applies secondary FCE instances to these scalar features:<br>- Layer 1 (Trajectory FCE): The primary curvature-based correction in state space<br>- Layer 2 (Energy FCE): A scalar FCE tracking $E(t) - E_0$ as a 1D waveform<br>- Layer 2b (Distance FCE): A scalar FCE tracking $d_{\min}(t)$ for multi-body systems<br>- Layer 2c (Lobe-Transition FCE): A scalar FCE tracking $x(t)$ for the Lorenz system</p> <p>The outputs of the secondary layers modulate the primary correction --- dampening it during sensitive periods (close encounters, lobe transitions) and strengthening energy projection when drift accelerates. This architecture was motivated by the dual-correction approach successfully deployed in prior five-body celestial simulations [7], where a primary 3-body correction and a secondary pair correction were applied to different gravitational subgroups detected via hierarchical clustering.</p> <p>### 1.4 Paper Structure</p> <p>Section 2 presents the complete mathematical formulation of the FCE and its multi-layer extension. Section 3 describes the three test systems and their equations of motion. Section 4 details the computational implementation. Section 5 presents numerical results with quantitative comparisons. Section 6 discusses implications and limitations. Section 7 concludes with future directions.</p> <p>---</p> <p>## 2. Mathematical Formulation</p> <p>### 2.1 Governing Equation</p> <p>Consider a dynamical system evolving in $\mathbb{R}^n$:</p> <p>$$\frac{d\mathbf{X}}{dt} = \mathbf{F}(\mathbf{X}, t)$$</p> <p>where $\mathbf{X}(t) \in \mathbb{R}^n$ is the state vector and $\mathbf{F}$ encodes the intrinsic dynamics (gravitational forces, Lagrangian mechanics, etc.). The FCE augments this with a correction term:</p> <p>$$\frac{d\mathbf{X}}{dt} = \mathbf{F}(\mathbf{X}, t) + \hat{\mathbf{C}}(t)$$</p> <p>where $\hat{\mathbf{C}}(t) \in \mathbb{R}^n$ is the fractal correction vector. The remainder of this section specifies how $\hat{\mathbf{C}}(t)$ is computed from the trajectory's own curvature structure.</p> <p>### 2.2 Local Curvature Extraction</p> <p>Let $\mathbf{v}(t) = d\mathbf{X}/dt$ be the velocity and $\mathbf{a}(t) = d^2\mathbf{X}/dt^2$ the acceleration. The signed curvature (2D) and unsigned curvature (3D) are:</p> <p>$$\kappa_{2D}(t) = \frac{v_x \, a_y - v_y \, a_x}{\\|\mathbf{v}\\|^3}$$</p> <p>$$\kappa_{3D}(t) = \frac{\\|\mathbf{v} \times \mathbf{a}\\|}{\\|\mathbf{v}\\|^3}$$</p> <p>These are the standard Frenet--Serret curvatures of the trajectory as a space curve. The curvature $\kappa(t)$ encodes the instantaneous rate of turning per unit arc length and is invariant under reparameterization of time --- it is a geometric property of the path itself.</p> <p>### 2.3 Pi-Scaled Phase Variable</p> <p>The cumulative turning angle is defined as:</p> <p>$$\theta(t) = \int_0^t \|\kappa(\tau)\| \, \\|\mathbf{v}(\tau)\\| \, d\tau$$</p> <p>This integral computes the total angle swept by the tangent vector along the trajectory. We define the fractal phase as:</p> <p>$$\varphi(t) = \frac{\theta(t)}{\pi} \bmod 2 \quad \in [0, 2)$$</p> <p>The normalization by $\pi$ is fundamental: $\pi$ radians represents a half-turn, so $\varphi$ completes one full cycle ($0 \to 2$) for every full revolution ($2\pi$ radians) of the tangent vector. For quasi-periodic orbits this phase variable advances nearly uniformly. For chaotic orbits the advance rate fluctuates, but the phase-curvature relationship retains statistical regularity due to the fractal self-similarity of the attractor.</p> <p>### 2.4 Phase-Binned Curvature Prediction</p> <p>The phase domain $[0, 2)$ is discretized into $N_b$ bins of width $\Delta\varphi = 2/N_b$. At each timestep, the current curvature $\kappa(t)$ is deposited into the bin corresponding to $\varphi(t)$, with exponential decay weighting:</p> <p>$$B_k^{(\text{sum})} \leftarrow \gamma \, B_k^{(\text{sum})} + \kappa(t) \, \delta_{k, b(t)}$$</p> <p>$$B_k^{(\text{weight})} \leftarrow \gamma \, B_k^{(\text{weight})} + \delta_{k, b(t)}$$</p> <p>where $b(t) = \lfloor \varphi(t) \cdot N_b / 2 \rfloor$ is the bin index, $\gamma = 0.995$ is the decay factor, and $\delta_{k,b}$ is the Kronecker delta. The predicted curvature at phase $\varphi$ is:</p> <p>$$\kappa_{\text{pred}}(\varphi) = \frac{B_{b(\varphi)}^{(\text{sum})}}{B_{b(\varphi)}^{(\text{weight})}}$$</p> <p>If the primary bin has insufficient samples ($B_k^{(\text{count})} < 2$), neighboring bins are queried with offsets $\{0, \pm 1, \pm 2\}$. If no bin has sufficient data, an exponentially-weighted moving average of the curvature history is used as fallback:</p> <p>$$\kappa_{\text{pred}} = \frac{\sum_{i=1}^{W} w_i \, \kappa_{t-i}}{\sum_{i=1}^{W} w_i}, \quad w_i = e^{-3(W-i)/W}$$</p> <p>where $W$ is the curvature window length.</p> <p>### 2.5 Perpendicular Correction Vector</p> <p>The correction is applied in the direction perpendicular to the velocity. Define the unit normal vector:</p> <p>In 2D:</p> <p>$$\hat{\mathbf{n}} = \frac{1}{\\|\mathbf{v}\\|} \begin{pmatrix} -v_y \\ v_x \end{pmatrix}$$</p> <p>In 3D:</p> <p>$$\hat{\mathbf{n}} = \frac{\mathbf{a}_\perp}{\\|\mathbf{a}_\perp\\|}, \quad \mathbf{a}_\perp = \mathbf{a} - (\mathbf{a} \cdot \hat{\mathbf{v}})\,\hat{\mathbf{v}}$$</p> <p>where $\hat{\mathbf{v}} = \mathbf{v}/\\|\mathbf{v}\\|$. If $\\|\mathbf{a}_\perp\\| < \varepsilon$, a reference vector is used to construct $\hat{\mathbf{n}}$ via cross product.</p> <p>The correction vector is:</p> <p>$$\hat{\mathbf{C}}(t) = \alpha \, \Delta\kappa \, \\|\mathbf{v}\\|^2 \, \hat{\mathbf{n}}$$</p> <p>where $\alpha$ is the correction strength and $\Delta\kappa = \kappa_{\text{pred}} - \kappa$ is the curvature mismatch.</p> <p>Perpendicularity and energy preservation. Since $\hat{\mathbf{C}} \perp \mathbf{v}$, the instantaneous power delivered by the correction is:</p> <p>$$P = \hat{\mathbf{C}} \cdot \mathbf{v} = \alpha \, \Delta\kappa \, \\|\mathbf{v}\\|^2 \, (\hat{\mathbf{n}} \cdot \mathbf{v}) = 0$$</p> <p>The correction changes the direction of the velocity (and hence the trajectory curvature) without directly injecting or removing kinetic energy. For Hamiltonian systems this is essential: it means the FCE modifies the orbit geometry while leaving the energy surface largely intact.</p> <p>### 2.6 Recursive Refinement</p> <p>The correction is computed iteratively over $M$ recursive iterations:</p> <p>$$\hat{\mathbf{C}}^{(0)} = \mathbf{0}$$</p> <p>$$\hat{\mathbf{C}}^{(j+1)} = \hat{\mathbf{C}}^{(j)} + \alpha \, (\kappa_{\text{pred}} - \kappa^{(j)}) \, \\|\mathbf{v}\\|^2 \, \hat{\mathbf{n}}$$</p> <p>where $\kappa^{(j)} = \kappa(\mathbf{v}, \mathbf{a} + \hat{\mathbf{C}}^{(j)})$ is the curvature recomputed with the accumulated correction. This iterative scheme converges the actual curvature toward the predicted value.</p> <p>### 2.7 Correction Clamping</p> <p>To prevent numerical instability, the correction magnitude is clamped:</p> <p>$$\\|\hat{\mathbf{C}}\\| \leq \min\left(r \\|\mathbf{a}\\|, \; C_{\max}\right)$$</p> <p>where $r$ is the maximum correction-to-acceleration ratio (typically 0.05) and $C_{\max}$ is an absolute ceiling. This ensures corrections remain small perturbations to the natural dynamics.</p> <p>### 2.8 Scalar FCE for 1D Waveform Tracking</p> <p>For scalar signals $s(t)$ (energy, distance, etc.), we define the 1D curvature as:</p> <p>$$\kappa_{1D}(t) = \frac{\|s''(t)\|}{(1 + s'(t)^2)^{3/2}}$$</p> <p>This is the standard curvature of the graph $(t, s(t))$ as a plane curve. The phase variable, bin accumulation, and prediction follow identically to the vector case. The output is a scalar correction factor $c_s \in [-c_{\max}, c_{\max}]$:</p> <p>$$c_s(t) = \alpha_s \, (\kappa_{\text{pred}} - \kappa_{1D}) \, \text{sgn}(s'(t))$$</p> <p>The sign of $s'(t)$ ensures the correction opposes the predicted curvature deviation in the correct direction.</p> <p>### 2.9 Multi-Layer Architecture</p> <p>The full multi-layer correction at each timestep proceeds as:</p> <p>Step 1 --- Layer 2b (Distance modulation, N-body only):</p> <p>$$\mu_d = \text{clip}\!\left(1 + c_{\text{dist}},\; 0.2,\; 1.5\right)$$</p> <p>where $c_{\text{dist}}$ is the scalar FCE output tracking $d_{\min}(t)$.</p> <p>Step 2 --- Layer 1 (Primary trajectory correction):</p> <p>$$\hat{\mathbf{C}}_i = \mu_d \cdot \text{FCE}_i(\mathbf{v}_i, \mathbf{a}_i, \Delta t) \quad \text{for each body } i$$</p> <p>For multi-body systems, momentum conservation is enforced:</p> <p>$$\hat{\mathbf{C}}_i \leftarrow \hat{\mathbf{C}}_i - \frac{1}{M_{\text{total}}} \sum_j m_j \hat{\mathbf{C}}_j$$</p> <p>ensuring $\sum_i m_i \hat{\mathbf{C}}_i = \mathbf{0}$, which prevents center-of-mass drift.</p> <p>Step 3 --- Integration:</p> <p>$$\mathbf{X}(t + \Delta t) = \text{RK4}\!\left[\mathbf{F}(\mathbf{X}) + \hat{\mathbf{C}},\; \Delta t\right]$$</p> <p>Step 4 --- Layer 2 (Energy sub-FCE):</p> <p>$$e_{\text{hint}} = \text{ScalarFCE}\!\left(E(t+\Delta t) - E_0,\; \Delta t\right)$$</p> <p>The energy hint adaptively controls the projection blend:</p> <p>$$\beta = \beta_0 + \beta_1 \, \min\!\left(\frac{\|e_{\text{hint}}\|}{e_{\text{threshold}}},\; 1\right)$$</p> <p>Step 5 --- Energy projection (Hamiltonian systems only):</p> <p>The velocities are gently rescaled toward energy conservation:</p> <p>$$T_{\text{target}} = E_0 - V(\mathbf{X}(t+\Delta t))$$</p> <p>$$\mathbf{v} \leftarrow \mathbf{v} \cdot \left[1 + \beta\left(\sqrt{\frac{T_{\text{target}}}{T_{\text{current}}}} - 1\right)\right]$$</p> <p>where $T$ denotes kinetic energy and $\beta \in [0.02, 0.10]$ is the blend factor. The gentle blending avoids discrete jumps while steadily correcting energy drift.</p> <p>Step 6 --- Layer 2c (Lobe-transition dampening, Lorenz only):</p> <p>When the Lorenz $x$-component crosses zero (lobe transition), the primary correction is dampened for a cooldown window:</p> <p>$$\mu_{\text{lobe}} = \begin{cases} 0.3 + 0.7(1 - k/K) & \text{if within cooldown window} \\ \text{clip}(1 - \|c_{\text{lobe}}\|, 0.2, 1.0) & \text{otherwise} \end{cases}$$</p> <p>where $k$ counts steps since the last transition and $K = 50$ is the cooldown duration.</p> <p>### 2.10 Warmup Phase</p> <p>During an initial warmup period of $N_w$ timesteps, the FCE accumulates curvature and phase data without applying any corrections ($\hat{\mathbf{C}} = \mathbf{0}$). This allows the phase bins to populate with representative curvature values before prediction begins. The warmup length depends on the system's characteristic timescale: for the three-body figure-eight orbit (period $T \approx 6.3$), we use $N_w = 8000$ steps at $\Delta t = 0.001$ (approximately 1.27 orbital periods).</p> <p>---</p> <p>## 3. Test Systems</p> <p>### 3.1 Three-Body Gravitational Problem</p> <p>The gravitational three-body problem governs the motion of three point masses $m_1, m_2, m_3$ interacting via Newtonian gravity. The state vector $\mathbf{X} = (\mathbf{r}_1, \mathbf{r}_2, \mathbf{r}_3, \mathbf{v}_1, \mathbf{v}_2, \mathbf{v}_3) \in \mathbb{R}^{12}$ evolves according to:</p> <p>$$\frac{d\mathbf{r}_i}{dt} = \mathbf{v}_i, \quad \frac{d\mathbf{v}_i}{dt} = \sum_{j \neq i} \frac{G \, m_j \, (\mathbf{r}_j - \mathbf{r}_i})}{(\\|\mathbf{r}_j - \mathbf{r}_i\\|^2 + \epsilon^2)^{3/2}}$$</p> <p>where $G$ is the gravitational constant and $\epsilon = 0.02$ is a Plummer softening parameter that prevents numerical divergence during close encounters.</p> <p>The total energy is conserved:</p> <p>$$E = \sum_i \frac{1}{2} m_i \\|\mathbf{v}_i\\|^2 - \sum_{i < j} \frac{G \, m_i \, m_j}{\\|\mathbf{r}_j - \mathbf{r}_i\\|}$$</p> <p>as is the total angular momentum:</p> <p>$$L = \sum_i m_i \, (\mathbf{r}_i \times \mathbf{v}_i)$$</p> <p>Initial conditions. We use the Chenciner--Montgomery figure-eight solution [6] with equal masses $m_1 = m_2 = m_3 = 1$:</p> <p>\| Body \| Position $(x, y)$ \| Velocity $(v_x, v_y)$ \|<br>\|------\|-------------------\|----------------------\|<br>\| 1 \| $(-0.97000436, \; 0.24308753)$ \| $(0.4662036850, \; 0.4323657300)$ \|<br>\| 2 \| $(0, \; 0)$ \| $(-0.93240737, \; -0.86473146)$ \|<br>\| 3 \| $(0.97000436, \; -0.24308753)$ \| $(0.4862036850, \; 0.4123657300)$ \|</p> <p>Body 3's velocity is perturbed by $(+0.02, -0.02)$ from the exact figure-eight solution to seed chaotic divergence. The simulation runs for 80,000 timesteps at $\Delta t = 0.001$ (80 time units, approximately 12.7 figure-eight periods).</p> <p>### 3.2 Double Pendulum</p> <p>The double pendulum consists of two point masses $m_1, m_2$ on rigid massless rods of lengths $l_1, l_2$, subject to gravitational acceleration $g$. The state vector is $\mathbf{X} = (\theta_1, \theta_2, \dot\theta_1, \dot\theta_2)$ with equations of motion derived from the Lagrangian $\mathcal{L} = T - V$:</p> <p>$$T = \frac{1}{2}(m_1 + m_2) l_1^2 \dot\theta_1^2 + \frac{1}{2} m_2 l_2^2 \dot\theta_2^2 + m_2 l_1 l_2 \dot\theta_1 \dot\theta_2 \cos(\theta_1 - \theta_2)$$</p> <p>$$V = -(m_1 + m_2) g l_1 \cos\theta_1 - m_2 g l_2 \cos\theta_2$$</p> <p>The angular accelerations are obtained from the Euler--Lagrange equations:</p> <p>$$\ddot\theta_1 = \frac{-g(2m_1 + m_2)\sin\theta_1 - m_2 g \sin(\theta_1 - 2\theta_2) - 2\sin(\theta_1 - \theta_2)\,m_2(\dot\theta_2^2 l_2 + \dot\theta_1^2 l_1 \cos(\theta_1 - \theta_2))}{l_1(2m_1 + m_2 - m_2 \cos(2(\theta_1 - \theta_2)))}$$</p> <p>$$\ddot\theta_2 = \frac{2\sin(\theta_1 - \theta_2)\left(\dot\theta_1^2 l_1 (m_1 + m_2) + g(m_1 + m_2)\cos\theta_1 + \dot\theta_2^2 l_2 m_2 \cos(\theta_1 - \theta_2)\right)}{l_2(2m_1 + m_2 - m_2\cos(2(\theta_1 - \theta_2)))}$$</p> <p>Parameters: $m_1 = m_2 = 1$ kg, $l_1 = l_2 = 1$ m, $g = 9.81$ m/s$^2$, $\theta_1(0) = \pi/2$, $\theta_2(0) = \pi$, $\dot\theta_1(0) = \dot\theta_2(0) = 0$. The simulation runs 60,000 steps at $\Delta t = 0.001$ s.</p> <p>The system is Hamiltonian with conserved total energy $E = T + V$. For this initial condition (near-vertical pendulum 2), the Lyapunov exponent is positive and the motion is strongly chaotic [8, 9].</p> <p>### 3.3 Lorenz Attractor</p> <p>The Lorenz system [2] is the canonical example of a dissipative chaotic system:</p> <p>$$\dot{x} = \sigma(y - x)$$</p> <p>$$\dot{y} = x(\rho - z) - y$$</p> <p>$$\dot{z} = xy - \beta z$$</p> <p>with standard parameters $\sigma = 10$, $\rho = 28$, $\beta = 8/3$ and initial condition $(x_0, y_0, z_0) = (1, 1, 1)$.</p> <p>Unlike the Hamiltonian systems above, the Lorenz system is dissipative: the phase-space volume contracts at rate $\nabla \cdot \mathbf{F} = -(\sigma + 1 + \beta) = -13.67$. There is no conserved energy. Instead, the trajectory asymptotically approaches the Lorenz attractor, a fractal set with Hausdorff dimension $D \approx 2.06$ [10].</p> <p>For exact second derivatives (needed by the FCE), we compute the acceleration analytically via the Jacobian:</p> <p>$$\frac{d^2\mathbf{X}}{dt^2} = J(\mathbf{X}) \cdot \mathbf{F}(\mathbf{X}), \quad J = \begin{pmatrix} -\sigma & \sigma & 0 \\ \rho - z & -1 & -x \\ y & x & -\beta \end{pmatrix}$$</p> <p>This avoids the numerical noise inherent in finite-difference approximations of the acceleration and is essential for stable curvature computation.</p> <p>The simulation runs 30,000 steps at $\Delta t = 0.005$ (150 time units).</p> <p>---</p> <p>## 4. Computational Implementation</p> <p>### 4.1 Integration Scheme</p> <p>All systems are integrated with the classical 4th-order Runge--Kutta (RK4) method. For the corrected trajectory, the correction vector is held constant across the four RK4 substeps:</p> <p>$$k_1 = \mathbf{F}(\mathbf{X}_n) + \hat{\mathbf{C}}_n$$<br>$$k_2 = \mathbf{F}(\mathbf{X}_n + \tfrac{1}{2}\Delta t \, k_1) + \hat{\mathbf{C}}_n$$<br>$$k_3 = \mathbf{F}(\mathbf{X}_n + \tfrac{1}{2}\Delta t \, k_2) + \hat{\mathbf{C}}_n$$<br>$$k_4 = \mathbf{F}(\mathbf{X}_n + \Delta t \, k_3) + \hat{\mathbf{C}}_n$$<br>$$\mathbf{X}_{n+1} = \mathbf{X}_n + \frac{\Delta t}{6}(k_1 + 2k_2 + 2k_3 + k_4)$$</p> <p>Holding $\hat{\mathbf{C}}$ constant over the substeps ensures the correction acts as a smooth external forcing rather than a state-dependent feedback within the integration step.</p> <p>### 4.2 FCE Parameters</p> <p>The FCE parameters were tuned for each system based on characteristic timescales and dynamics:</p> <p>Table 1: FCE Configuration Parameters</p> <p>\| Parameter \| Three-Body \| Double Pendulum \| Lorenz \|<br>\|-----------\|-----------\|----------------\|--------\|<br>\| Dimension \| 2 (per body) \| 2 (angle space) \| 3 \|<br>\| Correction strength $\alpha$ \| 0.005 \| 0.03 \| 0.008 \|<br>\| Phase bins $N_b$ \| 64 \| 32 \| 48 \|<br>\| Warmup steps $N_w$ \| 8000 \| 200 \| 300 \|<br>\| Curvature window $W$ \| 1000 \| 500 \| 800 \|<br>\| Recursive iterations $M$ \| 3 \| 3 \| 3 \|<br>\| Max absolute correction $C_{\max}$ \| 0.015 \| 2.0 \| 0.5 \|</p> <p>Table 2: Scalar Sub-FCE Configuration</p> <p>\| Sub-FCE Layer \| System \| Strength \| Bins \| Warmup \| Max Output \|<br>\|--------------\|--------\|----------\|------\|--------\|------------\|<br>\| Energy (Layer 2) \| Three-Body \| 0.003 \| 32 \| 8000 \| 0.01 \|<br>\| Energy (Layer 2) \| Double Pendulum \| 0.008 \| 32 \| 200 \| 0.03 \|<br>\| Distance (Layer 2b) \| Three-Body \| 0.008 \| 32 \| 8000 \| 0.4 \|<br>\| Lobe (Layer 2c) \| Lorenz \| 0.01 \| 32 \| 300 \| 0.8 \|</p> <p>### 4.3 Three-Body Momentum Conservation</p> <p>For the three-body system, each body has its own FCE instance operating in 2D position space. After computing per-body corrections $\hat{\mathbf{C}}_i$, momentum conservation is enforced:</p> <p>$$\hat{\mathbf{C}}_i \leftarrow \hat{\mathbf{C}}_i - \frac{\sum_j m_j \hat{\mathbf{C}}_j}{M_{\text{total}}}$$</p> <p>This subtraction ensures $\sum_i m_i \hat{\mathbf{C}}_i = \mathbf{0}$, preventing secular drift of the center of mass. Without this constraint, even small per-body corrections accumulate into a net momentum that causes the entire system to translate, eventually leading to apparent escape.</p> <p>### 4.4 Double Pendulum Angle-Space Correction</p> <p>For the double pendulum, the FCE operates in angle space $(\theta_1, \theta_2)$ rather than Cartesian space. The velocity vector is $(\dot\theta_1, \dot\theta_2)$ and the acceleration vector is $(\ddot\theta_1, \ddot\theta_2)$. The correction is added to the angular acceleration components:</p> <p>$$\text{correction vector} = (0, \; 0, \; C_1, \; C_2)$$</p> <p>applied to the state vector $(\theta_1, \theta_2, \dot\theta_1, \dot\theta_2)$.</p> <p>### 4.5 Lorenz Analytical Acceleration</p> <p>For the Lorenz system, finite-difference approximation of $\ddot{\mathbf{X}}$ introduces numerical noise at the level of $O(\Delta t^{-1})$, which contaminates curvature estimates. Instead, we compute the exact acceleration:</p> <p>$$\ddot{\mathbf{X}} = \frac{d\mathbf{F}}{dt} = \frac{\partial \mathbf{F}}{\partial \mathbf{X}} \cdot \frac{d\mathbf{X}}{dt} = J(\mathbf{X}) \cdot \mathbf{F}(\mathbf{X})$$</p> <p>This analytical approach eliminates numerical differentiation entirely and is critical for stable curvature computation in the Lorenz system.</p> <p>### 4.6 Lyapunov Exponent Estimation</p> <p>The maximal Lyapunov exponent is estimated via the standard twin-trajectory method [11]:</p> <p>1. Initialize two trajectories $\mathbf{X}_1(0)$ and $\mathbf{X}_2(0) = \mathbf{X}_1(0) + \delta_0 \hat{\mathbf{e}}_1$ with $\delta_0 = 10^{-8}$<br>2. Integrate both trajectories for $R$ steps (renormalization interval)<br>3. Compute the divergence: $d = \\|\mathbf{X}_2 - \mathbf{X}_1\\|$<br>4. Accumulate: $\Lambda \leftarrow \Lambda + \ln(d / \delta_0)$<br>5. Renormalize: $\mathbf{X}_2 \leftarrow \mathbf{X}_1 + \delta_0 (\mathbf{X}_2 - \mathbf{X}_1) / d$<br>6. Repeat; the Lyapunov exponent converges to $\lambda_{\max} = \Lambda / (N \cdot R \cdot \Delta t)$</p> <p>For corrected trajectories of the double pendulum and Lorenz systems, each twin has its own independent FCE instance. For the three-body system, an energy-ratio proxy is used because the FCE's internal state (phase bin history) becomes inconsistent after the renormalization step that moves $\mathbf{X}_2$ to a new location while leaving its FCE state unchanged.</p> <p>---</p> <p>## 5. Results</p> <p>### 5.1 Overview</p> <p>Each system was simulated twice: once with the uncorrected dynamics $\mathbf{F}(\mathbf{X})$ alone, and once with the full multi-layer FCE correction $\mathbf{F}(\mathbf{X}) + \hat{\mathbf{C}}(t)$. Both simulations use identical initial conditions, timesteps, and integration methods. All results are deterministic and fully reproducible from the companion code.</p> <p>Table 3: Summary of Multi-Layer FCE Results</p> <p>\| System \| Metric \| Uncorrected \| FCE Corrected \| Improvement \|<br>\|--------\|--------\|------------\|--------------\|-------------\|<br>\| Three-Body \| Energy max drift \| $1.137 \times 10^{1}$ \| $7.791 \times 10^{0}$ \| +31.5% \|<br>\| \| Energy std / $\|E_0\|$ \| $6.521 \times 10^{-2}$ \| $5.559 \times 10^{-2}$ \| +14.7% \|<br>\| \| Lyapunov $\lambda_{\max}$ \| 0.1908 \| 0.1626 \| +14.7% \|<br>\| Double Pendulum \| Energy max drift \| $7.209 \times 10^{0}$ \| $5.122 \times 10^{-1}$ \| +92.9% \|<br>\| \| Energy std / $\|E_0\|$ \| $2.000 \times 10^{0}$ \| $1.508 \times 10^{-1}$ \| +92.5% \|<br>\| \| Lyapunov $\lambda_{\max}$ \| 0.2262 \| 0.1388 \| +38.6% \|<br>\| Lorenz \| Lyapunov $\lambda_{\max}$ \| 0.8415 \| 0.7186 \| +14.6% \|<br>\| \| $\\|X\\|^2$ std / $\\|X_0\\|^2$ \| $1.662 \times 10^2$ \| $1.640 \times 10^2$ \| +1.3% \|</p> <p>### 5.2 Three-Body Problem Results</p> <p>The perturbed figure-eight orbit is a stringent test: the perturbation of $(+0.02, -0.02)$ on body 3 drives the system away from the exact periodic solution into chaotic three-body dynamics over roughly 80 time units.</p> <p>Energy conservation. The uncorrected trajectory exhibits a maximum relative energy drift of $1.137 \times 10^{1}$ (the system eventually undergoes a close encounter that produces a large energy excursion due to the softening approximation). The FCE-corrected trajectory reduces this to $7.791 \times 10^{0}$, a 31.5% improvement in the worst-case energy excursion. The energy standard deviation (normalized by $\|E_0\|$) improves from $6.521 \times 10^{-2}$ to $5.559 \times 10^{-2}$, a 14.7% improvement.</p> <p>Active FCE layers:<br>- Layer 1: Per-body trajectory FCE (3 independent 2D engines) with momentum-conserving projection<br>- Layer 2: Energy waveform sub-FCE tracking $E(t) - E_0$, adaptively controlling projection blend ($\beta \in [0.02, 0.10]$)<br>- Layer 2b: Distance sub-FCE tracking $d_{\min}(t)$, modulating primary correction strength ($\mu_d \in [0.2, 1.5]$)<br>- Energy projection: velocity rescaling toward $E_0$ with adaptive blend</p> <p>Correction statistics: Mean correction magnitude $1.148 \times 10^{-3}$, maximum $1.435 \times 10^{-2}$. The corrections are 2--3 orders of magnitude smaller than the gravitational accelerations, confirming that the FCE operates as a small perturbation.</p> <p>### 5.3 Double Pendulum Results</p> <p>The double pendulum at $(\theta_1, \theta_2) = (\pi/2, \pi)$ exhibits vigorous chaotic motion with frequent over-the-top rotations. This system produces the most dramatic FCE improvement.</p> <p>Energy conservation. The uncorrected RK4 integration exhibits substantial energy drift: the maximum relative drift reaches $7.209$ (a 720% excursion from initial energy) and the energy standard deviation is $2.000 \times \|E_0\|$. The multi-layer FCE reduces these to $5.122 \times 10^{-1}$ and $1.508 \times 10^{-1}$ respectively --- improvements of 92.9% and 92.5%.</p> <p>This dramatic improvement results from the synergy of three mechanisms:<br>1. Layer 1 (trajectory FCE) steers the angular velocities toward the phase-predicted curvature<br>2. Layer 2 (energy sub-FCE) detects when the energy waveform is about to spike and pre-emptively rescales angular velocities<br>3. Energy projection gently nudges velocities toward the energy surface $E = E_0$ at each step</p> <p>The energy projection alone would force conservation but could distort the dynamics. By combining it with the curvature-based trajectory correction and the energy sub-FCE's predictive modulation, the system maintains both energy conservation and trajectory fidelity.</p> <p>Lyapunov exponent. The maximal Lyapunov exponent is reduced from $0.2262$ to $0.1388$, a 38.6% reduction. While the corrected system remains chaotic ($\lambda_{\max} > 0$), the reduced Lyapunov exponent means the prediction horizon is extended by a factor of approximately $0.2262 / 0.1388 = 1.63$ --- a 63% longer window of useful prediction.</p> <p>### 5.4 Lorenz Attractor Results</p> <p>The Lorenz system is dissipative, so energy conservation is not applicable as a metric. Instead, we focus on Lyapunov exponent reduction and trajectory structure preservation.</p> <p>Lyapunov exponent. The maximal Lyapunov exponent is reduced from $0.8415$ to $0.7186$, a 14.6% reduction. The prediction horizon is extended by a factor of $0.8415 / 0.7186 = 1.17$.</p> <p>Lobe-transition awareness. The Lorenz attractor consists of two lobes centered at the fixed points $(\pm\sqrt{\beta(\rho-1)}, \pm\sqrt{\beta(\rho-1)}, \rho-1) \approx (\pm 8.49, \pm 8.49, 27)$. Transitions between lobes (when $x$ crosses zero) are the most sensitive moments in the dynamics --- small perturbations here determine which lobe the trajectory enters next. The Layer 2c sub-FCE detects these transitions and dampens the primary correction during a 50-step cooldown window. This prevents the FCE from artificially forcing the trajectory into a particular lobe, preserving the natural attractor structure while still reducing overall divergence.</p> <p>Correction statistics: Mean correction magnitude $3.749 \times 10^{-1}$, maximum $5.000 \times 10^{-1}$ (hitting the absolute ceiling). The Lorenz system's rapid dynamics require relatively larger corrections compared to the gentler orbital mechanics of the three-body problem.</p> <p>### 5.5 Comparison: Single-Layer vs. Multi-Layer FCE</p> <p>To quantify the contribution of the layered architecture, we compare against the single-layer (Layer 1 only) FCE from our earlier work:</p> <p>Table 4: Single-Layer vs. Multi-Layer FCE Improvement</p> <p>\| System \| Metric \| Single-Layer \| Multi-Layer \| Gain \|<br>\|--------\|--------\|-------------\|-------------\|------\|<br>\| Three-Body \| Energy std improvement \| 12.8% \| 14.7% \| +15% relative \|<br>\| Double Pendulum \| Energy std improvement \| 28.8% \| 92.5% \| +221% relative \|<br>\| Double Pendulum \| Lyapunov reduction \| 38.6% \| 38.6% \| same \|<br>\| Lorenz \| Lyapunov reduction \| 14.6% \| 14.6% \| same \|</p> <p>The dominant gain is in the double pendulum's energy conservation, where the addition of the energy sub-FCE and energy projection produces a 3.2-fold improvement over the single-layer approach. The Lyapunov reductions are maintained but not significantly enhanced by the secondary layers, suggesting that trajectory divergence is primarily governed by the Layer 1 curvature correction, while energy drift is a separate phenomenon best addressed by dedicated monitoring.</p> <p>---</p> <p>## 6. Discussion</p> <p>### 6.1 Why Perpendicular Corrections Work</p> <p>The perpendicular correction is the core design principle of the FCE. By applying force normal to the velocity:</p> <p>1. Kinetic energy is preserved instantaneously: $\hat{\mathbf{C}} \cdot \mathbf{v} = 0$, so no work is done<br>2. Curvature is directly controlled: A normal force changes the direction of motion, which is precisely the curvature<br>3. Natural dynamics are minimally perturbed: The correction does not fight the system's own forces; it redirects them</p> <p>This is analogous to how a magnetic field deflects charged particles without doing work: $\mathbf{F}_B = q\mathbf{v} \times \mathbf{B} \perp \mathbf{v}$.</p> <p>### 6.2 The Role of Pi in Phase Normalization</p> <p>The choice to normalize by $\pi$ (rather than $2\pi$ or an arbitrary constant) is deliberate. The curvature phase $\varphi = \theta/\pi \bmod 2$ completes one cycle for every full $2\pi$-radian revolution of the tangent vector. This means:</p> <p>- Closed orbits have integer-valued $\varphi$ periods<br>- Quasi-periodic orbits produce dense but structured phase-curvature maps<br>- Chaotic orbits produce fractal phase-curvature distributions with self-similar structure</p> <p>The $\pi$-normalization connects the phase variable to the fundamental geometry of curves: $\pi$ radians is a half-turn, the minimal unit of directional reversal.</p> <p>### 6.3 Hierarchical Application of the FCE</p> <p>The multi-layer architecture demonstrates a key principle: the FCE can be applied recursively to any feature that exhibits waveform structure. The energy $E(t)$ of a chaotic Hamiltonian system is not constant (due to numerical integration error) but oscillates in a quasi-periodic pattern. The distance $d_{\min}(t)$ between bodies varies as they orbit and undergo close encounters. The Lorenz $x(t)$ traces out an irregular oscillation between positive and negative values.</p> <p>Each of these scalar signals has its own curvature structure that the ScalarFCE can learn and predict. The outputs then feed back into the primary correction loop, creating a hierarchy:</p> <p>$$\text{Scalar features} \xrightarrow{\text{ScalarFCE}} \text{Modulation signals} \xrightarrow{\text{multiply/blend}} \text{Primary correction} \xrightarrow{\text{RK4}} \text{Next state}$$</p> <p>This is directly analogous to the dual-correction architecture used in five-body celestial simulations [7], where hierarchical clustering identifies gravitationally bound subgroups and applies targeted corrections to each --- a 3:2 grouping might switch to 2:3 as the bodies rearrange, with the FCE predicting the transition.</p> <p>### 6.4 Limitations</p> <p>1. Warmup requirement. The FCE requires a learning phase before corrections begin. For the three-body system, this is 8000 steps --- 10% of the total simulation. During this time, the system evolves uncorrected and may already accumulate significant numerical error.</p> <p>2. Parameter sensitivity. The correction strength $\alpha$, phase bins $N_b$, and warmup length $N_w$ require system-specific tuning. While the general algorithm is universal, optimal parameters depend on the characteristic timescales and curvature magnitudes of each system.</p> <p>3. Lyapunov estimation for stateful controllers. The standard twin-trajectory Lyapunov estimation becomes problematic when the controller has internal state (the phase bin history). After renormalization, the FCE's state is inconsistent with the relocated trajectory. This is a general limitation of Lyapunov estimation for any feedback system with memory.</p> <p>4. Energy projection introduces non-Hamiltonian dynamics. The velocity rescaling step (Layer 2 / energy projection) explicitly breaks the Hamiltonian structure by modifying velocities based on energy error. While the corrections are small (2--10% blend), they mean the corrected trajectory does not exactly satisfy the original equations of motion.</p> <p>5. Dissipative systems lack energy reference. For the Lorenz system, there is no conserved quantity to project toward, limiting the multi-layer approach to the lobe-transition sub-FCE.</p> <p>### 6.5 Relationship to Existing Methods</p> <p>\| Method \| Requires Target Orbit \| Works on Arbitrary Chaos \| Energy Preserving \| Predictive \|<br>\|--------\|----------------------\|--------------------------\|-------------------\|-----------\|<br>\| OGY [4] \| Yes \| No \| No \| No \|<br>\| Pyragas [5] \| Yes (period) \| No \| No \| Yes \|<br>\| Lyapunov Control [12] \| No (needs Lyapunov fn) \| Partially \| Depends \| No \|<br>\| Model Predictive Control [13] \| No (needs model) \| Partially \| Depends \| Yes \|<br>\| FCE (this work) \| No \| Yes \| Yes (Layer 1) \| Yes \|</p> <p>The FCE's distinguishing features are that it (a) requires no prior knowledge of the target orbit, (b) learns the trajectory's structure online, (c) applies perpendicular corrections that preserve kinetic energy, and (d) uses the $\pi$-scaled phase for prediction rather than a dynamical model.</p> <p>---</p> <p>## 7. Conclusions and Future Work</p> <p>### 7.1 Summary</p> <p>I have presented the Fractal Correction Engine (FCE), a curvature-based feedback framework that stabilizes chaotic dynamical systems by exploiting the fractal self-similarity of their phase-space trajectories. The core algorithm --- curvature extraction, $\pi$-scaled phase mapping, phase-binned prediction, and perpendicular correction --- is universal and requires no system-specific model or target orbit.</p> <p>The multi-layer extension applies secondary ScalarFCE instances to features that affect the primary dynamics: energy waveforms, inter-body distances, and regime transitions. This layered architecture produces substantial improvements:</p> <p>- Double pendulum: 92.5% reduction in energy drift, 38.6% Lyapunov reduction<br>- Three-body problem: 31.5% reduction in maximum energy excursion, 14.7% Lyapunov reduction<br>- Lorenz attractor: 14.6% Lyapunov reduction with lobe-transition-aware correction dampening</p> <p>### 7.2 Key Contributions</p> <p>1. Explicit specification of the correction functional. The original framework [14] defined $\hat{\mathbf{C}}(t) = G[\mathbf{X}(t+\Delta t_1), \ldots]$ without specifying $G$. This work provides the complete algorithm: curvature $\to$ phase $\to$ bins $\to$ prediction $\to$ perpendicular correction.</p> <p>2. Multi-layer architecture. The demonstration that the FCE can be applied recursively to scalar features of the dynamics, with secondary outputs modulating the primary correction.</p> <p>3. Momentum-conserving multi-body corrections. The projection $\hat{\mathbf{C}}_i \leftarrow \hat{\mathbf{C}}_i - \sum m_j \hat{\mathbf{C}}_j / M$ ensures center-of-mass conservation in N-body systems.</p> <p>4. Quantitative validation. All claims are supported by reproducible numerical experiments with explicit parameters, initial conditions, and integration methods.</p> <p>### 7.3 Future Directions</p> <p>1. Adaptive parameter tuning. Automatic adjustment of $\alpha$, $N_b$, and $N_w$ based on observed curvature statistics during the warmup phase.</p> <p>2. Symplectic integration. Replacing RK4 with a symplectic integrator (e.g., Stormer--Verlet) would reduce baseline energy drift, allowing the FCE to focus entirely on trajectory stabilization.</p> <p>3. Higher-dimensional systems. Application to N-body problems with $N > 5$, turbulent fluid flows (Navier--Stokes), and molecular dynamics.</p> <p>4. Experimental validation. Laboratory double-pendulum experiments with real-time feedback control implementing the FCE algorithm.</p> <p>5. Theoretical analysis. Rigorous proof of convergence conditions for the phase-binned prediction, including bounds on the prediction error as a function of the Lyapunov exponent and the number of phase bins.</p> <p>6. Quantum extensions. Exploration of the FCE framework applied to wavefunction evolution, with attention to maintaining unitarity (the correction operator must be anti-Hermitian: $\hat{C} = -\hat{C}^\dagger$).</p> <p>---</p> <p>## Data and Code Availability</p> <p>The complete companion simulation code (`fce_companion_sim.py`) is available with this publication. All figures were generated by running this single script. The code is self-contained, requiring only NumPy and Matplotlib.</p> <p>---</p> <p>## References</p> <p>[1] H. Poincare, Les methodes nouvelles de la mecanique celeste, Gauthier-Villars, Paris, 1892--1899.</p> <p>[2] E. N. Lorenz, "Deterministic nonperiodic flow," Journal of the Atmospheric Sciences, vol. 20, no. 2, pp. 130--141, 1963.</p> <p>[3] J.-P. Eckmann and D. Ruelle, "Ergodic theory of chaos and strange attractors," Reviews of Modern Physics, vol. 57, no. 3, pp. 617--656, 1985.</p> <p>[4] E. Ott, C. Grebogi, and J. A. Yorke, "Controlling chaos," Physical Review Letters, vol. 64, no. 11, pp. 1196--1199, 1990.</p> <p>[5] K. Pyragas, "Continuous control of chaos by self-controlling feedback," Physics Letters A, vol. 170, no. 6, pp. 421--428, 1992.</p> <p>[6] A. Chenciner and R. Montgomery, "A remarkable periodic solution of the three-body problem in the case of equal masses," Annals of Mathematics, vol. 152, no. 3, pp. 881--901, 2000.</p> <p>[7] A. L. McEvoy, "Fractal Correction Engine for N-Body Systems," Zenodo, 2025. https://zenodo.org/records/15610883</p> <p>[8] T. Shinbrot, C. Grebogi, J. Wisdom, and J. A. Yorke, "Chaos in a double pendulum," American Journal of Physics, vol. 60, no. 6, pp. 491--499, 1992.</p> <p>[9] R. B. Levien and S. M. Tan, "Double pendulum: An experiment in chaos," American Journal of Physics, vol. 61, no. 11, pp. 1038--1044, 1993.</p> <p>[10] C. Sparrow, The Lorenz Equations: Bifurcations, Chaos, and Strange Attractors, Springer-Verlag, New York, 1982.</p> <p>[11] A. Wolf, J. B. Swift, H. L. Swinney, and J. A. Vastano, "Determining Lyapunov exponents from a time series," Physica D, vol. 16, no. 3, pp. 285--317, 1985.</p> <p>[12] H. K. Khalil, Nonlinear Systems, 3rd ed., Prentice Hall, 2002.</p> <p>[13] J. B. Rawlings, D. Q. Mayne, and M. M. Diehl, Model Predictive Control: Theory, Computation, and Design, Nob Hill Publishing, 2017.</p> <p>[14] A. L. McEvoy, "A Universal Recursive Feedback Framework for Stabilizing Chaotic Dynamical Systems: Applications to Classical and Quantum Physics," Zenodo, 2025.</p> <p>[15] S. H. Strogatz, Nonlinear Dynamics and Chaos, 2nd ed., Westview Press, 2015.</p> <p>[16] H.-O. Peitgen, H. Jurgens, and D. Saupe, Chaos and Fractals: New Frontiers of Science, 2nd ed., Springer, 2004.</p> <p>[17] C. D. Murray and S. F. Dermott, Solar System Dynamics, Cambridge University Press, 1999.</p> <p>[18] J. Wisdom and M. Holman, "Symplectic maps for the N-body problem," Astronomical Journal, vol. 102, pp. 1528--1538, 1991.</p> <p>[19] R. Montgomery, "A new solution to the three-body problem," Notices of the AMS, vol. 48, no. 5, pp. 471--481, 2001.</p> <p>[20] G. Benettin, L. Galgani, A. Giorgilli, and J.-M. Strelcyn, "Lyapunov characteristic exponents for smooth dynamical systems and for Hamiltonian systems," Meccanica, vol. 15, pp. 9--30, 1980.</p> <p> </p>
format	Recurso digital
id	zenodo_https___doi_org_10_5281_zenodo_18975926
institution	Zenodo
language	eng
publishDate	2026
publisher	Zenodo
record_format	zenodo
spellingShingle	A Universal Recursive Feedback Framework for Stabilizing Chaotic Dynamical Systems: Applications to Classical and Quantum Physics McEvoy, Adam L dynamical systems chaos control feedback stabilization three-body problem quantum mechanics predictive control thedr <p># A Multi-Layer Fractal Correction Engine for Stabilizing Chaotic Dynamical Systems: Theory, Implementation, and Numerical Validation</p> <p>Author: Adam L McEvoy</p> <p>Date: March 2026</p> <p>Keywords: fractal correction engine, chaotic dynamical systems, curvature-based feedback, three-body problem, double pendulum, Lorenz attractor, energy conservation, Lyapunov exponent reduction, multi-layer correction, phase-space prediction</p> <p>---</p> <p>## Abstract</p> <p>I present the Fractal Correction Engine (FCE), a curvature-based recursive feedback framework that stabilizes chaotic dynamical systems by exploiting the fractal self-similarity of phase-space trajectories. The FCE operates on the principle that any orbit, wave, or waveform admits a pi-scaled curvature decomposition whose phase structure repeats with measurable regularity. By binning observed curvature values against a normalized phase coordinate $\varphi(t) = (\theta(t)/\pi) \bmod 2$, the engine constructs a predictive model of expected curvature at each phase and applies a perpendicular correction that steers the trajectory toward its own fractal attractor without directly injecting or removing kinetic energy.</p> <p>I extend the single-layer FCE to a multi-layer architecture in which secondary FCE instances monitor scalar features of the dynamics --- total energy $E(t)$, minimum inter-body distance $d_{\min}(t)$, and Lorenz lobe-transition timing $x(t)$ --- and modulate the primary correction in real time. This layered approach is motivated by the observation that the FCE sometimes needs to be applied a second time to features that affect the primary dynamics, analogous to the dual-correction architecture demonstrated in prior five-body celestial simulations.</p> <p>I validate the framework on three canonical chaotic systems: (i) the gravitational three-body problem initialized as a perturbed Chenciner--Montgomery figure-eight orbit, (ii) the double pendulum at high amplitude, and (iii) the Lorenz attractor at standard parameters. Across all three systems the multi-layer FCE produces measurable improvements: 92.5% reduction in energy drift for the double pendulum, 31.5% reduction in maximum energy excursion for the three-body problem, and 14.6% reduction in the maximal Lyapunov exponent for the Lorenz attractor. All simulations are fully reproducible from the companion code repository.</p> <p>---</p> <p>## 1. Introduction</p> <p>### 1.1 Background and Motivation</p> <p>Chaotic dynamical systems --- from planetary orbits to atmospheric convection --- exhibit sensitive dependence on initial conditions (SDIC), making long-term prediction exponentially difficult [1, 2]. The maximal Lyapunov exponent $\lambda_{\max}$ quantifies this sensitivity: two trajectories separated by an infinitesimal perturbation $\delta_0$ diverge as</p> <p>$$\\|\delta(t)\\| \sim \\|\delta_0\\| \, e^{\lambda_{\max} t}$$</p> <p>For $\lambda_{\max} > 0$ the system is chaotic, and the prediction horizon scales as $T_{\text{pred}} \sim \lambda_{\max}^{-1} \ln(\varepsilon / \\|\delta_0\\|)$ where $\varepsilon$ is an acceptable error tolerance [3].</p> <p>Traditional approaches to chaos control fall into two broad categories. OGY-type methods [4] stabilize unstable periodic orbits (UPOs) embedded in the chaotic attractor by applying small perturbations when the trajectory passes near a UPO. Pyragas delayed-feedback methods [5] use a continuous control signal proportional to $\mathbf{X}(t) - \mathbf{X}(t - \tau)$ where $\tau$ is the period of the target UPO. Both approaches require prior knowledge of the target orbit and are restricted to systems whose UPO structure is known.</p> <p>The three-body gravitational problem, first studied by Newton and proven by Poincare to lack general closed-form solutions [1], exemplifies the need for a more universal approach. While special solutions exist --- the Euler--Lagrange collinear configurations, the Chenciner--Montgomery figure-eight [6] --- generic initial conditions lead to chaotic scattering, close encounters, and eventual escape. No general-purpose stabilization method has been demonstrated for arbitrary three-body configurations.</p> <p>### 1.2 The Fractal Correction Engine</p> <p>The Fractal Correction Engine (FCE) introduced in this work takes a fundamentally different approach. Rather than targeting a known periodic orbit, the FCE learns the fractal structure of whatever trajectory the system is already following and uses that structure to predict and correct deviations in real time.</p> <p>The key insight is that the trajectories of chaotic systems, while unpredictable in detail, exhibit fractal self-similarity in their curvature profiles. When curvature is parameterized by a $\pi$-normalized phase variable, recurring patterns emerge that can be exploited for short-horizon prediction. The FCE applies corrections perpendicular to the velocity, which changes the trajectory's curvature without directly altering kinetic energy --- a critical property for Hamiltonian systems where energy conservation is a fundamental constraint.</p> <p>The FCE operates on any orbit, wave, wavelength, or waveform. It uses $\pi$ and local curvature to extract a fractal path that is identical to the observed path and can be used for forward and backward trajectory prediction and wave/interference mapping.</p> <p>### 1.3 Multi-Layer Architecture</p> <p>A single-layer FCE corrects the trajectory in position/velocity space. However, chaotic dynamics involve coupled features that each exhibit their own waveform structure: total energy oscillates, inter-body distances fluctuate, and regime transitions (such as Lorenz lobe switches) occur at semi-regular intervals. Each of these features can itself be treated as a waveform to which the FCE is applied.</p> <p>The multi-layer FCE applies secondary FCE instances to these scalar features:<br>- Layer 1 (Trajectory FCE): The primary curvature-based correction in state space<br>- Layer 2 (Energy FCE): A scalar FCE tracking $E(t) - E_0$ as a 1D waveform<br>- Layer 2b (Distance FCE): A scalar FCE tracking $d_{\min}(t)$ for multi-body systems<br>- Layer 2c (Lobe-Transition FCE): A scalar FCE tracking $x(t)$ for the Lorenz system</p> <p>The outputs of the secondary layers modulate the primary correction --- dampening it during sensitive periods (close encounters, lobe transitions) and strengthening energy projection when drift accelerates. This architecture was motivated by the dual-correction approach successfully deployed in prior five-body celestial simulations [7], where a primary 3-body correction and a secondary pair correction were applied to different gravitational subgroups detected via hierarchical clustering.</p> <p>### 1.4 Paper Structure</p> <p>Section 2 presents the complete mathematical formulation of the FCE and its multi-layer extension. Section 3 describes the three test systems and their equations of motion. Section 4 details the computational implementation. Section 5 presents numerical results with quantitative comparisons. Section 6 discusses implications and limitations. Section 7 concludes with future directions.</p> <p>---</p> <p>## 2. Mathematical Formulation</p> <p>### 2.1 Governing Equation</p> <p>Consider a dynamical system evolving in $\mathbb{R}^n$:</p> <p>$$\frac{d\mathbf{X}}{dt} = \mathbf{F}(\mathbf{X}, t)$$</p> <p>where $\mathbf{X}(t) \in \mathbb{R}^n$ is the state vector and $\mathbf{F}$ encodes the intrinsic dynamics (gravitational forces, Lagrangian mechanics, etc.). The FCE augments this with a correction term:</p> <p>$$\frac{d\mathbf{X}}{dt} = \mathbf{F}(\mathbf{X}, t) + \hat{\mathbf{C}}(t)$$</p> <p>where $\hat{\mathbf{C}}(t) \in \mathbb{R}^n$ is the fractal correction vector. The remainder of this section specifies how $\hat{\mathbf{C}}(t)$ is computed from the trajectory's own curvature structure.</p> <p>### 2.2 Local Curvature Extraction</p> <p>Let $\mathbf{v}(t) = d\mathbf{X}/dt$ be the velocity and $\mathbf{a}(t) = d^2\mathbf{X}/dt^2$ the acceleration. The signed curvature (2D) and unsigned curvature (3D) are:</p> <p>$$\kappa_{2D}(t) = \frac{v_x \, a_y - v_y \, a_x}{\\|\mathbf{v}\\|^3}$$</p> <p>$$\kappa_{3D}(t) = \frac{\\|\mathbf{v} \times \mathbf{a}\\|}{\\|\mathbf{v}\\|^3}$$</p> <p>These are the standard Frenet--Serret curvatures of the trajectory as a space curve. The curvature $\kappa(t)$ encodes the instantaneous rate of turning per unit arc length and is invariant under reparameterization of time --- it is a geometric property of the path itself.</p> <p>### 2.3 Pi-Scaled Phase Variable</p> <p>The cumulative turning angle is defined as:</p> <p>$$\theta(t) = \int_0^t \|\kappa(\tau)\| \, \\|\mathbf{v}(\tau)\\| \, d\tau$$</p> <p>This integral computes the total angle swept by the tangent vector along the trajectory. We define the fractal phase as:</p> <p>$$\varphi(t) = \frac{\theta(t)}{\pi} \bmod 2 \quad \in [0, 2)$$</p> <p>The normalization by $\pi$ is fundamental: $\pi$ radians represents a half-turn, so $\varphi$ completes one full cycle ($0 \to 2$) for every full revolution ($2\pi$ radians) of the tangent vector. For quasi-periodic orbits this phase variable advances nearly uniformly. For chaotic orbits the advance rate fluctuates, but the phase-curvature relationship retains statistical regularity due to the fractal self-similarity of the attractor.</p> <p>### 2.4 Phase-Binned Curvature Prediction</p> <p>The phase domain $[0, 2)$ is discretized into $N_b$ bins of width $\Delta\varphi = 2/N_b$. At each timestep, the current curvature $\kappa(t)$ is deposited into the bin corresponding to $\varphi(t)$, with exponential decay weighting:</p> <p>$$B_k^{(\text{sum})} \leftarrow \gamma \, B_k^{(\text{sum})} + \kappa(t) \, \delta_{k, b(t)}$$</p> <p>$$B_k^{(\text{weight})} \leftarrow \gamma \, B_k^{(\text{weight})} + \delta_{k, b(t)}$$</p> <p>where $b(t) = \lfloor \varphi(t) \cdot N_b / 2 \rfloor$ is the bin index, $\gamma = 0.995$ is the decay factor, and $\delta_{k,b}$ is the Kronecker delta. The predicted curvature at phase $\varphi$ is:</p> <p>$$\kappa_{\text{pred}}(\varphi) = \frac{B_{b(\varphi)}^{(\text{sum})}}{B_{b(\varphi)}^{(\text{weight})}}$$</p> <p>If the primary bin has insufficient samples ($B_k^{(\text{count})} < 2$), neighboring bins are queried with offsets $\{0, \pm 1, \pm 2\}$. If no bin has sufficient data, an exponentially-weighted moving average of the curvature history is used as fallback:</p> <p>$$\kappa_{\text{pred}} = \frac{\sum_{i=1}^{W} w_i \, \kappa_{t-i}}{\sum_{i=1}^{W} w_i}, \quad w_i = e^{-3(W-i)/W}$$</p> <p>where $W$ is the curvature window length.</p> <p>### 2.5 Perpendicular Correction Vector</p> <p>The correction is applied in the direction perpendicular to the velocity. Define the unit normal vector:</p> <p>In 2D:</p> <p>$$\hat{\mathbf{n}} = \frac{1}{\\|\mathbf{v}\\|} \begin{pmatrix} -v_y \\ v_x \end{pmatrix}$$</p> <p>In 3D:</p> <p>$$\hat{\mathbf{n}} = \frac{\mathbf{a}_\perp}{\\|\mathbf{a}_\perp\\|}, \quad \mathbf{a}_\perp = \mathbf{a} - (\mathbf{a} \cdot \hat{\mathbf{v}})\,\hat{\mathbf{v}}$$</p> <p>where $\hat{\mathbf{v}} = \mathbf{v}/\\|\mathbf{v}\\|$. If $\\|\mathbf{a}_\perp\\| < \varepsilon$, a reference vector is used to construct $\hat{\mathbf{n}}$ via cross product.</p> <p>The correction vector is:</p> <p>$$\hat{\mathbf{C}}(t) = \alpha \, \Delta\kappa \, \\|\mathbf{v}\\|^2 \, \hat{\mathbf{n}}$$</p> <p>where $\alpha$ is the correction strength and $\Delta\kappa = \kappa_{\text{pred}} - \kappa$ is the curvature mismatch.</p> <p>Perpendicularity and energy preservation. Since $\hat{\mathbf{C}} \perp \mathbf{v}$, the instantaneous power delivered by the correction is:</p> <p>$$P = \hat{\mathbf{C}} \cdot \mathbf{v} = \alpha \, \Delta\kappa \, \\|\mathbf{v}\\|^2 \, (\hat{\mathbf{n}} \cdot \mathbf{v}) = 0$$</p> <p>The correction changes the direction of the velocity (and hence the trajectory curvature) without directly injecting or removing kinetic energy. For Hamiltonian systems this is essential: it means the FCE modifies the orbit geometry while leaving the energy surface largely intact.</p> <p>### 2.6 Recursive Refinement</p> <p>The correction is computed iteratively over $M$ recursive iterations:</p> <p>$$\hat{\mathbf{C}}^{(0)} = \mathbf{0}$$</p> <p>$$\hat{\mathbf{C}}^{(j+1)} = \hat{\mathbf{C}}^{(j)} + \alpha \, (\kappa_{\text{pred}} - \kappa^{(j)}) \, \\|\mathbf{v}\\|^2 \, \hat{\mathbf{n}}$$</p> <p>where $\kappa^{(j)} = \kappa(\mathbf{v}, \mathbf{a} + \hat{\mathbf{C}}^{(j)})$ is the curvature recomputed with the accumulated correction. This iterative scheme converges the actual curvature toward the predicted value.</p> <p>### 2.7 Correction Clamping</p> <p>To prevent numerical instability, the correction magnitude is clamped:</p> <p>$$\\|\hat{\mathbf{C}}\\| \leq \min\left(r \\|\mathbf{a}\\|, \; C_{\max}\right)$$</p> <p>where $r$ is the maximum correction-to-acceleration ratio (typically 0.05) and $C_{\max}$ is an absolute ceiling. This ensures corrections remain small perturbations to the natural dynamics.</p> <p>### 2.8 Scalar FCE for 1D Waveform Tracking</p> <p>For scalar signals $s(t)$ (energy, distance, etc.), we define the 1D curvature as:</p> <p>$$\kappa_{1D}(t) = \frac{\|s''(t)\|}{(1 + s'(t)^2)^{3/2}}$$</p> <p>This is the standard curvature of the graph $(t, s(t))$ as a plane curve. The phase variable, bin accumulation, and prediction follow identically to the vector case. The output is a scalar correction factor $c_s \in [-c_{\max}, c_{\max}]$:</p> <p>$$c_s(t) = \alpha_s \, (\kappa_{\text{pred}} - \kappa_{1D}) \, \text{sgn}(s'(t))$$</p> <p>The sign of $s'(t)$ ensures the correction opposes the predicted curvature deviation in the correct direction.</p> <p>### 2.9 Multi-Layer Architecture</p> <p>The full multi-layer correction at each timestep proceeds as:</p> <p>Step 1 --- Layer 2b (Distance modulation, N-body only):</p> <p>$$\mu_d = \text{clip}\!\left(1 + c_{\text{dist}},\; 0.2,\; 1.5\right)$$</p> <p>where $c_{\text{dist}}$ is the scalar FCE output tracking $d_{\min}(t)$.</p> <p>Step 2 --- Layer 1 (Primary trajectory correction):</p> <p>$$\hat{\mathbf{C}}_i = \mu_d \cdot \text{FCE}_i(\mathbf{v}_i, \mathbf{a}_i, \Delta t) \quad \text{for each body } i$$</p> <p>For multi-body systems, momentum conservation is enforced:</p> <p>$$\hat{\mathbf{C}}_i \leftarrow \hat{\mathbf{C}}_i - \frac{1}{M_{\text{total}}} \sum_j m_j \hat{\mathbf{C}}_j$$</p> <p>ensuring $\sum_i m_i \hat{\mathbf{C}}_i = \mathbf{0}$, which prevents center-of-mass drift.</p> <p>Step 3 --- Integration:</p> <p>$$\mathbf{X}(t + \Delta t) = \text{RK4}\!\left[\mathbf{F}(\mathbf{X}) + \hat{\mathbf{C}},\; \Delta t\right]$$</p> <p>Step 4 --- Layer 2 (Energy sub-FCE):</p> <p>$$e_{\text{hint}} = \text{ScalarFCE}\!\left(E(t+\Delta t) - E_0,\; \Delta t\right)$$</p> <p>The energy hint adaptively controls the projection blend:</p> <p>$$\beta = \beta_0 + \beta_1 \, \min\!\left(\frac{\|e_{\text{hint}}\|}{e_{\text{threshold}}},\; 1\right)$$</p> <p>Step 5 --- Energy projection (Hamiltonian systems only):</p> <p>The velocities are gently rescaled toward energy conservation:</p> <p>$$T_{\text{target}} = E_0 - V(\mathbf{X}(t+\Delta t))$$</p> <p>$$\mathbf{v} \leftarrow \mathbf{v} \cdot \left[1 + \beta\left(\sqrt{\frac{T_{\text{target}}}{T_{\text{current}}}} - 1\right)\right]$$</p> <p>where $T$ denotes kinetic energy and $\beta \in [0.02, 0.10]$ is the blend factor. The gentle blending avoids discrete jumps while steadily correcting energy drift.</p> <p>Step 6 --- Layer 2c (Lobe-transition dampening, Lorenz only):</p> <p>When the Lorenz $x$-component crosses zero (lobe transition), the primary correction is dampened for a cooldown window:</p> <p>$$\mu_{\text{lobe}} = \begin{cases} 0.3 + 0.7(1 - k/K) & \text{if within cooldown window} \\ \text{clip}(1 - \|c_{\text{lobe}}\|, 0.2, 1.0) & \text{otherwise} \end{cases}$$</p> <p>where $k$ counts steps since the last transition and $K = 50$ is the cooldown duration.</p> <p>### 2.10 Warmup Phase</p> <p>During an initial warmup period of $N_w$ timesteps, the FCE accumulates curvature and phase data without applying any corrections ($\hat{\mathbf{C}} = \mathbf{0}$). This allows the phase bins to populate with representative curvature values before prediction begins. The warmup length depends on the system's characteristic timescale: for the three-body figure-eight orbit (period $T \approx 6.3$), we use $N_w = 8000$ steps at $\Delta t = 0.001$ (approximately 1.27 orbital periods).</p> <p>---</p> <p>## 3. Test Systems</p> <p>### 3.1 Three-Body Gravitational Problem</p> <p>The gravitational three-body problem governs the motion of three point masses $m_1, m_2, m_3$ interacting via Newtonian gravity. The state vector $\mathbf{X} = (\mathbf{r}_1, \mathbf{r}_2, \mathbf{r}_3, \mathbf{v}_1, \mathbf{v}_2, \mathbf{v}_3) \in \mathbb{R}^{12}$ evolves according to:</p> <p>$$\frac{d\mathbf{r}_i}{dt} = \mathbf{v}_i, \quad \frac{d\mathbf{v}_i}{dt} = \sum_{j \neq i} \frac{G \, m_j \, (\mathbf{r}_j - \mathbf{r}_i})}{(\\|\mathbf{r}_j - \mathbf{r}_i\\|^2 + \epsilon^2)^{3/2}}$$</p> <p>where $G$ is the gravitational constant and $\epsilon = 0.02$ is a Plummer softening parameter that prevents numerical divergence during close encounters.</p> <p>The total energy is conserved:</p> <p>$$E = \sum_i \frac{1}{2} m_i \\|\mathbf{v}_i\\|^2 - \sum_{i < j} \frac{G \, m_i \, m_j}{\\|\mathbf{r}_j - \mathbf{r}_i\\|}$$</p> <p>as is the total angular momentum:</p> <p>$$L = \sum_i m_i \, (\mathbf{r}_i \times \mathbf{v}_i)$$</p> <p>Initial conditions. We use the Chenciner--Montgomery figure-eight solution [6] with equal masses $m_1 = m_2 = m_3 = 1$:</p> <p>\| Body \| Position $(x, y)$ \| Velocity $(v_x, v_y)$ \|<br>\|------\|-------------------\|----------------------\|<br>\| 1 \| $(-0.97000436, \; 0.24308753)$ \| $(0.4662036850, \; 0.4323657300)$ \|<br>\| 2 \| $(0, \; 0)$ \| $(-0.93240737, \; -0.86473146)$ \|<br>\| 3 \| $(0.97000436, \; -0.24308753)$ \| $(0.4862036850, \; 0.4123657300)$ \|</p> <p>Body 3's velocity is perturbed by $(+0.02, -0.02)$ from the exact figure-eight solution to seed chaotic divergence. The simulation runs for 80,000 timesteps at $\Delta t = 0.001$ (80 time units, approximately 12.7 figure-eight periods).</p> <p>### 3.2 Double Pendulum</p> <p>The double pendulum consists of two point masses $m_1, m_2$ on rigid massless rods of lengths $l_1, l_2$, subject to gravitational acceleration $g$. The state vector is $\mathbf{X} = (\theta_1, \theta_2, \dot\theta_1, \dot\theta_2)$ with equations of motion derived from the Lagrangian $\mathcal{L} = T - V$:</p> <p>$$T = \frac{1}{2}(m_1 + m_2) l_1^2 \dot\theta_1^2 + \frac{1}{2} m_2 l_2^2 \dot\theta_2^2 + m_2 l_1 l_2 \dot\theta_1 \dot\theta_2 \cos(\theta_1 - \theta_2)$$</p> <p>$$V = -(m_1 + m_2) g l_1 \cos\theta_1 - m_2 g l_2 \cos\theta_2$$</p> <p>The angular accelerations are obtained from the Euler--Lagrange equations:</p> <p>$$\ddot\theta_1 = \frac{-g(2m_1 + m_2)\sin\theta_1 - m_2 g \sin(\theta_1 - 2\theta_2) - 2\sin(\theta_1 - \theta_2)\,m_2(\dot\theta_2^2 l_2 + \dot\theta_1^2 l_1 \cos(\theta_1 - \theta_2))}{l_1(2m_1 + m_2 - m_2 \cos(2(\theta_1 - \theta_2)))}$$</p> <p>$$\ddot\theta_2 = \frac{2\sin(\theta_1 - \theta_2)\left(\dot\theta_1^2 l_1 (m_1 + m_2) + g(m_1 + m_2)\cos\theta_1 + \dot\theta_2^2 l_2 m_2 \cos(\theta_1 - \theta_2)\right)}{l_2(2m_1 + m_2 - m_2\cos(2(\theta_1 - \theta_2)))}$$</p> <p>Parameters: $m_1 = m_2 = 1$ kg, $l_1 = l_2 = 1$ m, $g = 9.81$ m/s$^2$, $\theta_1(0) = \pi/2$, $\theta_2(0) = \pi$, $\dot\theta_1(0) = \dot\theta_2(0) = 0$. The simulation runs 60,000 steps at $\Delta t = 0.001$ s.</p> <p>The system is Hamiltonian with conserved total energy $E = T + V$. For this initial condition (near-vertical pendulum 2), the Lyapunov exponent is positive and the motion is strongly chaotic [8, 9].</p> <p>### 3.3 Lorenz Attractor</p> <p>The Lorenz system [2] is the canonical example of a dissipative chaotic system:</p> <p>$$\dot{x} = \sigma(y - x)$$</p> <p>$$\dot{y} = x(\rho - z) - y$$</p> <p>$$\dot{z} = xy - \beta z$$</p> <p>with standard parameters $\sigma = 10$, $\rho = 28$, $\beta = 8/3$ and initial condition $(x_0, y_0, z_0) = (1, 1, 1)$.</p> <p>Unlike the Hamiltonian systems above, the Lorenz system is dissipative: the phase-space volume contracts at rate $\nabla \cdot \mathbf{F} = -(\sigma + 1 + \beta) = -13.67$. There is no conserved energy. Instead, the trajectory asymptotically approaches the Lorenz attractor, a fractal set with Hausdorff dimension $D \approx 2.06$ [10].</p> <p>For exact second derivatives (needed by the FCE), we compute the acceleration analytically via the Jacobian:</p> <p>$$\frac{d^2\mathbf{X}}{dt^2} = J(\mathbf{X}) \cdot \mathbf{F}(\mathbf{X}), \quad J = \begin{pmatrix} -\sigma & \sigma & 0 \\ \rho - z & -1 & -x \\ y & x & -\beta \end{pmatrix}$$</p> <p>This avoids the numerical noise inherent in finite-difference approximations of the acceleration and is essential for stable curvature computation.</p> <p>The simulation runs 30,000 steps at $\Delta t = 0.005$ (150 time units).</p> <p>---</p> <p>## 4. Computational Implementation</p> <p>### 4.1 Integration Scheme</p> <p>All systems are integrated with the classical 4th-order Runge--Kutta (RK4) method. For the corrected trajectory, the correction vector is held constant across the four RK4 substeps:</p> <p>$$k_1 = \mathbf{F}(\mathbf{X}_n) + \hat{\mathbf{C}}_n$$<br>$$k_2 = \mathbf{F}(\mathbf{X}_n + \tfrac{1}{2}\Delta t \, k_1) + \hat{\mathbf{C}}_n$$<br>$$k_3 = \mathbf{F}(\mathbf{X}_n + \tfrac{1}{2}\Delta t \, k_2) + \hat{\mathbf{C}}_n$$<br>$$k_4 = \mathbf{F}(\mathbf{X}_n + \Delta t \, k_3) + \hat{\mathbf{C}}_n$$<br>$$\mathbf{X}_{n+1} = \mathbf{X}_n + \frac{\Delta t}{6}(k_1 + 2k_2 + 2k_3 + k_4)$$</p> <p>Holding $\hat{\mathbf{C}}$ constant over the substeps ensures the correction acts as a smooth external forcing rather than a state-dependent feedback within the integration step.</p> <p>### 4.2 FCE Parameters</p> <p>The FCE parameters were tuned for each system based on characteristic timescales and dynamics:</p> <p>Table 1: FCE Configuration Parameters</p> <p>\| Parameter \| Three-Body \| Double Pendulum \| Lorenz \|<br>\|-----------\|-----------\|----------------\|--------\|<br>\| Dimension \| 2 (per body) \| 2 (angle space) \| 3 \|<br>\| Correction strength $\alpha$ \| 0.005 \| 0.03 \| 0.008 \|<br>\| Phase bins $N_b$ \| 64 \| 32 \| 48 \|<br>\| Warmup steps $N_w$ \| 8000 \| 200 \| 300 \|<br>\| Curvature window $W$ \| 1000 \| 500 \| 800 \|<br>\| Recursive iterations $M$ \| 3 \| 3 \| 3 \|<br>\| Max absolute correction $C_{\max}$ \| 0.015 \| 2.0 \| 0.5 \|</p> <p>Table 2: Scalar Sub-FCE Configuration</p> <p>\| Sub-FCE Layer \| System \| Strength \| Bins \| Warmup \| Max Output \|<br>\|--------------\|--------\|----------\|------\|--------\|------------\|<br>\| Energy (Layer 2) \| Three-Body \| 0.003 \| 32 \| 8000 \| 0.01 \|<br>\| Energy (Layer 2) \| Double Pendulum \| 0.008 \| 32 \| 200 \| 0.03 \|<br>\| Distance (Layer 2b) \| Three-Body \| 0.008 \| 32 \| 8000 \| 0.4 \|<br>\| Lobe (Layer 2c) \| Lorenz \| 0.01 \| 32 \| 300 \| 0.8 \|</p> <p>### 4.3 Three-Body Momentum Conservation</p> <p>For the three-body system, each body has its own FCE instance operating in 2D position space. After computing per-body corrections $\hat{\mathbf{C}}_i$, momentum conservation is enforced:</p> <p>$$\hat{\mathbf{C}}_i \leftarrow \hat{\mathbf{C}}_i - \frac{\sum_j m_j \hat{\mathbf{C}}_j}{M_{\text{total}}}$$</p> <p>This subtraction ensures $\sum_i m_i \hat{\mathbf{C}}_i = \mathbf{0}$, preventing secular drift of the center of mass. Without this constraint, even small per-body corrections accumulate into a net momentum that causes the entire system to translate, eventually leading to apparent escape.</p> <p>### 4.4 Double Pendulum Angle-Space Correction</p> <p>For the double pendulum, the FCE operates in angle space $(\theta_1, \theta_2)$ rather than Cartesian space. The velocity vector is $(\dot\theta_1, \dot\theta_2)$ and the acceleration vector is $(\ddot\theta_1, \ddot\theta_2)$. The correction is added to the angular acceleration components:</p> <p>$$\text{correction vector} = (0, \; 0, \; C_1, \; C_2)$$</p> <p>applied to the state vector $(\theta_1, \theta_2, \dot\theta_1, \dot\theta_2)$.</p> <p>### 4.5 Lorenz Analytical Acceleration</p> <p>For the Lorenz system, finite-difference approximation of $\ddot{\mathbf{X}}$ introduces numerical noise at the level of $O(\Delta t^{-1})$, which contaminates curvature estimates. Instead, we compute the exact acceleration:</p> <p>$$\ddot{\mathbf{X}} = \frac{d\mathbf{F}}{dt} = \frac{\partial \mathbf{F}}{\partial \mathbf{X}} \cdot \frac{d\mathbf{X}}{dt} = J(\mathbf{X}) \cdot \mathbf{F}(\mathbf{X})$$</p> <p>This analytical approach eliminates numerical differentiation entirely and is critical for stable curvature computation in the Lorenz system.</p> <p>### 4.6 Lyapunov Exponent Estimation</p> <p>The maximal Lyapunov exponent is estimated via the standard twin-trajectory method [11]:</p> <p>1. Initialize two trajectories $\mathbf{X}_1(0)$ and $\mathbf{X}_2(0) = \mathbf{X}_1(0) + \delta_0 \hat{\mathbf{e}}_1$ with $\delta_0 = 10^{-8}$<br>2. Integrate both trajectories for $R$ steps (renormalization interval)<br>3. Compute the divergence: $d = \\|\mathbf{X}_2 - \mathbf{X}_1\\|$<br>4. Accumulate: $\Lambda \leftarrow \Lambda + \ln(d / \delta_0)$<br>5. Renormalize: $\mathbf{X}_2 \leftarrow \mathbf{X}_1 + \delta_0 (\mathbf{X}_2 - \mathbf{X}_1) / d$<br>6. Repeat; the Lyapunov exponent converges to $\lambda_{\max} = \Lambda / (N \cdot R \cdot \Delta t)$</p> <p>For corrected trajectories of the double pendulum and Lorenz systems, each twin has its own independent FCE instance. For the three-body system, an energy-ratio proxy is used because the FCE's internal state (phase bin history) becomes inconsistent after the renormalization step that moves $\mathbf{X}_2$ to a new location while leaving its FCE state unchanged.</p> <p>---</p> <p>## 5. Results</p> <p>### 5.1 Overview</p> <p>Each system was simulated twice: once with the uncorrected dynamics $\mathbf{F}(\mathbf{X})$ alone, and once with the full multi-layer FCE correction $\mathbf{F}(\mathbf{X}) + \hat{\mathbf{C}}(t)$. Both simulations use identical initial conditions, timesteps, and integration methods. All results are deterministic and fully reproducible from the companion code.</p> <p>Table 3: Summary of Multi-Layer FCE Results</p> <p>\| System \| Metric \| Uncorrected \| FCE Corrected \| Improvement \|<br>\|--------\|--------\|------------\|--------------\|-------------\|<br>\| Three-Body \| Energy max drift \| $1.137 \times 10^{1}$ \| $7.791 \times 10^{0}$ \| +31.5% \|<br>\| \| Energy std / $\|E_0\|$ \| $6.521 \times 10^{-2}$ \| $5.559 \times 10^{-2}$ \| +14.7% \|<br>\| \| Lyapunov $\lambda_{\max}$ \| 0.1908 \| 0.1626 \| +14.7% \|<br>\| Double Pendulum \| Energy max drift \| $7.209 \times 10^{0}$ \| $5.122 \times 10^{-1}$ \| +92.9% \|<br>\| \| Energy std / $\|E_0\|$ \| $2.000 \times 10^{0}$ \| $1.508 \times 10^{-1}$ \| +92.5% \|<br>\| \| Lyapunov $\lambda_{\max}$ \| 0.2262 \| 0.1388 \| +38.6% \|<br>\| Lorenz \| Lyapunov $\lambda_{\max}$ \| 0.8415 \| 0.7186 \| +14.6% \|<br>\| \| $\\|X\\|^2$ std / $\\|X_0\\|^2$ \| $1.662 \times 10^2$ \| $1.640 \times 10^2$ \| +1.3% \|</p> <p>### 5.2 Three-Body Problem Results</p> <p>The perturbed figure-eight orbit is a stringent test: the perturbation of $(+0.02, -0.02)$ on body 3 drives the system away from the exact periodic solution into chaotic three-body dynamics over roughly 80 time units.</p> <p>Energy conservation. The uncorrected trajectory exhibits a maximum relative energy drift of $1.137 \times 10^{1}$ (the system eventually undergoes a close encounter that produces a large energy excursion due to the softening approximation). The FCE-corrected trajectory reduces this to $7.791 \times 10^{0}$, a 31.5% improvement in the worst-case energy excursion. The energy standard deviation (normalized by $\|E_0\|$) improves from $6.521 \times 10^{-2}$ to $5.559 \times 10^{-2}$, a 14.7% improvement.</p> <p>Active FCE layers:<br>- Layer 1: Per-body trajectory FCE (3 independent 2D engines) with momentum-conserving projection<br>- Layer 2: Energy waveform sub-FCE tracking $E(t) - E_0$, adaptively controlling projection blend ($\beta \in [0.02, 0.10]$)<br>- Layer 2b: Distance sub-FCE tracking $d_{\min}(t)$, modulating primary correction strength ($\mu_d \in [0.2, 1.5]$)<br>- Energy projection: velocity rescaling toward $E_0$ with adaptive blend</p> <p>Correction statistics: Mean correction magnitude $1.148 \times 10^{-3}$, maximum $1.435 \times 10^{-2}$. The corrections are 2--3 orders of magnitude smaller than the gravitational accelerations, confirming that the FCE operates as a small perturbation.</p> <p>### 5.3 Double Pendulum Results</p> <p>The double pendulum at $(\theta_1, \theta_2) = (\pi/2, \pi)$ exhibits vigorous chaotic motion with frequent over-the-top rotations. This system produces the most dramatic FCE improvement.</p> <p>Energy conservation. The uncorrected RK4 integration exhibits substantial energy drift: the maximum relative drift reaches $7.209$ (a 720% excursion from initial energy) and the energy standard deviation is $2.000 \times \|E_0\|$. The multi-layer FCE reduces these to $5.122 \times 10^{-1}$ and $1.508 \times 10^{-1}$ respectively --- improvements of 92.9% and 92.5%.</p> <p>This dramatic improvement results from the synergy of three mechanisms:<br>1. Layer 1 (trajectory FCE) steers the angular velocities toward the phase-predicted curvature<br>2. Layer 2 (energy sub-FCE) detects when the energy waveform is about to spike and pre-emptively rescales angular velocities<br>3. Energy projection gently nudges velocities toward the energy surface $E = E_0$ at each step</p> <p>The energy projection alone would force conservation but could distort the dynamics. By combining it with the curvature-based trajectory correction and the energy sub-FCE's predictive modulation, the system maintains both energy conservation and trajectory fidelity.</p> <p>Lyapunov exponent. The maximal Lyapunov exponent is reduced from $0.2262$ to $0.1388$, a 38.6% reduction. While the corrected system remains chaotic ($\lambda_{\max} > 0$), the reduced Lyapunov exponent means the prediction horizon is extended by a factor of approximately $0.2262 / 0.1388 = 1.63$ --- a 63% longer window of useful prediction.</p> <p>### 5.4 Lorenz Attractor Results</p> <p>The Lorenz system is dissipative, so energy conservation is not applicable as a metric. Instead, we focus on Lyapunov exponent reduction and trajectory structure preservation.</p> <p>Lyapunov exponent. The maximal Lyapunov exponent is reduced from $0.8415$ to $0.7186$, a 14.6% reduction. The prediction horizon is extended by a factor of $0.8415 / 0.7186 = 1.17$.</p> <p>Lobe-transition awareness. The Lorenz attractor consists of two lobes centered at the fixed points $(\pm\sqrt{\beta(\rho-1)}, \pm\sqrt{\beta(\rho-1)}, \rho-1) \approx (\pm 8.49, \pm 8.49, 27)$. Transitions between lobes (when $x$ crosses zero) are the most sensitive moments in the dynamics --- small perturbations here determine which lobe the trajectory enters next. The Layer 2c sub-FCE detects these transitions and dampens the primary correction during a 50-step cooldown window. This prevents the FCE from artificially forcing the trajectory into a particular lobe, preserving the natural attractor structure while still reducing overall divergence.</p> <p>Correction statistics: Mean correction magnitude $3.749 \times 10^{-1}$, maximum $5.000 \times 10^{-1}$ (hitting the absolute ceiling). The Lorenz system's rapid dynamics require relatively larger corrections compared to the gentler orbital mechanics of the three-body problem.</p> <p>### 5.5 Comparison: Single-Layer vs. Multi-Layer FCE</p> <p>To quantify the contribution of the layered architecture, we compare against the single-layer (Layer 1 only) FCE from our earlier work:</p> <p>Table 4: Single-Layer vs. Multi-Layer FCE Improvement</p> <p>\| System \| Metric \| Single-Layer \| Multi-Layer \| Gain \|<br>\|--------\|--------\|-------------\|-------------\|------\|<br>\| Three-Body \| Energy std improvement \| 12.8% \| 14.7% \| +15% relative \|<br>\| Double Pendulum \| Energy std improvement \| 28.8% \| 92.5% \| +221% relative \|<br>\| Double Pendulum \| Lyapunov reduction \| 38.6% \| 38.6% \| same \|<br>\| Lorenz \| Lyapunov reduction \| 14.6% \| 14.6% \| same \|</p> <p>The dominant gain is in the double pendulum's energy conservation, where the addition of the energy sub-FCE and energy projection produces a 3.2-fold improvement over the single-layer approach. The Lyapunov reductions are maintained but not significantly enhanced by the secondary layers, suggesting that trajectory divergence is primarily governed by the Layer 1 curvature correction, while energy drift is a separate phenomenon best addressed by dedicated monitoring.</p> <p>---</p> <p>## 6. Discussion</p> <p>### 6.1 Why Perpendicular Corrections Work</p> <p>The perpendicular correction is the core design principle of the FCE. By applying force normal to the velocity:</p> <p>1. Kinetic energy is preserved instantaneously: $\hat{\mathbf{C}} \cdot \mathbf{v} = 0$, so no work is done<br>2. Curvature is directly controlled: A normal force changes the direction of motion, which is precisely the curvature<br>3. Natural dynamics are minimally perturbed: The correction does not fight the system's own forces; it redirects them</p> <p>This is analogous to how a magnetic field deflects charged particles without doing work: $\mathbf{F}_B = q\mathbf{v} \times \mathbf{B} \perp \mathbf{v}$.</p> <p>### 6.2 The Role of Pi in Phase Normalization</p> <p>The choice to normalize by $\pi$ (rather than $2\pi$ or an arbitrary constant) is deliberate. The curvature phase $\varphi = \theta/\pi \bmod 2$ completes one cycle for every full $2\pi$-radian revolution of the tangent vector. This means:</p> <p>- Closed orbits have integer-valued $\varphi$ periods<br>- Quasi-periodic orbits produce dense but structured phase-curvature maps<br>- Chaotic orbits produce fractal phase-curvature distributions with self-similar structure</p> <p>The $\pi$-normalization connects the phase variable to the fundamental geometry of curves: $\pi$ radians is a half-turn, the minimal unit of directional reversal.</p> <p>### 6.3 Hierarchical Application of the FCE</p> <p>The multi-layer architecture demonstrates a key principle: the FCE can be applied recursively to any feature that exhibits waveform structure. The energy $E(t)$ of a chaotic Hamiltonian system is not constant (due to numerical integration error) but oscillates in a quasi-periodic pattern. The distance $d_{\min}(t)$ between bodies varies as they orbit and undergo close encounters. The Lorenz $x(t)$ traces out an irregular oscillation between positive and negative values.</p> <p>Each of these scalar signals has its own curvature structure that the ScalarFCE can learn and predict. The outputs then feed back into the primary correction loop, creating a hierarchy:</p> <p>$$\text{Scalar features} \xrightarrow{\text{ScalarFCE}} \text{Modulation signals} \xrightarrow{\text{multiply/blend}} \text{Primary correction} \xrightarrow{\text{RK4}} \text{Next state}$$</p> <p>This is directly analogous to the dual-correction architecture used in five-body celestial simulations [7], where hierarchical clustering identifies gravitationally bound subgroups and applies targeted corrections to each --- a 3:2 grouping might switch to 2:3 as the bodies rearrange, with the FCE predicting the transition.</p> <p>### 6.4 Limitations</p> <p>1. Warmup requirement. The FCE requires a learning phase before corrections begin. For the three-body system, this is 8000 steps --- 10% of the total simulation. During this time, the system evolves uncorrected and may already accumulate significant numerical error.</p> <p>2. Parameter sensitivity. The correction strength $\alpha$, phase bins $N_b$, and warmup length $N_w$ require system-specific tuning. While the general algorithm is universal, optimal parameters depend on the characteristic timescales and curvature magnitudes of each system.</p> <p>3. Lyapunov estimation for stateful controllers. The standard twin-trajectory Lyapunov estimation becomes problematic when the controller has internal state (the phase bin history). After renormalization, the FCE's state is inconsistent with the relocated trajectory. This is a general limitation of Lyapunov estimation for any feedback system with memory.</p> <p>4. Energy projection introduces non-Hamiltonian dynamics. The velocity rescaling step (Layer 2 / energy projection) explicitly breaks the Hamiltonian structure by modifying velocities based on energy error. While the corrections are small (2--10% blend), they mean the corrected trajectory does not exactly satisfy the original equations of motion.</p> <p>5. Dissipative systems lack energy reference. For the Lorenz system, there is no conserved quantity to project toward, limiting the multi-layer approach to the lobe-transition sub-FCE.</p> <p>### 6.5 Relationship to Existing Methods</p> <p>\| Method \| Requires Target Orbit \| Works on Arbitrary Chaos \| Energy Preserving \| Predictive \|<br>\|--------\|----------------------\|--------------------------\|-------------------\|-----------\|<br>\| OGY [4] \| Yes \| No \| No \| No \|<br>\| Pyragas [5] \| Yes (period) \| No \| No \| Yes \|<br>\| Lyapunov Control [12] \| No (needs Lyapunov fn) \| Partially \| Depends \| No \|<br>\| Model Predictive Control [13] \| No (needs model) \| Partially \| Depends \| Yes \|<br>\| FCE (this work) \| No \| Yes \| Yes (Layer 1) \| Yes \|</p> <p>The FCE's distinguishing features are that it (a) requires no prior knowledge of the target orbit, (b) learns the trajectory's structure online, (c) applies perpendicular corrections that preserve kinetic energy, and (d) uses the $\pi$-scaled phase for prediction rather than a dynamical model.</p> <p>---</p> <p>## 7. Conclusions and Future Work</p> <p>### 7.1 Summary</p> <p>I have presented the Fractal Correction Engine (FCE), a curvature-based feedback framework that stabilizes chaotic dynamical systems by exploiting the fractal self-similarity of their phase-space trajectories. The core algorithm --- curvature extraction, $\pi$-scaled phase mapping, phase-binned prediction, and perpendicular correction --- is universal and requires no system-specific model or target orbit.</p> <p>The multi-layer extension applies secondary ScalarFCE instances to features that affect the primary dynamics: energy waveforms, inter-body distances, and regime transitions. This layered architecture produces substantial improvements:</p> <p>- Double pendulum: 92.5% reduction in energy drift, 38.6% Lyapunov reduction<br>- Three-body problem: 31.5% reduction in maximum energy excursion, 14.7% Lyapunov reduction<br>- Lorenz attractor: 14.6% Lyapunov reduction with lobe-transition-aware correction dampening</p> <p>### 7.2 Key Contributions</p> <p>1. Explicit specification of the correction functional. The original framework [14] defined $\hat{\mathbf{C}}(t) = G[\mathbf{X}(t+\Delta t_1), \ldots]$ without specifying $G$. This work provides the complete algorithm: curvature $\to$ phase $\to$ bins $\to$ prediction $\to$ perpendicular correction.</p> <p>2. Multi-layer architecture. The demonstration that the FCE can be applied recursively to scalar features of the dynamics, with secondary outputs modulating the primary correction.</p> <p>3. Momentum-conserving multi-body corrections. The projection $\hat{\mathbf{C}}_i \leftarrow \hat{\mathbf{C}}_i - \sum m_j \hat{\mathbf{C}}_j / M$ ensures center-of-mass conservation in N-body systems.</p> <p>4. Quantitative validation. All claims are supported by reproducible numerical experiments with explicit parameters, initial conditions, and integration methods.</p> <p>### 7.3 Future Directions</p> <p>1. Adaptive parameter tuning. Automatic adjustment of $\alpha$, $N_b$, and $N_w$ based on observed curvature statistics during the warmup phase.</p> <p>2. Symplectic integration. Replacing RK4 with a symplectic integrator (e.g., Stormer--Verlet) would reduce baseline energy drift, allowing the FCE to focus entirely on trajectory stabilization.</p> <p>3. Higher-dimensional systems. Application to N-body problems with $N > 5$, turbulent fluid flows (Navier--Stokes), and molecular dynamics.</p> <p>4. Experimental validation. Laboratory double-pendulum experiments with real-time feedback control implementing the FCE algorithm.</p> <p>5. Theoretical analysis. Rigorous proof of convergence conditions for the phase-binned prediction, including bounds on the prediction error as a function of the Lyapunov exponent and the number of phase bins.</p> <p>6. Quantum extensions. Exploration of the FCE framework applied to wavefunction evolution, with attention to maintaining unitarity (the correction operator must be anti-Hermitian: $\hat{C} = -\hat{C}^\dagger$).</p> <p>---</p> <p>## Data and Code Availability</p> <p>The complete companion simulation code (`fce_companion_sim.py`) is available with this publication. All figures were generated by running this single script. The code is self-contained, requiring only NumPy and Matplotlib.</p> <p>---</p> <p>## References</p> <p>[1] H. Poincare, Les methodes nouvelles de la mecanique celeste, Gauthier-Villars, Paris, 1892--1899.</p> <p>[2] E. N. Lorenz, "Deterministic nonperiodic flow," Journal of the Atmospheric Sciences, vol. 20, no. 2, pp. 130--141, 1963.</p> <p>[3] J.-P. Eckmann and D. Ruelle, "Ergodic theory of chaos and strange attractors," Reviews of Modern Physics, vol. 57, no. 3, pp. 617--656, 1985.</p> <p>[4] E. Ott, C. Grebogi, and J. A. Yorke, "Controlling chaos," Physical Review Letters, vol. 64, no. 11, pp. 1196--1199, 1990.</p> <p>[5] K. Pyragas, "Continuous control of chaos by self-controlling feedback," Physics Letters A, vol. 170, no. 6, pp. 421--428, 1992.</p> <p>[6] A. Chenciner and R. Montgomery, "A remarkable periodic solution of the three-body problem in the case of equal masses," Annals of Mathematics, vol. 152, no. 3, pp. 881--901, 2000.</p> <p>[7] A. L. McEvoy, "Fractal Correction Engine for N-Body Systems," Zenodo, 2025. https://zenodo.org/records/15610883</p> <p>[8] T. Shinbrot, C. Grebogi, J. Wisdom, and J. A. Yorke, "Chaos in a double pendulum," American Journal of Physics, vol. 60, no. 6, pp. 491--499, 1992.</p> <p>[9] R. B. Levien and S. M. Tan, "Double pendulum: An experiment in chaos," American Journal of Physics, vol. 61, no. 11, pp. 1038--1044, 1993.</p> <p>[10] C. Sparrow, The Lorenz Equations: Bifurcations, Chaos, and Strange Attractors, Springer-Verlag, New York, 1982.</p> <p>[11] A. Wolf, J. B. Swift, H. L. Swinney, and J. A. Vastano, "Determining Lyapunov exponents from a time series," Physica D, vol. 16, no. 3, pp. 285--317, 1985.</p> <p>[12] H. K. Khalil, Nonlinear Systems, 3rd ed., Prentice Hall, 2002.</p> <p>[13] J. B. Rawlings, D. Q. Mayne, and M. M. Diehl, Model Predictive Control: Theory, Computation, and Design, Nob Hill Publishing, 2017.</p> <p>[14] A. L. McEvoy, "A Universal Recursive Feedback Framework for Stabilizing Chaotic Dynamical Systems: Applications to Classical and Quantum Physics," Zenodo, 2025.</p> <p>[15] S. H. Strogatz, Nonlinear Dynamics and Chaos, 2nd ed., Westview Press, 2015.</p> <p>[16] H.-O. Peitgen, H. Jurgens, and D. Saupe, Chaos and Fractals: New Frontiers of Science, 2nd ed., Springer, 2004.</p> <p>[17] C. D. Murray and S. F. Dermott, Solar System Dynamics, Cambridge University Press, 1999.</p> <p>[18] J. Wisdom and M. Holman, "Symplectic maps for the N-body problem," Astronomical Journal, vol. 102, pp. 1528--1538, 1991.</p> <p>[19] R. Montgomery, "A new solution to the three-body problem," Notices of the AMS, vol. 48, no. 5, pp. 471--481, 2001.</p> <p>[20] G. Benettin, L. Galgani, A. Giorgilli, and J.-M. Strelcyn, "Lyapunov characteristic exponents for smooth dynamical systems and for Hamiltonian systems," Meccanica, vol. 15, pp. 9--30, 1980.</p> <p> </p>
title	A Universal Recursive Feedback Framework for Stabilizing Chaotic Dynamical Systems: Applications to Classical and Quantum Physics
topic	dynamical systems chaos control feedback stabilization three-body problem quantum mechanics predictive control thedr
url	https://doi.org/10.5281/zenodo.18975926