arxiv: 2605.13461 · v1 · submitted 2026-05-13 · 🌌 astro-ph.EP

Recognition: 2 theorem links

· Lean Theorem

Self-gravity in thin protoplanetary discs: 2. Numerical convergence solved and revealing the overestimation in mass of formed planets with softening

S. Rendon Restrepo

Authors on Pith no claims yet

Pith reviewed 2026-05-14 18:13 UTC · model grok-4.3

classification 🌌 astro-ph.EP

keywords protoplanetary discsgravitational instabilityself-gravityplanet formationnumerical simulationsBessel kernelsoftening length

0 comments

The pith

The 2D Bessel kernel for self-gravity resolves numerical convergence in thin disc simulations and shows that softening overestimates planet masses by a factor of two to three.

A machine-rendered reading of the paper's core claim, the machinery that carries it, and where it could break.

This paper introduces and tests a first-principles 2D Bessel kernel for calculating self-gravity in thin protoplanetary discs. The kernel incorporates a natural transition scale from three-dimensional to two-dimensional gravitational behavior. Simulations using this kernel converge properly at high resolution, unlike those with standard softening lengths. Small softening values produce too many fragments whose masses are inflated by a factor of two to three, while large softening suppresses fragmentation entirely. The Bessel approach yields bound fragments that match expected physical behavior better than softening methods.

Core claim

The 2D Bessel formalism of gravity effectively resolves the convergence issues encountered in 2D simulations of gravitational instability. Simulations with the Bessel kernel produce fewer fragments than those with small softening parameters, and the final fragment masses are overestimated by a factor of 2-3 when softening is used. A softening parameter of 0.6 H fails to keep fragments bound, unlike the Bessel method.

What carries the argument

The 2D Bessel kernel prescription for self-gravity, which introduces a characteristic length below which gravity transitions smoothly from 3D to 2D scaling.

If this is right

2D simulations of gravitational instability must use the Bessel kernel to achieve numerical convergence.
Standard softening prescriptions overestimate the number of fragments and their masses by a factor of 2-3.
Large softening parameters inhibit gravitational effects and prevent bound fragment formation.
Adopting the Bessel prescription ensures consistent and accurate treatment of gravity in thin discs.

Where Pith is reading between the lines

These are editorial extensions of the paper, not claims the author makes directly.

Previous 2D simulations using softening may have produced unrealistically high rates of planet formation through gravitational instability.
Observational estimates of early planet masses in massive discs could be revised downward if Bessel gravity is the correct approach.
Hybrid 2D-3D models of disc evolution might benefit from incorporating the Bessel transition scale for improved accuracy.

Load-bearing premise

The first-principles Bessel kernel accurately captures the vertical structure and the 3D-to-2D transition in real thin discs without direct comparison to full three-dimensional simulations.

What would settle it

Running equivalent high-resolution simulations of gravitational instability in full 3D and comparing the number and masses of formed fragments to the 2D Bessel results would test whether the overestimation is real.

Figures

Figures reproduced from arXiv: 2605.13461 by S. Rendon Restrepo.

**Figure 1.** Figure 1: Numerical convergence of 2D Gravitational instability simulations using the gravitational Bessel Kernel. The solid line represents a sigmoid fit to the critical cooling threshold that separates fragmentation from gravito-turbulence. Above β = 5 none of my simulations fragments. 3. Numerical convergence of 2D global simulations with the Bessel kernel As it is understood from the literature, the numerical … view at source ↗

**Figure 2.** Figure 2: Numerical convergence of the properties of the fragmentation regime when using the Bessel kernel for β = 1, 2, 3. The black dash line indicates the fragmentation threshold. 4. Characterisation of GI for different gravity prescriptions Now that I have demonstrated the convergence of GI simulations using the Bessel kernel, my next objective is to compare the results obtained with this approach to those usin… view at source ↗

**Figure 3.** Figure 3: Numerical convergence of the time and space averaged Reynolds and gravitational stresses in the regime of gravito-turbulence for different cooling rates. remain within the box subsequently undergo successive mergers, ultimately forming a single, extremely massive stellar-type object with a mass in the range of 0.15-0.2 M⊙, which is more than 50% of the initial mass of the disc. By the end of the simulat… view at source ↗

**Figure 4.** Figure 4: Comparison of the fragmentation regime (β = 2) for the different gravity prescriptions at different times. Fragments are highlighted with green circles. For ϵ/Hrms = 1.2 (not shown here), no fragments form. Article number, page 8 of 14 [PITH_FULL_IMAGE:figures/full_fig_p008_4.png] view at source ↗

**Figure 5.** Figure 5: Properties of the fragmentation regime for the different gravity prescriptions. The black dash line indicate the fragmentation threshold. for all the gravity prescriptions under consideration. The mass of the formed fragment in the unsoftened simulation is shown in [PITH_FULL_IMAGE:figures/full_fig_p009_5.png] view at source ↗

**Figure 6.** Figure 6: Evolution of gravitationally bound objects under different gravity prescriptions for β = 2. The solid line represents the total mass, while cross markers denote individual fragments. The brown dashed lines indicate the lower and upper mass limits for brown dwarf formation. For ϵ/Hrms = 0.0, no fragments form, and for ϵ/Hrms = 0.6, the total mass may be overestimated beyond ∼ 6.7 kyr (see Sect. 4.1). For … view at source ↗

**Figure 7.** Figure 7: Detailed comparison of fragment evolution using the Bessel kernel (left) and the Plummer potential paradigm with ϵ/Hrms = 0.6 (right). The snapshots, from top to bottom, correspond to t = 6.7, 7.7, 8.05, 8.1 and 8.6 kyr. All plots are centered at the peak density. Article number, page 12 of 14 [PITH_FULL_IMAGE:figures/full_fig_p012_7.png] view at source ↗

**Figure 8.** Figure 8: Comparison of the gravito-turbulent regime (β = 8) for the different gravity prescriptions [PITH_FULL_IMAGE:figures/full_fig_p013_8.png] view at source ↗

**Figure 9.** Figure 9: Evolution of the gravitationally bound objects for ϵ/Hrmrms = 0 and β = 8. Article number, page 13 of 14 [PITH_FULL_IMAGE:figures/full_fig_p013_9.png] view at source ↗

read the original abstract

The Gravitational Instability (GI) is a leading theory for explaining early planet formation in massive discs. In the early 2010s, 3D SPH simulations of GI failed to converge, initially attributed to resolution-dependent viscosity but later appearing in 2D SPH and grid-based simulations, suggesting a numerical artifact inherent to the 2D approximation of gravity. Recently, we derived from first principles a much improved prescription for gravity in 2D discs (via a Bessel kernel). This prescription introduces a characteristic length below which gravity smoothly transitions from a 3D to a 2D scaling. This cannot be captured by standard smoothing length approaches, widely used in 2D simulations. We employ this new prescription to resolve the convergence issue of GI in 2D, and compare the outcomes of GI in runs using the Bessel kernel with those obtained using softening prescriptions at high resolution. We conducted numerical simulations with the FargoCPT code, where the Bessel prescription was implemented. The 2D Bessel formalism of gravity effectively resolves the convergence issues encountered in 2D simulations. When compared to simulations employing softened or unsoftened potentials, I observe that a small softening parameter tends to overestimate gravitational effects. This results in an artificially high number of fragments, leading to final fragment masses that are overestimated by a factor of 2-3. Conversely, employing large softening parameters inhibits gravitational effects. Although our analysis initially suggests that a softening parameter of 0.6 H might offer the best compromise, in reality, the resulting fragments fail to remain gravitationally bound-a behavior not observed when using the Bessel kernel. Our findings strongly suggest that the Bessel prescription should be adopted to ensure a consistent and accurate treatment of gravity in thin discs.

Editorial analysis

A structured set of objections, weighed in public.

Desk editor's note, referee report, simulated authors' rebuttal, and a circularity audit. Tearing a paper down is the easy half of reading it; the pith above is the substance, this is the friction.

Desk Editor's Note private letter to a colleague

Bessel kernel fixes 2D GI convergence and shows softening overestimates fragment masses by 2-3x, but lacks 3D benchmarks.

read the letter

This paper shows that a Bessel kernel derived from first principles fixes convergence in 2D gravitational instability simulations of thin discs. It also indicates that common softening prescriptions overestimate the final masses of formed planets by a factor of two to three. The new element is the kernel itself, which builds in a length scale for the 3D-to-2D gravity transition that softening cannot match. They put this into FargoCPT and ran high-resolution comparisons. The Bessel runs converge nicely, while small softening lengths create excess fragments and thus higher masses. Big softening lengths stop fragments from staying bound. The numerical evidence for convergence and the mass difference looks like a real improvement over prior work. It gives modelers a better default for 2D gravity. The soft spot is the missing link to 3D. The kernel assumes an idealized vertical profile, and real discs have more complexity. Without 3D reference runs, the claim that softening overestimates could flip if the Bessel underestimates instead. Details on exact resolutions and measurement methods are also thin in the abstract. This is for anyone simulating massive discs in 2D to study early planet formation. It offers a practical fix for a known numerical problem. I would send it for peer review. The idea is grounded enough and the results are specific enough to warrant referee input on the derivation and the comparisons.

Referee Report

1 major / 2 minor

Summary. The paper claims that a first-principles 2D Bessel kernel for self-gravity in thin protoplanetary discs resolves long-standing convergence problems in gravitational instability (GI) simulations, and that standard softening prescriptions (especially small softening lengths) produce an artificially high number of fragments, overestimating final planet masses by a factor of 2-3 compared to the Bessel results; large softening inhibits bound fragments, and the Bessel approach is recommended for accurate 2D modeling.

Significance. If the central claim holds, the work is significant for the field: it supplies a physically derived, parameter-free gravity prescription that addresses a decade-old numerical artifact in 2D GI simulations, potentially yielding more reliable fragment-mass predictions for early planet formation and standardizing how self-gravity is treated in thin-disc codes.

major comments (1)

[Abstract] Abstract and results comparison: the claim that small softening overestimates fragment masses by a factor of 2-3 is load-bearing for the recommendation to adopt the Bessel kernel, yet the manuscript provides no high-resolution 3D simulations to anchor which 2D result is closer to physical reality; the derivation assumes an idealized vertical structure and instantaneous 3D-to-2D transition, so the reported overestimate cannot yet be distinguished from possible underestimation by the Bessel kernel itself.

minor comments (2)

[Abstract] The abstract states convergence is achieved but supplies no error bars, resolution values, or quantitative metrics (e.g., L2 norms or fragment-mass histograms) to support the factor-of-2-3 claim; these must be added in the methods and results sections.
[Methods] Implementation details for the Bessel kernel in FargoCPT (e.g., how the characteristic length scale is computed on the grid and any additional numerical parameters) are not specified; this information is required for reproducibility.

Simulated Author's Rebuttal

1 responses · 0 unresolved

We thank the referee for their constructive review and for acknowledging the potential significance of resolving convergence issues in 2D GI simulations. We address the major comment below and will revise the manuscript to incorporate appropriate caveats while preserving the core findings.

read point-by-point responses

Referee: [Abstract] Abstract and results comparison: the claim that small softening overestimates fragment masses by a factor of 2-3 is load-bearing for the recommendation to adopt the Bessel kernel, yet the manuscript provides no high-resolution 3D simulations to anchor which 2D result is closer to physical reality; the derivation assumes an idealized vertical structure and instantaneous 3D-to-2D transition, so the reported overestimate cannot yet be distinguished from possible underestimation by the Bessel kernel itself.

Authors: We agree that high-resolution 3D simulations would be required to definitively determine which 2D prescription yields masses closer to physical reality. Our work derives the Bessel kernel from first principles for thin discs and demonstrates that it eliminates the long-standing non-convergence artifact seen with softening prescriptions. The reported factor of 2-3 is a relative difference between converged Bessel results and non-converged softening outcomes in our 2D runs. We will revise the abstract to explicitly state that the overestimate is relative to the Bessel kernel and add a new discussion paragraph on the idealized assumptions in the derivation together with the need for future 3D benchmarking. This revision clarifies the scope without changing our recommendation that the Bessel kernel is the appropriate choice for accurate 2D thin-disc modeling. revision: yes

Circularity Check

0 steps flagged

No significant circularity; Bessel kernel from independent first-principles derivation

full rationale

The paper states that the 2D Bessel kernel was derived from first principles in prior work and is applied here without refitting or redefinition to the current simulation results. The central claims (convergence resolution and factor-of-2-3 mass overestimation under softening) rest on direct numerical comparison between the kernel runs and softening runs, which are independent of the kernel's internal construction. No equation or claim reduces by construction to a self-definition, fitted parameter renamed as prediction, or load-bearing self-citation chain; the derivation is treated as an external mathematical input whose validity is not re-proven inside this manuscript. Self-citation of the kernel derivation is present but does not create circularity under the stated criteria.

Axiom & Free-Parameter Ledger

0 free parameters · 1 axioms · 0 invented entities

The central claim rests on the validity of the first-principles derivation of the Bessel kernel for thin discs and the assumption that 2D simulations with this kernel faithfully represent 3D gravitational instability outcomes.

axioms (1)

domain assumption Thin-disc approximation permits a 2D treatment of gravity with a smooth 3D-to-2D transition at a characteristic length scale
Invoked in the derivation of the Bessel kernel as stated in the abstract

pith-pipeline@v0.9.0 · 5629 in / 1240 out tokens · 42688 ms · 2026-05-14T18:13:45.099541+00:00 · methodology

discussion (0)

Lean theorems connected to this paper

Citations machine-checked in the Pith Canon. Every link opens the source theorem in the public Lean library.

IndisputableMonolith/Foundation/AlexanderDuality.lean alexander_duality_circle_linking unclear

?

unclear
Relation between the paper passage and the cited Recognition theorem.

Recently, we derived from first principles a much improved prescription for gravity in 2D discs (via a Bessel kernel). This prescription introduces a characteristic length, H_rms, below which gravity smoothly transitions from a 3D to a 2D scaling.
IndisputableMonolith/Cost/FunctionalEquation.lean washburn_uniqueness_aczel unclear

?

unclear
Relation between the paper passage and the cited Recognition theorem.

The 2D Bessel formalism of gravity effectively resolves the convergence issues encountered in 2D simulations... final fragment masses that are overestimated by a factor of 2-3.

What do these tags mean?

matches: The paper's claim is directly supported by a theorem in the formal canon.
supports: The theorem supports part of the paper's argument, but the paper may add assumptions or extra steps.
extends: The paper goes beyond the formal theorem; the theorem is a base layer rather than the whole result.
uses: The paper appears to rely on the theorem as machinery.
contradicts: The paper's claim conflicts with a theorem or certificate in the canon.
unclear: Pith found a possible connection, but the passage is too broad, indirect, or ambiguous to say the theorem truly supports the claim.

Reference graph

Works this paper leans on

57 extracted references · 3 canonical work pages

[1]

C., Ruden, S

Adams, F. C., Ruden, S. P., & Shu, F. H. 1989, ApJ, 347, 959

1989
[2]

2011, in Advances in Imaging and Electron Physics, V ol

Baddour, N. 2011, in Advances in Imaging and Electron Physics, V ol. 165, Ad- vances in Imaging and Electron Physics, ed. P. W. Hawkes (Elsevier), 1–45

2011
[3]

& Klahr, H

Baehr, H. & Klahr, H. 2015, The Astrophysical Journal, 814, 155

2015
[4]

Baehr, H., Klahr, H., & Kratter, K. M. 2017, ApJ, 848, 40

2017
[5]

& Masset, F

Baruteau, C. & Masset, F. 2008, ApJ, 678, 483

2008
[6]

2011, MNRAS, 416, 1971

Baruteau, C., Meru, F., & Paardekooper, S.-J. 2011, MNRAS, 416, 1971

2011
[7]

& Zhu, Z

Baruteau, C. & Zhu, Z. 2016, Monthly Notices of the Royal Astronomical Soci- ety, 458, 3927

2016
[8]

Beckman, P. G. & O’Neil, M. 2024, arXiv e-prints, arXiv:2411.09583 Benítez-Llambay, P. & Masset, F. S. 2016, ApJS, 223, 11 Béthune, W., Latter, H., & Kley, W. 2021, A&A, 650, A49

work page arXiv 2024
[9]

Boss, A. P. 1997, Science, 276, 1836

1997
[10]

Brown, J. J. & Ogilvie, G. I. 2024, MNRAS, 534, 39 Article number, page 11 of 14 A&A proofs:manuscript no. aanda ´Calovi´c, A., Nayakshin, S., Casewell, S., & Miret-Roig, N. 2026, MNRAS, 545, staf2097

2024
[11]

J., Ziampras, A., Brown, J

Cordwell, A. J., Ziampras, A., Brown, J. J., & Rafikov, R. R. 2025, MNRAS, 543, 4198

2025
[12]

& Bhatta, D

Debnath, L. & Bhatta, D. 2014, Integral Transforms and Their Applications, 3rd edn. (Chapman and Hall/CRC)

2014
[13]

2017, ApJ, 847, 43

Deng, H., Mayer, L., & Meru, F. 2017, ApJ, 847, 43

2017
[14]

C., & Klessen, R

Federrath, C., Banerjee, R., Clark, P. C., & Klessen, R. S. 2010, ApJ, 713, 269

2010
[15]

& Rice, K

Forgan, D. & Rice, K. 2013, MNRAS, 432, 3168

2013
[16]

Gammie, C. F. 1996, ApJ, 457, 355

1996
[17]

Gammie, C. F. 2001, ApJ, 553, 174

2001
[18]

G., Mamatsashvili, G

Gibbons, P. G., Mamatsashvili, G. R., & Rice, W. K. M. 2014, MNRAS, 442, 361

2014
[19]

& Eastwood, J

Hockney, R. & Eastwood, J. 2021, Computer Simulation Using Particles (CRC Press)

2021
[20]

& Schreiber, A

Klahr, H. & Schreiber, A. 2020, ApJ, 901, 54

2020
[21]

F., Jung, M., & Duschl, W

Klee, J., Illenseer, T. F., Jung, M., & Duschl, W. J. 2017, A&A, 606, A70

2017
[22]

F., Jung, M., & Duschl, W

Klee, J., Illenseer, T. F., Jung, M., & Duschl, W. J. 2019, A&A, 632, A35

2019
[23]

& Lodato, G

Kratter, K. & Lodato, G. 2016, ARA&A, 54, 271

2016
[24]

M., Murray-Clay, R

Kratter, K. M., Murray-Clay, R. A., & Youdin, A. N. 2010, ApJ, 710, 1375

2010
[25]

2014, MNRAS, 440, 683

Lega, E., Crida, A., Bitsch, B., & Morbidelli, A. 2014, MNRAS, 440, 683

2014
[26]

Lin, C. C. & Shu, F. H. 1964, ApJ, 140, 646

1964
[27]

& Kratter, K

Lin, M.-K. & Kratter, K. M. 2016, ApJ, 824, 91

2016
[28]

& Clarke, C

Lodato, G. & Clarke, C. J. 2011, MNRAS, 413, 2735

2011
[29]

Longarini, C., Lodato, G., Bertin, G., & Armitage, P. J. 2023b, MNRAS, 519, 2017

2017
[30]

J., Kratter, K

Longarini, C., Price, D. J., Kratter, K. M., Lodato, G., & Clarke, C. J. 2025, MNRAS, 541, 1145

2025
[31]

Lovelace, R. V . E. & Hohlfeld, R. G. 2013, MNRAS, 429, 529

2013
[32]

Lovelace, R. V . E., Li, H., Colgate, S. A., & Nelson, A. F. 1999, ApJ, 513, 805

1999
[33]

& Pringle, J

Lynden-Bell, D. & Pringle, J. E. 1974, MNRAS, 168, 603

1974
[34]

2000, A&AS, 141, 165

Masset, F. 2000, A&AS, 141, 165

2000
[35]

& Bate, M

Meru, F. & Bate, M. R. 2012, MNRAS, 427, 2022 Müller, T. W. A., Kley, W., & Meru, F. 2012, A&A, 541, A123

2012
[36]

2026, MNRAS, 546, stag043

Nayakshin, S., Zhang, L., ´Calovi´c, A., et al. 2026, MNRAS, 546, stag043

2026
[37]

2025, ApJ, 995, 96

Ni, Y ., Deng, H., & Bai, X.-N. 2025, ApJ, 995, 96

2025
[38]

2012, MNRAS, 421, 3286

Paardekooper, S.-J. 2012, MNRAS, 421, 3286

2012
[39]

2011, MNRAS, 416, L65 Regály, Z

Paardekooper, S.-J., Baruteau, C., & Meru, F. 2011, MNRAS, 416, L65 Regály, Z. & V orobyov, E. 2017, MNRAS, 471, 2204 Rendon Restrepo, S. & Barge, P. 2022, A&A, 666, A92 Rendon Restrepo, S., Rometsch, T., Ziegler, U., & Gressel, O. 2025, arXiv e- prints, arXiv:2506.10812

work page arXiv 2011
[40]

2015, MNRAS, 454, 1940

Rice, K., Lopez, E., Forgan, D., & Biller, B. 2015, MNRAS, 454, 1940

2015
[41]

Rice, W. K. M., Armitage, P. J., Bate, M. R., & Bonnell, I. A. 2003, MNRAS, 339, 1025

2003
[42]

Rice, W. K. M., Forgan, D. H., & Armitage, P. J. 2012, MNRAS, 420, 1640

2012
[43]

Rice, W. K. M., Lodato, G., Pringle, J. E., Armitage, P. J., & Bonnell, I. A. 2004, MNRAS, 355, 543

2004
[44]

Rice, W. K. M., Lodato, G., Pringle, J. E., Armitage, P. J., & Bonnell, I. A. 2006, MNRAS, 372, L9

2006
[45]

Rice, W. K. M., Paardekooper, S.-J., Forgan, D. H., & Armitage, P. J. 2014, MNRAS, 438, 1593

2014
[46]

M., Moldenhauer, T

Rometsch, T., Jordan, L. M., Moldenhauer, T. W., et al. 2024, A&A, 684, A192

2024
[47]

Safronov, V . S. 1960, Annales d’Astrophysique, 23, 979

1960
[48]

& Whitworth, A

Stamatellos, D. & Whitworth, A. P. 2009, MNRAS, 392, 413

2009
[49]

Z., Tsukamoto, Y ., & Inutsuka, S

Takahashi, S. Z., Tsukamoto, Y ., & Inutsuka, S. 2016, MNRAS, 458, 3597

2016
[50]

1964, ApJ, 139, 1217

Toomre, A. 1964, ApJ, 139, 1217

1964
[51]

K., Klein, R

Truelove, J. K., Klein, R. I., McKee, C. F., et al. 1997, ApJ, 489, L179

1997
[52]

Tscharnuter, W. M. & Winkler, K.-H. A. 1979, Computer Physics Communica- tions, 18, 171 van der Walt, S., Schönberger, J. L., Nunez-Iglesias, J., et al. 2014, PeerJ, 2, e453 V on Neumann, J. & Richtmyer, R. D. 1950, Journal of Applied Physics, 21, 232 V orobyov, E. I. & Elbakyan, V . G. 2018, A&A, 618, A7

1979
[53]

1992, Earth Moon and Planets, 56, 173

Willerding, E. 1992, Earth Moon and Planets, 56, 173

1992
[54]

Young, M. D. & Clarke, C. J. 2015, MNRAS, 451, 3987

2015
[55]

Zhang, L., Nayakshin, S., Baruteau, C., Thebault, P., & V orobyov, E. I. 2026, arXiv e-prints, arXiv:2603.02395

work page arXiv 2026
[56]

& Baruteau, C

Zhu, Z. & Baruteau, C. 2016, MNRAS, 458, 3918

2016
[57]

& Springel, V

Zier, O. & Springel, V . 2023, MNRAS, 520, 3097 Besselϵ/H rms =0.6 Fig. 7.Detailed comparison of fragment evolution using the Bessel ker- nel (left) and the Plummer potential paradigm withϵ/Hrms =0.6 (right). The snapshots, from top to bottom, correspond tot=6.7,7.7,8.05,8.1 and 8.6 kyr. All plots are centered at the peak density. Article number, page 12 ...

2023