{"id":"be9fd5de-be2c-4561-accf-a0fe1cdbbba6","arxiv_id":"2502.01496","paper_version":1,"verdict":"CONDITIONAL","confidence":"MODERATE","novelty_score":5.0,"correctness_risk":"medium","formal_verification":"none","parameter_count":2,"one_line_summary":"PHANTOM's fast multipole gravity solver violates total linear momentum conservation by computing leaf-cell-to-supercell forces without reverse pairs; in neutron star simulations this shifts the center of mass by up to about one stellar radius.","lead":"This paper shows that the fast multipole gravity solver in the widely used PHANTOM code does not conserve linear momentum: forces are computed from leaf cells to larger tree cells but not back. In neutron star simulations this causes the star's center of mass to drift by up to about one stellar radius within tens of orbital periods.","discovery_kind":"extension","skeptic_critique":{"model":"deepseek-v4-flash","headline":"The central non-conservation claim rests on an unverified assertion about PHANTOM's tree traversal: only \"leaf cell←super-cell\" pairs are considered, producing the one-sided force in Eq. (30). If PHANTOM symmetrizes cell interactions, the mechanism vanishes.","rationale":"The reader's weakest_assumption identifies exactly the same premise: the one-sided traversal asserted in Section 4 and encoded in Eq. (30) is the load-bearing model. I agree that this premise is unverified against the PHANTOM source and that the quantitative claims (sigma scaling, T^{3/2} law, theta > 0.2 irreducibility) collapse if the traversal is actually symmetric. This is a stronger concern than the random-walk model details or the single-realization statistics: those only affect the quantitative law, whereas the traversal asymmetry is required for there to be any non-conservation at all from the FMM mechanism described. The paper's Section 3 proof that symmetric cell-cell pairs satisfy Newton's third law is an independent and plausible analytic result, but it does not establish the central conclusion without the Section 4 premise. The proposed concrete test is feasible and decisive: read the relevant PHANTOM gravity-tree code and check whether reverse pairs are considered, or equivalently run a small configuration and measure net force. If the premise holds, the paper's central claim is supported and the remaining quantitative issues are secondary; if it fails, the paper's mechanism is absent. Since the reader already issued CONDITIONAL, this stress-test does not change the verdict; it sharpens the condition to a specific, falsifiable check on the code traversal.","tokens_in":17543,"tokens_out":2239,"duration_ms":23922,"concrete_test":"Pin the exact PHANTOM version used for the runs (e.g., the commit hash or release tag) and inspect the FMM gravity traversal source, typically the routine that walks the kd-tree and accumulates gravitational accelerations in src/tree/grav.F90 or equivalent. Verify directly whether, for a leaf cell alpha and a supercell beta, the code computes the beta-to-alpha force and applies it to alpha's particles without ever applying the reciprocal alpha-to-beta force to beta's leaf children. If the code instead considers both directions of each cell-cell pair, or if it uses a symmetric tree walk, then the residual force in Eq. (30) is identically zero and the sigma(theta) scaling is moot.","verdict_should_be":"UNCHANGED","load_bearing_attack":"The load-bearing premise is stated, not demonstrated: Section 4 asserts that in PHANTOM \"such a pair is searched for only for leaf cells, but not for super-cells,\" so that \"direct pairs 'leaf cell←super-cell' will be encountered, but the reverse pairs will never be considered.\" Equation (30) and Fig. 2 build the entire residual-force model on this asymmetry. The subsequent scaling law (47), the exponent sigma fitted in Fig. 8, and the conclusion that theta > 0.2 makes the error irreducible all depend on Eq. (31)-(33) being the actual net force in PHANTOM. The paper never quotes or inspects the relevant code path (the FMM gravity walk in PHANTOM's tree/grav routines); it is a claim about a specific code's traversal, and it is exactly the kind of claim that can be wrong if the implementation symmetrizes interactions, uses a dual tree walk, or applies forces in both directions when accumulating accelerations. If the asymmetry is absent, the residual force is zero at machine precision, the random-walk model does not describe a real effect, and the numerical drift in Fig. 4 must have a different origin. The theta = 0 control shown in Fig. 7 does not rescue this premise: it only shows that long-range interactions are involved in whatever causes the drift, not that the cause is the particular one-sided traversal assumed in Eq. (30). Thus the central claim is conditional on an unverified fact about PHANTOM internals.","agreement_with_reader":"agree"},"referee_report":{"model":"deepseek-v4-flash","summary":"The paper examines the fast multipole method (FMM) implementation in the PHANTOM smoothed particle hydrodynamics code for self-gravity. The authors prove that for any pair of well-separated kd-tree cells, the force computed from PHANTOM's multipole acceleration formulas satisfies Newton's third law. They then argue that, despite this pair-level symmetry, the global linear momentum is not conserved because PHANTOM's tree traversal only considers 'leaf cell←super-cell' pairs in one direction and never the reverse pairs. This asymmetric treatment is modeled as generating an uncompensated random force on the system, leading to a center-of-mass displacement that grows as T^{3/2} (a 'random memory walk'). Numerical simulations of a single neutron star with N from 5000 to 150000 at tree opening parameter θ=0.5 show a displacement of up to ~1.75 stellar radii over 4×10^4 time units, while a θ=0 (direct summation) control shows no such drift. The authors derive scaling laws for the displacement as a function of mass, radius, particle number, and time, and conclude that for θ≳0.2 increasing N does not reduce the error for neutron-star-scale objects.","tokens_in":17851,"tokens_out":7694,"duration_ms":70087,"significance":"If the premise about PHANTOM's traversal asymmetry is correct, this paper identifies a serious and actionable issue for a widely used astrophysics code: spurious center-of-mass migration in neutron-star simulations, with the error growing with time and not decreasing with particle number for typical opening angles. The pair-level Newton's-third-law proof is a clean analytical result, and the θ=0 control convincingly ties the observed drift to the tree-based long-range force rather than to SPH discretization. The scaling formulas provide a practical guide for error estimation in other objects (white dwarfs, stars). However, the central claim rests on an unverified assertion about the internal traversal logic of a specific code, and the abstract overstates the evidence by referring to a binary-neutron-star simulation that was not performed. The paper therefore has the potential to be an important warning to the PHANTOM user community, but the key premise and supporting models need to be substantiated before publication.","major_comments":[{"comment":"The premise that PHANTOM's gravity tree walk only considers 'leaf cell←super-cell' pairs and never the reverse pairs is asserted in the text but not demonstrated against the actual code. The entire residual-force model in Eq. (30) and the scaling relation (47) depend on this one-sided traversal. The authors should cite the specific PHANTOM version and quote or inspect the relevant tree-walk routine (or provide a simple instrumented test that counts interaction pairs). Without this verification, the mechanism remains a hypothesis about the code's behavior, and the simulated drift could have an alternative cause, such as truncation error in the multipole expansion.","section":"Section 4, Eq. (30)"},{"comment":"The abstract and conclusion state that for a pair of neutron stars the center-of-mass displacement is comparable to the object radii within a few tens of Keplerian revolutions, but no binary simulation is presented. The reported simulations are of a single neutron star. The extrapolation from a single star to a binary system requires justification (e.g., showing that the spurious force on a binary's center of mass has the same statistics), or the abstract should be rephrased to state explicitly that the effect is demonstrated for an isolated star and expected to be similar in a binary.","section":"Abstract and Section 7"},{"comment":"The proof that PHANTOM's acceleration formula (8) is equivalent to a symmetric second-order Taylor expansion relies on the claim that the added terms (n,m)=(1,2),(2,1),(2,2) lie inside the residual term R_3(r,x) of the p=3 expansion. No explicit bound for R_3 is provided. Without such a bound, the pair-level conservation proof is incomplete, and the alternative explanation that the observed momentum drift originates from the truncation of the multipole expansion is not ruled out. The authors should either supply a rigorous error estimate or validate the equivalence numerically for the cell configurations used in the simulations.","section":"Section 3, Eqs. (27)-(29)"},{"comment":"The random-walk model assumes that the uncompensated force at each step is isotropic, independent of previous steps, and uniformly distributed in modulus up to F. These assumptions are not tested. Because the kd-tree is rebuilt deterministically as the particles move, successive residual forces may be correlated, which would invalidate the T^{3/2} scaling. The authors should measure the autocorrelation of the force sequence and compare the predicted displacement distribution to the simulations without fitting F from the same data (e.g., by estimating F from a separate calibration run and then predicting the drift).","section":"Section 5, Eq. (35)"},{"comment":"The exponent sigma, which controls the N-dependence of the displacement through N^{sigma-1/6}, is fitted from the same set of runs that exhibit the drift. The conclusion that 'for theta ≳ 0.2, no changes in N can eliminate the error' hinges on this fit. The paper should report fit uncertainties and, ideally, validate the predicted scaling by comparing the N^{sigma-1/6} prediction directly against the measured displacement amplitudes in Fig. 4, or by running additional N values. The direct observation in Fig. 4 shows no improvement over the limited range 5000-150000, but that does not justify the strong claim 'for any N'.","section":"Section 6, Fig. 8 and Eq. (51)"}],"minor_comments":[{"comment":"The PHANTOM version and configuration used for the simulations are not stated; this information is essential for reproducibility of the claim about the tree traversal.","section":"Throughout"},{"comment":"The notation 'a^{β→α,i}(r,x)' conflates the acceleration of a particle with the acceleration of a cell; consider renaming the cell-level force or clarifying the distinction in the text.","section":"Section 2, Eq. (8)"},{"comment":"The curves for the six N values are not labeled with distinct markers in the text description; a table of final displacements (or a legend with markers) would help the reader evaluate the N-dependence directly.","section":"Section 6, Fig. 4"},{"comment":"Reference [24] is a preprint; if it has appeared in a peer-reviewed venue, it should be updated.","section":"References"}],"recommendation":"major_revision","confidential_remarks":"The paper makes a strong claim about a specific, widely used code. The editor should consider whether the authors have adequately verified the internal behavior of PHANTOM; independent confirmation by a PHANTOM developer or a code-level inspection would substantially de-risk the paper. The title's phrase 'widely used hierarchical methods' is broader than the scope, which is limited to PHANTOM's specific traversal; this should be reflected either in the title or in the framing. If the traversal premise is confirmed, the paper would be a valuable contribution; as it stands, the central conclusion is conditional on an unverified implementation detail."},"author_rebuttal":null,"desk_editor":{"model":"deepseek-v4-flash","letter":"Dear [Colleague],\n\nThe short version: this paper claims that PHANTOM's FMM self-gravity implementation violates total linear momentum because its tree walk only considers 'leaf cell ← super-cell' pairs, never the reverse. The numerical test shows a single neutron star drifting by up to ~1.75 stellar radii over 4×10^4 time units at the default theta=0.5. If that diagnosis is correct, it is a real problem for compact-object simulations.\n\nWhat is new and good: the pairwise Newton's third law proof is clean; the theta=0 control correctly separates the FMM artifact from SPH discretization; and the residual-force formulas (31)-(33) give a concrete analytic handle on the effect. The paper also credits Dehnen (2000, 2002) for the general idea of momentum-conserving tree walks, so the contribution is framed as a specific code diagnosis rather than a new universal principle.\n\nThe soft spots are real, but not fatal. The load-bearing premise is the claim that PHANTOM's traversal is one-sided. That is asserted in Section 4 without a code citation or inspection. Barnes-Hut per-particle walks are typically asymmetric, so the claim is plausible, but the authors need to demonstrate it against a pinned PHANTOM version. If the implementation symmetrizes cell-cell interactions, the residual force formula does not describe the actual error. The random-walk model with independent isotropic kicks is a reasonable starting point, but it is not tested; the sigma exponent is fitted from the same runs, no error bars, and the binary extrapolation in the abstract comes from a single-star calculation.\n\nNone of this destroys the central argument; it is a call for more rigor. The paper is exactly the kind of work that deserves referee time: it is specific, testable, and potentially important for code users. My recommendation is to send it to peer review with a request that the referee verify the tree traversal, run repeated realizations, and test the binary case directly. I would not cite it as established fact until that is done.\n\nWorth bringing to the reading group, though, for the discussion of tree code symmetries.","headline":"Plausible, potentially important diagnosis of PHANTOM's FMM momentum drift, but the key traversal claim is unverified.","tokens_in":18395,"tokens_out":3593,"would_cite":false,"duration_ms":31506,"reading_group":"yes","serious_thinker":"yes","would_accept_peer_review":true},"rs_alignment":null,"lean_confirmation":null,"pith_extraction":{"msc":[],"pacs":[],"model":"deepseek-v4-flash","headline":"PHANTOM's fast multipole method conserves momentum pairwise but not globally, causing neutron-star centers of mass to drift by up to about 1.75 stellar radii.","keywords":["fast multipole method","PHANTOM code","kd-tree","momentum conservation","self-gravity","neutron star","tree opening parameter","random walk"],"falsifier":"Instrument the PHANTOM force loop to accumulate the net force on every kd-tree cell pair that interacts in a single time step; if the sum over all pairs is zero to round-off, the paper's claimed residual-force mechanism is absent. Alternatively, run the same neutron-star relaxation at $\\theta=0.5$ with a version of the traversal that adds reverse \"super-cell → leaf cell\" pairs; if the center-of-mass displacement persists at the same level, the explanation would need revision.","tokens_in":17283,"feed_emoji":"⭐","tokens_out":9487,"duration_ms":71408,"temperature":0.7,"pith_summary":"This paper argues that the fast multipole method (FMM) as implemented in the PHANTOM code for self-gravitating gas dynamics fails to conserve total linear momentum, even though Newton's third law holds for every pair of kd-tree cells engaged in mutual interaction. The cause is a one-sided tree traversal: PHANTOM computes long-range forces only for \"leaf cell ← super-cell\" pairs, never the reverse pairs, so the forces acting on a super-cell's children are not balanced by forces acting on the leaf cell. The resulting uncompensated force acts like a random kick applied at every time step, producing a center-of-mass migration that grows as $T^{3/2}$ in time. For neutron-star-scale objects, the displacement becomes comparable to the stellar radius within a few tens of Keplerian revolutions, and for tree opening parameter $\\theta \\gtrsim 0.2$ it cannot be reduced by increasing particle number. The paper therefore calls for replacing PHANTOM's traversal with a symmetric tree-walk scheme.","feed_headline":"Neutron-star centers drift under PHANTOM's FMM gravity","feed_subtitle":"A one-way kd-tree traversal creates a random force; displacement grows as time^{3/2} and more particles don't help.","key_machinery":"The central mechanism is the asymmetric pair \"leaf cell ← super-cell\" in the kd-tree traversal. For a leaf cell $\\alpha$ and a supercell $\\beta$ made of leaf children $\\beta_1,\\beta_2$, the total force on the system is taken as $F^{\\beta\\to\\alpha}$ plus $F^{\\alpha\\to\\beta_1}+F^{\\alpha\\to\\beta_2}$, with no reciprocal $F^{\\alpha\\to\\beta}$ to balance $F^{\\beta\\to\\alpha}$. Expanding the multipole forces in the displacement $x$ between the supercell's center of mass and its children gives a leading residual force $F\\simeq \\frac{4}{r^5} M_\\alpha M_\\beta \\frac{M_{\\beta_1}}{M_{\\beta_2}}\\left(1-\\frac{M_{\\beta_1}}{M_{\\beta_2}}\\right) l^3$, independent of internal quadrupole moments. This residual force is a random variable whose kicks accumulate as a \"random memory walk\", yielding $\\sqrt{\\langle r^2\\rangle}\\propto T^{3/2}$.","core_discovery":"The paper's central discovery is that the FMM force calculation in PHANTOM, although pairwise symmetric between any two cells, is not symmetric at the system level. In PHANTOM's traversal, a leaf cell receives long-range acceleration from a supercell, but the supercell does not receive the reciprocal acceleration; instead, the leaf children of the supercell receive their own accelerations from the leaf cell. Because the supercell's children are displaced from the supercell's center of mass, these two sets of forces do not cancel. Section 4 derives the residual force $F \\simeq \\frac{4}{r^5} M_\\alpha M_\\beta \\frac{M_{\\beta_1}}{M_{\\beta_2}}\\left(1-\\frac{M_{\\beta_1}}{M_{\\beta_2}}\\right) l^3$ for the elementary asymmetric pair, independent of the cells' quadrupole moments. The net effect is a random \"kick\" on the whole system, whose mean-square displacement satisfies $\\sqrt{\\langle r^2\\rangle} = \\frac{F}{M}\\sqrt{\\frac{\\Delta t}{3}}\\,T^{3/2}$. Numerical relaxation of a single neutron star with $N=5000$ to $150000$ particles at $\\theta=0.5$ shows center-of-mass displacements up to about $1.75$ stellar radii in $4\\times 10^4$ time units, while direct summation at $\\theta=0$ produces no such drift. The paper concludes that the code must be corrected, for example by adopting the symmetric tree-walk scheme of [9;10].","pith_inferences":["Inference: The same residual-force mechanism should also make the reported angular-momentum non-conservation grow with time; an isolated star with zero initial velocities would acquire a slow net spin, which could masquerade as physical rotation in merger studies.","Inference: Because the kick magnitude depends on cell geometry through the mass imbalance $M_{\\beta_1}/M_{\\beta_2}$, the drift is configuration-dependent, so unequal-mass or offset-density binaries would drift at different rates and binary orbital elements inferred from PHANTOM relaxation runs may carry a systematic bias.","Inference: The equal-time-step assumption in the random-walk model is an idealization; with PHANTOM's adaptive CFL timesteps the $\\Delta t$ entering (35) varies, so the clean $T^{3/2}$ scaling is a limiting law rather than an exact prediction. Runs with artificially fixed timesteps could test how strongly the scaling changes."],"forward_implications":["Current PHANTOM FMM results for neutron-star-scale objects carry a center-of-mass drift on the order of the stellar radius over tens of Keplerian orbits at $\\theta=0.5$.","Increasing the number of particles does not reduce the drift for $\\theta \\gtrsim 0.2$; the displacement scales with particle number as $N^{\\sigma-1/6}$, with $\\sigma\\approx 0.55$ at $\\theta=0.5$.","For smaller $\\theta\\lesssim 0.2$ the drift decreases only as $N^{-1/6}$, so cutting the error by a factor of ten would require six orders of magnitude more particles, and the increased short-range region makes the calculation significantly slower.","White dwarfs, hot stars, and nebulae are largely unaffected over the same timescale, but white-dwarf binaries on very long runs would accumulate a comparable error.","The proposed correction is a symmetric tree-walk that considers all well-separated cell pairs with a symmetric opening criterion, while noting that FMM schemes in general still have known angular-momentum non-conservation."],"supporting_citations":[{"why":"Supplies the symmetric tree-walk algorithm that conserves momentum and is proposed as the corrective standard against which PHANTOM's one-sided traversal is compared.","marker":"[9; 10]"},{"why":"Defines the PHANTOM SPH self-gravity force calculation, equations of motion, and time-stepping that the paper analyzes and tests.","marker":"[17]"},{"why":"The authors' earlier preprint on momentum non-conservation in PHANTOM FMM, which this article complements and revises.","marker":"[24]"},{"why":"Provides the adaptive gravitational softening formalism and the symmetric Green's function used in the proof that Newton's third law holds for a single cell pair.","marker":"[25]"},{"why":"Gives the random-walk 'random memory walk' model that yields the $T^{3/2}$ displacement scaling from velocity kicks.","marker":"[32]"},{"why":"Provides the polytropic equation of state used to build the neutron-star models in the numerical tests.","marker":"[35]"}],"fun_headline_variants":["PHANTOM's FMM violates momentum conservation","Neutron stars drift from FMM's asymmetric forces","FMM gravity kicks neutron stars: momentum not conserved","Hierarchical gravity breaks Newton's third law","PHANTOM FMM: no momentum conservation, random drift"],"cache_read_input_tokens":3200,"weakest_assumption_plain":"The load-bearing premise is that PHANTOM's actual tree traversal applies long-range forces in only one direction, \"leaf cell ← super-cell\", and never the reverse; if the code in fact symmetrizes those interactions, the derived residual-force formulas and the predicted drift would not describe its real momentum error.","fun_headline_variants_meta":{"raw":{"variants":["PHANTOM's FMM violates momentum conservation","Neutron stars drift from FMM's asymmetric forces","FMM gravity kicks neutron stars: momentum not conserved","Hierarchical gravity breaks Newton's third law","PHANTOM FMM: no momentum conservation, random drift"]},"model":"deepseek-v4-flash","effort":"low","cost_usd":0.000266,"raw_usage":{"total_tokens":1663,"prompt_tokens":1053,"completion_tokens":610,"prompt_tokens_details":{"cached_tokens":384},"prompt_cache_hit_tokens":384,"prompt_cache_miss_tokens":669,"completion_tokens_details":{"reasoning_tokens":534}},"tokens_in":669,"tokens_out":610,"duration_ms":6045,"temperature":1.0,"reasoning_tokens":534,"cache_read_input_tokens":384,"cache_creation_input_tokens":0},"cache_creation_input_tokens":0},"created_at":"2026-08-09T15:08:35.681083+00:00","model_set":{"reader":"deepseek-v4-flash"},"falsifier":"Instrument the PHANTOM force loop to accumulate the net force on every kd-tree cell pair that interacts in a single time step; if the sum over all pairs is zero to round-off, the paper's claimed residual-force mechanism is absent. Alternatively, run the same neutron-star relaxation at $\\theta=0.5$ with a version of the traversal that adds reverse \"super-cell → leaf cell\" pairs; if the center-of-mass displacement persists at the same level, the explanation would need revision.","supporting_citations":[{"cited_title":null,"cited_arxiv_id":null,"evidence_quote":"The authors' earlier preprint on momentum non-conservation in PHANTOM FMM, which this article complements and revises."},{"cited_title":"J., Monaghan J","cited_arxiv_id":null,"evidence_quote":"Provides the adaptive gravitational softening formalism and the symmetric Green's function used in the proof that Newton's third law holds for a single cell pair."},{"cited_title":"Diffusion","cited_arxiv_id":null,"evidence_quote":"Gives the random-walk 'random memory walk' model that yields the $T^{3/2}$ displacement scaling from velocity kicks."},{"cited_title":"K., Hebeler K., Lattimer J","cited_arxiv_id":null,"evidence_quote":"Provides the polytropic equation of state used to build the neutron-star models in the numerical tests."}],"review_version":1}