Skip to main content
Sandbox Physics
Kozai–Lidov cycles from torque to planetary orbits

AC09 · Tutorial

Measure a secular period from frozen Cartesian trajectories

Sandbox Physics ·

The downloadable experiment measures eccentricity cycles from actual saved Cartesian integrations. A Python standard-library script verifies each input hash, extracts orbital elements, integrates the averaged equations and measures peaks. The files are synthetic dynamical benchmarks with fully specified units and initial conditions.

Download one self-contained experiment

Use the public materials link below and extract the archive into an empty directory. Keep the analysis script, eight CSV files and provenance JSON together. The archive also contains the reference-generation script, this README, fixed analysis outputs and checksums. None requires access to the private development repository.

Run the command below from that directory. The program uses only the Python standard library and writes a new analysis folder. It first verifies the integrator identifier and every raw-file SHA-256 against the provenance record. A failed hash stops the calculation instead of silently analyzing a different experiment.

Compare the newly written period-analysis JSON with the supplied expected file. Numerical entries should reproduce under the stated Python floating-point arithmetic. Platform details and hashes belong to the input provenance, while derived results retain their processing method.

python3 analyze-frozen-orbits.py --data-dir . --output-dir rerun-analysis

Reading: REBOUND 5.2.2 official Kozai example and IAS15 integrator · Tremaine & Yavetz, Why do Earth satellites stay up?, arXiv:1309.5244v1

Know exactly what each row contains

Each CSV row holds time, three position components, three velocity components and the companion’s three relative position components. All are measured in the nonrotating primary-centered frame. Initially the inner periapsis and companion lie on the positive first axis; the outer plane is tilted about that axis.

The gravitational constant, primary mass, companion mass and initial inner semimajor axis are one. The inner Kepler period is two pi. Divide time by the cube of the outer-to-inner semimajor-axis ratio to obtain the dimensionless secular time. The saved cadence is 0.005 in that time.

The provenance distinguishes full Newtonian and instantaneous-quadrupole files. The first group was generated by a three-particle IAS15 integration. The two quadrupole controls integrate the central attraction plus the explicitly stated tidal acceleration with a prescribed circular companion clock. Both groups save actual integrated positions and velocities.

The analysis forms angular momentum by a cross product and eccentricity from the instantaneous state. It independently integrates the averaged vector equations at a maximum step of 0.001 secular times. No fitted eccentricity waveform enters either history.

τ=tA3,A=a2a1,e=v×(r×v)−rr\tau=\frac{t}{A^3},\qquad A=\frac{a_2}{a_1},\qquad \boldsymbol e=\boldsymbol v\times(\boldsymbol r\times\boldsymbol v)-\frac{\boldsymbol r}{r}
The final identity uses the unit primary gravitational parameter.

Reading: REBOUND 5.2.2 official Kozai example and IAS15 integrator · Tremaine & Yavetz, Why do Earth satellites stay up?, arXiv:1309.5244v1

Choose a period estimator that tolerates fast oscillations

The raw direct eccentricity contains orbital-timescale structure. A local maximum finder would count many of those wiggles as secular cycles. The script first takes a centered moving average over one outer period, then selects high maxima and groups candidates separated by less than one secular time. The largest maximum in each group represents that excursion.

Only curves whose smoothed eccentricity range is at least 0.2 enter this large-excursion estimator. The low-tilt branch therefore reports no large cycle period. A stopped near-polar record has no completed recurrence and also reports no period. These outcomes are explicit JSON nulls, with explanations.

The reported period is the mean spacing of the selected peaks. Repeating the measurement with half and twice the smoothing window gives a processing sensitivity range. It is not a statistical confidence interval: the trajectory has no random measurement noise, and the alternative windows probe a deliberate analysis choice.

The method is designed for these records. For a chaotic or slowly drifting series, report individual peak spacings and inspect the raw curve. A changing spacing can itself carry the physical information.

P^=1N−1∑k=1N−1(τk+1−τk)\widehat P=\frac{1}{N-1}\sum_{k=1}^{N-1}(\tau_{k+1}-\tau_k)
At least two identified peaks are needed. The script preserves their times alongside the average period.

Reading: REBOUND 5.2.2 official Kozai example and IAS15 integrator · Klein & Hamilton, Short-range forces beyond double averaging, arXiv:2606.19456v1

Reproduce the three-model comparison

For the default full-force orbit, the grouped peak times are 2.510 and 7.565 secular times, yielding 5.055. The independently averaged record yields 5.025. Its direct maximum eccentricity is about 0.843487, while the averaged maximum is 0.838909. These differences coexist with much smaller integrator-check residuals over the separately tested first three secular times.

Now use the matched seventy-five-degree cases. At an outer-to-inner ratio of six, read full Newtonian, instantaneous quadrupole and double-average periods of 5.005, 4.810 and 4.670. At a ratio of twenty, read 4.675, 4.665 and 4.670. The controlled-comparisons block in the output keeps those three models together.

Comparing instantaneous quadrupole with its double average isolates removal of fast phases at the same spatial order. Comparing full Newtonian with instantaneous quadrupole tests the spatial truncation. Long nonlinear phase differences and finite-record peak selection can mix their effects; the two residuals need not supply additive corrections.

Open the derived CSV files to plot the raw eccentricity, smoothed eccentricity, averaged eccentricity and common-time residual. The near-polar derived file ends at the shared surviving interval, while the summary retains each model’s own completed time.

Reading: Liu, Muñoz & Lai, Short-range forces in triple systems, arXiv:1409.6717v2 · REBOUND 5.2.2 official Kozai example and IAS15 integrator · Klein & Hamilton, Short-range forces beyond double averaging, arXiv:2606.19456v1

Change one analysis choice, then regenerate only if needed

A useful first change is to alter the smoothing window without touching the force model. Record which peaks moved and whether the period conclusion survives. Next compare a phase-sensitive position error with the slower eccentricity residual. This separates an analysis decision from an integration decision.

To regenerate the raw files, install the pinned REBOUND 5.2.2 release in a separate environment and run the optional generator below. Its provenance records the source commit, package license, platform, exact initial conditions, accuracy parameter and stopping rule. A different platform may produce different last-bit CSV hashes. Save its fresh hashes in the newly generated provenance.

The optional relativistic example in the output is averaged only, with strength 0.8 and maximum eccentricity about 0.629278. The script explicitly records that no direct relativistic integration was performed. Turning that into a direct relativistic validation would require adding and independently checking a consistent force model.

python3 generate-ias15-reference.py --output-dir regenerated-inputs
python3 analyze-frozen-orbits.py --data-dir regenerated-inputs --output-dir regenerated-analysis

Reading: REBOUND 5.2.2 official Kozai example and IAS15 integrator · Liu, Muñoz & Lai, Short-range forces in triple systems, arXiv:1409.6717v2

Try the orbit exchange

  1. Export a browser record for one named preset and compare its provenance with the full-resolution CSV. Browser preset states are downsampled actual IAS15 samples.
  2. Compare the compact and wide cases in the data row, then use the public three-model output to identify which approximation each pair tests.

Exercises with explanations

  1. The compact instantaneous-quadrupole period shifts when the smoothing window doubles. Is that an uncertainty in a measured exoplanet period?

    Read the explanation

    No. These are deterministic synthetic integrations. The shift describes sensitivity of this peak estimator to a processing window. It should be reported alongside the raw history, without treating it as an observational confidence interval.

  2. How would you extend this tutorial to compare two integrator tolerances fairly?

    Read the explanation

    Preserve the initial state, force model, output times and stop rules. Regenerate with a tighter tolerance and its own provenance, then compare Cartesian states, orbital elements and the appropriate invariant over their common interval.

Sources and calculation scope

  1. Tremaine & Yavetz, Why do Earth satellites stay up?, arXiv:1309.5244v1

    Fixed 20 September 2013 version, DOI 10.1119/1.4874853. Potential 18, vector equations 23, constraints 24 and section IV give the quadrupole model and near-circular critical-angle limit.

  2. Liu, Muñoz & Lai, Short-range forces in triple systems, arXiv:1409.6717v2

    Fixed 7 January 2015 version, DOI 10.1093/mnras/stu2396. Sections 2–4, quadrupole potential 7, octupole 8, GR 32–33/42 and Figures 3–4/7. The browser retains the circular-outer, test-particle quadrupole; the guide reads the wider model.

  3. REBOUND 5.2.2 official Kozai example and IAS15 integrator

    The official example motivates an independent direct reference. Frozen trajectories use REBOUND 5.2.2, source commit c34134716cd484ed9cdb9f7d54ef8d859bdf547d, IAS15 accuracy parameter one hundred-billionth, with separately stated initial conditions and stopping rules. REBOUND software is GPL-3.0-only; no publisher code is copied into the Lab.

  4. Klein & Hamilton, Short-range forces beyond double averaging, arXiv:2606.19456v1

    Fixed 17 June 2026 preprint version; an August v2 also exists. Figure 1 separates averaging/force choices, Figures 2–3 show cycle changes, and Appendix B compares approximations with direct N-body examples. This is a reading guide; the browser does not reproduce its single-averaged GR/tidal calculations.

FeedbackSomething broken, confusing, or missing on this page?Tell us
What kind of feedback?

Sent with this page’s address, its Lab ID, the site language and your device type. Please leave out personal details. Feedback is kept for up to 12 months.Privacy policy