[Paper Review] The high-order Euler method and the spin-orbit model A fast algorithm for solving differential equations with small, smooth nonlinearity

Michele V. Bartuccelli, Jonathan H. B. Deane+1 more|arXiv (Cornell University)|Oct 22, 2014

Quantum chaos and dynamical systems23 references4 citations

TL;DR

This paper presents a high-order Euler method tailored for fast, accurate numerical integration of smooth, time-periodic, dissipative ODEs with small nonlinearity—particularly the spin-orbit problem. By exploiting perturbation theory and Fourier series structure in the solution, the algorithm achieves speedups of at least 7× over standard solvers while maintaining double-precision accuracy, with excellent agreement to analytical predictions for periodic orbits and capture probabilities.

ABSTRACT

We present an algorithm for the rapid numerical integration of smooth, time-periodic differential equations with small nonlinearity, particularly suited to problems with small dissipation. The emphasis is on speed without compromising accuracy and we envisage applications in problems where integration over long time scales is required; for instance, orbit probability estimation via Monte Carlo simulation. We demonstrate the effectiveness of our algorithm by applying it to the spin-orbit problem, for which we have derived analytical results for comparison with those that we obtain numerically. Among other tests, we carry out a careful comparison of our numerical results with the analytically predicted set of periodic orbits that exists for given parameters. Further tests concern the long-term behaviour of solutions moving towards the quasi-periodic attractor, and capture probabilities for the periodic attractors computed from the formula of Goldreich and Peale. We implement the algorithm in standard double precision arithmetic and show that this is adequate to obtain an excellent measure of agreement between analytical predictions and the proposed fast algorithm.

Motivation & Objective

• To develop a fast, accurate numerical integration algorithm for smooth, time-periodic, dissipative ODEs with small nonlinearity, especially for long-time simulations.
• To address the computational bottleneck in Monte Carlo simulations of orbit capture probabilities, where long transient times and large numbers of initial conditions are required.
• To achieve high speed without sacrificing accuracy by leveraging perturbation theory and Fourier structure in the solution.
• To validate the algorithm against analytical results for the spin-orbit model, particularly periodic orbits and capture probabilities.
• To demonstrate that standard double-precision arithmetic is sufficient for high-fidelity numerical results when combined with the proposed method.

Proposed method

• The method uses a high-order Taylor expansion of the solution, exploiting the fact that derivatives in the ODE contain explicit powers of the small parameter ε.
• It expresses the solution at a time step h in Fourier series form: x(t₀ + h) = x₀ + Aᵢ,₀(y₀) + Σ εʲ [Aᵢ,ⱼ cos(2jx₀) + Bᵢ,ⱼ sin(2jx₀)], where coefficients Aᵢ,ⱼ and Bᵢ,ⱼ are polynomials in y₀, h, and system parameters.
• The algorithm is derived using computer algebra to systematically compute the coefficients Aᵢ,ⱼ and Bᵢ,ⱼ up to a desired order, ensuring all terms have the correct ε-power scaling.
• It leverages the structure that each derivative in the Taylor series contains a factor of ε^d, where d is the degree of trigonometric terms, enabling efficient computation.
• Only the first few Fourier coefficients (j=1,2,3) are needed to meet error criteria, reducing computational cost significantly.
• The method is implemented in standard double precision, avoiding arbitrary-precision arithmetic while maintaining high accuracy.

Experimental results

Research questions

• RQ1 Can a high-order Euler method be constructed that is both fast and accurate for smooth, time-periodic, dissipative ODEs with small nonlinearity?
• RQ2 Does the proposed algorithm achieve significant speedups over standard solvers like Runge-Kutta or symplectic methods without losing accuracy?
• RQ3 To what extent does the numerical solution agree with analytical predictions for periodic orbits in the spin-orbit model?
• RQ4 Can the algorithm reliably compute long-term behavior and capture probabilities for periodic attractors in the spin-orbit system?
• RQ5 Is double-precision arithmetic sufficient for high-fidelity results when using this structured, high-order approach?

Key findings

• The algorithm achieves a speedup of at least 7× compared to traditional solvers like Runge-Kutta and symplectic methods such as Yoshida’s algorithm.
• Numerical results for the spin-orbit model show excellent agreement with analytical predictions: the computed frequency shift Δω ≈ 2.28×10⁻¹² matches the perturbation theory value of 2.284×10⁻¹².
• The method accurately captures the long-term behavior of solutions approaching the quasi-periodic attractor, confirming convergence to expected dynamics.
• Capture probabilities for periodic attractors computed via the algorithm align well with the Goldreich-Peale formula, validating its use in Monte Carlo simulations.
• The first three Fourier coefficients (j=1,2,3) are sufficient to meet the error criterion, indicating that higher-order terms decay rapidly in magnitude.
• Numerical evidence shows that the coefficients Aᵢ,ⱼ and Bᵢ,ⱼ decrease faster than εʲ for j=2,3, supporting the efficiency of the method’s truncation strategy.