{"id":"d71549c5-79ff-4d87-9571-7cd7cf746618","arxiv_id":"2506.13745","paper_version":2,"verdict":"CONDITIONAL","confidence":"MODERATE","novelty_score":7.0,"correctness_risk":"medium","formal_verification":"none","parameter_count":5,"one_line_summary":"A FFT-based second-order secular perturbation method for planetary three-body systems reproduces n-body secular frequencies with sub-percent errors even near mean-motion resonances.","lead":"This paper presents a numerical implementation of second-order canonical perturbation theory for two planets around a star, using FFT to eliminate fast orbital angles without series expansions. It validates the method on the Sun-Jupiter-Saturn system and three exoplanet systems, showing improved secular precession frequencies near mean-motion resonances.","discovery_kind":"new_method","skeptic_critique":{"model":"deepseek-v4-flash","headline":"GJ876/TIC accuracy rests on an unverified K'=20 truncation of the resonant second-order sum; no convergence check is provided for these high-eccentricity systems.","rationale":"The paper's core construction, Lie-transform averaging with FFT-evaluated Fourier coefficients, appears internally consistent, and the SJS/WASP-148 tests give real evidence that second-order terms improve secular frequencies. The load-bearing weak point is the unverified convergence of the resonant H2 truncation for exactly the systems used to claim strong-MMR applicability. The text itself flags slower Fourier convergence for TIC (Sect. V), and K'=20 is an instability-avoidance choice, not a convergence choice. Because Eq. (23b) contains small divisors, the neglected high-order harmonics are not guaranteed to be negligible even if their raw Fourier amplitudes decay; the product with 1/(k·n'_0) can enhance them. The quoted residuals (0.3-1% for TIC, 2.6% for GJ876) are of the same order as the possible truncation error, so the improvement could be cutoff-dependent. The reader identified this same assumption; I agree. A convergence test along the lines of Figs. 1-2 for GJ876/TIC would settle it. Since the concern is addressable and the method has independent support on SJS/WASP-148, the existing CONDITIONAL verdict is appropriate; no change is needed.","tokens_in":23075,"tokens_out":9718,"duration_ms":104443,"concrete_test":"Recompute the per-harmonic contributions to the second-order resonant coefficients h_l^2 (Eq. 23b) for GJ876 and TIC at their filtered initial conditions, as done for SJS in Figs. 1-2, and compare the cumulative sum at K'=20 with cumulative sums up to K=64. Then re-integrate the resonant secular model with K'=24, 32, 40, and 48, using constant divisors from Eq. (27) to avoid small-divisor blow-ups, and measure g1 and g2 by frequency analysis. If g2 shifts by more than ~0.5% for GJ876 (or g1 by more than ~0.5% for TIC), the reported improvement is cutoff-dependent and the strong-MMR validation is not established.","verdict_should_be":"UNCHANGED","load_bearing_attack":"The central claim that second-order secular models sharply improve precession frequencies near strong MMRs is supported for SJS and WASP-148, but for the two resonant exoplanet benchmarks (TIC, GJ876) it depends on the Fourier-truncation assumption stated in Sect. IV.a. The second-order resonant Hamiltonian, Eq. (23b), sums non-resonant harmonics k with small denominators k·n'_0. For TIC and GJ876 the sum is cut at K'=20 solely to avoid small-divisor instabilities (Sect. IV.a and Sect. V). Equation (29) justifies truncation by exponential decay of h_k, but the only convergence demonstration is Figs. 1-2, for SJS. For TIC the text says convergence is slower due to higher eccentricities and masses; for GJ876 no check is reported. At high eccentricity the decay rate is smaller, and small divisors amplify the tail, so the neglected |k|>20 terms can contribute at a level comparable to the quoted residuals (1% on g1 for TIC, 2.6% on g2 for GJ876). The claimed improvement could therefore be an artifact of the cutoff. The use of filtered n-body initial conditions and a modified TIC orbit weakens external-validity claims but is not by itself fatal; the unverified truncation is the load-bearing element.","agreement_with_reader":"agree"},"referee_report":{"model":"deepseek-v4-flash","summary":"The manuscript presents a numerical implementation of second-order canonical perturbation theory for the planetary three-body problem. Using Lie transforms, the authors construct secular Hamiltonians to second order in the planet-to-star mass ratio for both non-resonant and resonant configurations, including a mean-motion frequency correction. The Fourier coefficients of the disturbing function and their derivatives are computed on a grid via FFT, and the secular equations are integrated with a high-order Adams method. The method is benchmarked against direct n-body integrations for the Sun-Jupiter-Saturn system, WASP-148, TIC 279401253, and GJ 876. For SJS, second-order terms reduce precession frequency errors from 14--20% to below 1%, and the resonant second-order model reduces the g6 error to 0.05%. For WASP-148, errors drop below 0.1% in the non-resonant second-order model. For the resonant systems TIC and GJ 876, the resonant second-order model improves the g2 frequency substantially (for GJ 876 from a large first-order error of reported 112% down to 2.6%), although for TIC the g1 error slightly increases from 0.2% to 1%.","tokens_in":23355,"tokens_out":21946,"duration_ms":187756,"significance":"If the results hold, the paper offers a useful and comparatively general numerical route to second-order secular dynamics, avoiding expansions in eccentricities, inclinations, and semi-major axis ratios. The derivation is self-contained and the SJS and WASP-148 validations are convincing, with clear algorithmic detail and quantified improvements. The main significance is the potential applicability to high-eccentricity, resonant exoplanetary systems where classical analytical expansions are problematic. However, the support for the central claim in exactly those regimes is incomplete: for TIC and GJ 876 the resonant second-order sum is truncated at K'=20 without any convergence check, and the benchmark initial conditions are taken from low-pass-filtered n-body solutions of the same systems. These issues weaken, but do not by themselves invalidate, the otherwise sound SJS/WASP-148 results.","major_comments":[{"comment":"The truncation of the resonant second-order sum in Eq. (23b) to K'=20 for TIC 279401253 and GJ 876 is introduced solely to avoid small-divisor instabilities, and no convergence check is provided for these systems. The exponential-decay justification in Eq. (29) is demonstrated only for SJS (Figs. 1 and 2), and the text itself states that convergence is slower for TIC because of higher eccentricities and masses. Since the claimed improvements for these systems (e.g., the g2 error reduction from 112% to 2.6% for GJ 876) are comparable to the expected magnitude of the truncated terms, the reported accuracy could be an artifact of the cutoff. I request a convergence study over K' (for example, 12, 16, 20, 24, 32) for both systems, reporting the secular frequencies g1 and g2 and, if possible, the eccentricity curves. If larger K' produces small-divisor instabilities, then the method is not validated for these high-eccentricity resonant configurations and the conclusions should be restricted accordingly.","section":"§IV.a, §V (TIC/GJ 876)"},{"comment":"The initial conditions for the secular integrations are obtained by low-pass filtering the very same n-body solution against which the secular models are compared (Eq. (33) and the text of Sect. IV.e). This makes the benchmark a test of internal consistency between the secular equations and filtered initial conditions, not a test of the method's ability to evolve a system from published osculating elements. I recommend adding at least one case (for example, SJS or WASP-148) in which the secular model is initialized directly from osculating elements, with the first- or second-order Lie-transform corrections applied, and compared to an n-body integration from the same osculating elements. Such a test would substantially strengthen the practical applicability of the method.","section":"§IV.e and §V"},{"comment":"For TIC, the second-order resonant model improves g2 (from 5.1% to 0.3%) but degrades g1 (from 0.2% to 1%). The paper's explanation in terms of high masses, high eccentricities, and the limited validity of the mean-motion correction is plausible, but the net improvement is less clear than for SJS and WASP-148. I suggest an error-budget analysis for the TIC case, separating the contributions of the K' truncation, the mean-motion correction, and the initial-condition offset, and a more guarded statement of the method's performance for strongly resonant, high-eccentricity systems.","section":"§V (TIC) and Table I"}],"minor_comments":[{"comment":"The title and running header contain a typo: 'perturbatio n theory' with a space.","section":"Title/running header"},{"comment":"The text reports a 112% relative error for g2 at first order, but the values in Table I (n-body -0.0202, resonant first order -0.0025) give a relative error of about 88%. Please correct this inconsistency.","section":"§V (GJ 876) and Table I"},{"comment":"The statement in the abstract that the method 'avoids the need for expansions in orbital elements' is too strong; the method still requires a finite truncation of the Fourier series (K and K') and relies on the exponential decay of the coefficients. Please qualify the claim.","section":"Abstract"},{"comment":"The statement in Eq. (29) that the neglected remainder of the Fourier series is 'of order ε²' should be clarified; the truncation error of H1 depends on the analyticity domain and is not automatically O(ε²), especially at high eccentricity.","section":"Eq. (29)"},{"comment":"The paper mentions an implementation in C but does not state whether the code will be made publicly available. Given the numerical nature of the work, a code availability statement would be useful for reproducibility.","section":"§VI"}],"recommendation":"major_revision","confidential_remarks":"The paper is a solid contribution from a well-known group, and the SJS and WASP-148 results are convincing. The main technical concern is the unverified K'=20 truncation for the resonant systems TIC and GJ 876, which is load-bearing for the claim that the method improves accuracy in strong MMR regimes. The benchmark's use of filtered n-body initial conditions is a weaker but not fatal limitation. The 112% vs 88% discrepancy for GJ 876 appears to be a simple numerical typo. I would encourage the editor to request a revision focused on the convergence analysis for K' and a test starting from osculating elements."},"author_rebuttal":null,"desk_editor":{"model":"deepseek-v4-flash","letter":"The paper deserves a serious referee. The genuinely new piece is the FFT-based numerical evaluation of the second-order Lie-transform secular Hamiltonian for the planetary 3-body problem, including a resonant variant, and it works. On Sun-Jupiter-Saturn the convergence test is explicit, and second-order plus resonant terms bring the g5/g6 frequency errors from 14-20% down to 0.05-1%. On WASP-148 the second-order model drops g1/g2 errors to below 0.1%. That is a real, reproducible improvement over first-order Gauss/Schubart-type averaging, and the derivation of Eqs. (17), (20), (23) is compact and internally consistent.\n\nThe soft spots are real but not fatal. The stress-test concern is on target: for TIC 279401253 and GJ 876, the second-order resonant sum is cut at K'=20 solely to avoid small-divisor instabilities, and no convergence check is shown for these high-eccentricity, high-mass systems. At high e the Fourier coefficients decay more slowly and the small denominators amplify the tail, so the quoted residuals (1% on g1 for TIC, 2.6% on g2 for GJ876) could partly be an artifact of the cutoff. This is an empirical gap, not a demonstrated error--the SJS convergence data supports the exponential-decay premise, and the authors are candid about slower convergence for TIC. But the headline claim, improved accuracy especially near strong MMRs, needs that check.\n\nThe benchmark design is also partially circular: initial conditions for the secular models come from low-pass filtering the same n-body run used as ground truth. The authors acknowledge this and justify it via Eq. (33), so it is not hidden, but it means the validations measure self-consistency rather than predictive power from raw orbital elements. A test starting from published astrometric or RV elements without filtering would strengthen the practical claim. The TIC test also uses a modified eccentricity pushed into the stable zone, again acknowledged, but it weakens external validity.\n\nThe citation pattern looks fine: Laskar 1985, Brouwer-van Woerkom, Locatelli-Giorgilli, and the Gauss/Schubart averaging literature are all there. I did not see a missing key reference.\n\nWho is this for: anyone working on secular dynamics of exoplanetary systems or on practical perturbation theory. It should go to peer review. The central derivation seems sound, the SJS and WASP-148 validations are convincing, and the method is a genuine step beyond first-order numerical averaging. The main required revision is a convergence check or quantitative bound for K' on the resonant benchmarks, and ideally one test from unfiltered initial conditions. My own verdict is conditional, but I would send it to a serious referee.","headline":"Second-order secular theory via FFT is genuinely new and the SJS/WASP-148 tests are convincing, but the resonant exoplanet benchmarks rest on an unverified truncation that needs a convergence check before the headline claim is fully secure.","tokens_in":23924,"tokens_out":2626,"would_cite":true,"duration_ms":24486,"reading_group":"yes","serious_thinker":"yes","would_accept_peer_review":true},"rs_alignment":null,"lean_confirmation":null,"pith_extraction":{"msc":["70F15","70H15"],"pacs":[],"model":"deepseek-v4-flash","headline":"Second-order secular theory reproduces exoplanet precession to 0.05-1%.","keywords":["second-order secular perturbation theory","Lie transform","mean-motion resonance","FFT secular averaging","three-body planetary problem","exoplanet long-term dynamics","orbital precession frequencies","small divisors"],"falsifier":"For GJ 876 or TIC 279401253, raise the resonant truncation from $K'=20$ to $K'=30$ or $40$, recompute the second-order precession frequencies $g_1$ and $g_2$, and compare with n-body values; if $g_2$ shifts by more than the reported few-percent accuracy, the claimed agreement is controlled by the unverified truncation rather than by the second-order theory.","tokens_in":22838,"feed_emoji":"🪐","tokens_out":9206,"duration_ms":80971,"temperature":0.7,"pith_summary":"The paper develops a numerical route to second-order canonical perturbation theory for two planets orbiting a star. Standard secular theories stop at first order in the planetary mass ratio, which distorts the slow precession frequencies of orbits near mean-motion resonances. The authors use a Lie transform to remove the fast orbital angles from the Hamiltonian and compute the required Fourier coefficients with a fast Fourier transform, avoiding expansions in eccentricity, inclination, or semi-major axis ratio. The second-order Hamiltonian reduces relative errors in the dominant precession frequencies from roughly 14-20% to 0.05-1% for the Sun-Jupiter-Saturn system, and from 112% to 2.6% for the resonant system GJ 876. If the results hold, the method extends accurate long-term secular modeling to compact, massive, and resonant exoplanetary systems.","feed_headline":"Second-order secular model cuts precession errors to under 1%","feed_subtitle":"A numerical Lie-transform method brings secular models of resonant exoplanets in line with n-body simulations.","key_machinery":"The load-bearing object is the second-order secular Hamiltonian of Eqs. (20) and (23): sums over non-zero Fourier harmonics $h^k_1$ of the disturbing function, divided by the small-divisor denominators $k\\cdot n'_0$ (or $k\\cdot \\tilde n_0$ in the linearized resonant version) and combined through Poisson brackets over the slow variables. It is built by a Lie transform that eliminates the fast mean longitudes; the Fourier coefficients and their first derivatives are computed by evaluating the disturbing function on a two-dimensional grid in the mean longitudes and applying a 2D FFT, which removes any need for expansions in eccentricities, inclinations, or semi-major axis ratios. The second derivatives needed for the equations of motion are evaluated with a five-point finite-difference stencil.","core_discovery":"The central claim is that a second-order secular Hamiltonian, obtained by Lie-transforming the three-body Hamiltonian and evaluating its Fourier sums numerically, reproduces the long-term secular dynamics of planetary systems near or inside mean-motion resonances much better than first-order theory. The second-order terms contain the small-divisor denominators $k\\cdot n'_0$ that encode near-resonant coupling, so they correct exactly what first-order averaging misses. In the validations, the relative error in the Saturn perihelion frequency $g_6$ drops from about 20% with the first-order model to 1% with the second-order secular model and to 0.05% when the 5:2 resonant harmonics are included; the error in the GJ 876 frequency $g_2$ drops from 112% to 2.6%; and the WASP-148 perihelion frequencies $g_1$ and $g_2$ improve from 12% and 6.5% error to below 0.1%. The paper concludes that second-order terms significantly improve the accuracy of orbital precession frequencies, most notably for systems with strong mean-motion resonances or large planetary masses.","pith_inferences":["A natural extension is to push the resonant truncation limit $K'$ upward for GJ 876 and TIC 279401253: a convergence scan would show whether their residual discrepancies come from the unverified truncation or from neglected third-order mass terms.","The same FFT-plus-Lie-transform construction should extend to more than two planets, since neither the Poisson-bracket sums nor the exponential-decay argument is specific to the three-body case.","For high-eccentricity systems, the mean-motion correction currently absorbs only the zero-eccentricity part of the perturbation; a version that folds in higher-degree terms could further improve systems like TIC.","If the method scales, it could serve as a fast screening tool for the long-term stability of newly discovered resonant exoplanet systems, replacing expensive n-body integrations for population-level studies."],"forward_implications":["Second-order secular models reproduce n-body precession frequencies with relative errors near or below 1% in systems where first-order models err by 10-100%.","Adding the resonant harmonics of a mean-motion resonance, as done for the 5:2 Jupiter-Saturn resonance, removes most of the remaining frequency error, showing that those residuals come from higher-order mass terms associated with the resonance.","Because the method avoids expansions in orbital elements, it can be applied to highly eccentric, inclined, or tightly packed planets where classical secular expansions fail to converge.","For systems locked in a 2:1 resonance, the non-resonant secular model is inadequate, and the second-order resonant secular model is required to track the evolution accurately.","The method isolates the main planetary interactions driving secular dynamics, making it a tool for studying long-term stability of exoplanetary architectures."],"supporting_citations":[{"why":"Supplies the Lie transform formalism the paper uses to eliminate the fast angles and build the secular Hamiltonian order by order.","marker":"[37, 38]"},{"why":"Supplies the technique of modifying the unperturbed Hamiltonian so the small divisors use a corrected mean-motion vector.","marker":"[8, 42]"},{"why":"Supplies the exponential decay of Fourier coefficients that justifies truncating the perturbation series at finite harmonic order.","marker":"[41, 43]"},{"why":"Supplies the fast Fourier transform algorithm used to compute the harmonic coefficients of the disturbing function on a grid.","marker":"[44]"},{"why":"Provides the second-order secular theory for the solar system and the frequency-analysis method used to measure precession frequencies.","marker":"[9]"},{"why":"Provides the initial conditions for the Sun-Jupiter-Saturn n-body model used in the central validation.","marker":"[56]"},{"why":"Supplies the symplectic n-body integrator used as the benchmark against which the secular models are compared.","marker":"[50]"},{"why":"Supplies the frequency-analysis technique used to extract the secular frequencies $g_5$, $g_6$, and $s_6$ from the orbital solutions.","marker":"[58]"}],"fun_headline_variants":["Second-order secular model cuts precession errors tenfold","Numerical Lie transform hits 0.1% error on exoplanet precession","Second-order theory tames resonant exoplanet orbits","Lie-transform method sharpens exoplanet precession predictions"],"cache_read_input_tokens":3200,"weakest_assumption_plain":"The load-bearing premise is that the Fourier coefficients of the perturbation decay exponentially with harmonic order, so truncating the series at $K$ (and at $K'$ for resonant systems) leaves the second-order secular dynamics essentially unchanged, even though convergence is demonstrated numerically only for the Sun-Jupiter-Saturn system.","fun_headline_variants_meta":{"raw":{"variants":["Second-order secular model cuts precession errors tenfold","Numerical Lie transform hits 0.1% error on exoplanet precession","Second-order theory tames resonant exoplanet orbits","Lie-transform method sharpens exoplanet precession predictions"]},"model":"deepseek-v4-flash","effort":"low","cost_usd":0.000187,"raw_usage":{"total_tokens":1372,"prompt_tokens":1029,"completion_tokens":343,"prompt_tokens_details":{"cached_tokens":384},"prompt_cache_hit_tokens":384,"prompt_cache_miss_tokens":645,"completion_tokens_details":{"reasoning_tokens":272}},"tokens_in":645,"tokens_out":343,"duration_ms":3596,"temperature":1.0,"reasoning_tokens":272,"cache_read_input_tokens":384,"cache_creation_input_tokens":0},"cache_creation_input_tokens":0},"created_at":"2026-08-15T19:58:47.120089+00:00","model_set":{"reader":"deepseek-v4-flash"},"falsifier":"For GJ 876 or TIC 279401253, raise the resonant truncation from $K'=20$ to $K'=30$ or $40$, recompute the second-order precession frequencies $g_1$ and $g_2$, and compare with n-body values; if $g_2$ shifts by more than the reported few-percent accuracy, the claimed agreement is controlled by the unverified truncation rather than by the second-order theory.","supporting_citations":[{"cited_title":"The computational cost involves N 2 evaluations per each function, followed by a 2D FFT, with complex- ity 2 N 2 log(N ) per function","cited_arxiv_id":null,"evidence_quote":"Supplies the fast Fourier transform algorithm used to compute the harmonic coefficients of the disturbing function on a grid."},{"cited_title":"Laskar, Astronomy and Astrophysics 198, 341 (1988)","cited_arxiv_id":null,"evidence_quote":"Provides the second-order secular theory for the solar system and the frequency-analysis method used to measure precession frequencies."},{"cited_title":"Varadi, M","cited_arxiv_id":null,"evidence_quote":"Provides the initial conditions for the Sun-Jupiter-Saturn n-body model used in the central validation."},{"cited_title":"Laskar, Astronomy and Astrophysics 287, L9 (1994)","cited_arxiv_id":null,"evidence_quote":"Supplies the symplectic n-body integrator used as the benchmark against which the secular models are compared."},{"cited_title":"Locatelli and A","cited_arxiv_id":null,"evidence_quote":"Supplies the frequency-analysis technique used to extract the secular frequencies $g_5$, $g_6$, and $s_6$ from the orbital solutions."}],"review_version":2}