© The Authors 2025

1. Introduction

As of the writing of this paper, dense and narrow rings have been detected around three minor bodies: the centaur object Chariklo (Braga-Ribas et al. 2014) and two dwarf planets, Haumea (Ortiz et al. 2017) and Quaoar (Morgado et al. 2023). Although these objects exhibit significant diversity in size (roughly between 127 and 712 km) and rotational period (between 4 and 17 hours), the rings show surprising similarities. All are located very close to spin-orbit resonances (SORs) with the central mass, mainly the 1/3, and they show no evidence of any shepherding companions to explain the lack of radial diffusion¹. Even more baffling, the location of the main ring around Quaoar lies outside the Roche radius, raising the question as to why such a dense population of bodies did not accrete into a single satellite. The existence of a small outer moon (Weywot) close to a 6/1 mean-motion resonance (MMR) with the ring may increase the eccentricities and inhibit accretion; however, there is little concrete evidence of such a scenario.

The proximity to resonances is particularly challenging and contrasts with the dynamical characteristics of ring systems around giant planets. To date no plausible mechanism has been presented to explain how SORs could have shaped the structure of rings around minor bodies or how orbital evolution could have led to such configurations. Resonance capture requires convergent migration, while the outward spiraling of ring particles due to collisions and/or diffusive effects gives rise to a divergent migration. These difficulties have resulted in a certain consensus that perhaps there is no causal relation between the observed population and SORs.

In this Letter we explore the possibility that resonance capture in SORs could have occurred as a result of the combination of three effects: (i) gravitational perturbations from the non-spherical component of the central mass, (ii) collisions in the proto-ring, and (iii) secular changes in the rotational period of the minor body due to angular momentum conservation of the full (body plus ring) dynamical system. While the roles of the first two interactions have been analyzed in the past, the back-reaction onto the central mass due to the orbital evolution of the proto-ring has not been considered and will prove crucial. For adequate values of system parameters, we will show how these interactions, particularly the spin-down of the central mass, can lead to resonance capture in a 1/3 SOR of a significant part of an initial particle disk and the formation of a stable narrow ring similar to those observed. Our analysis, however, is more oriented toward presenting a proof of concept of this mechanism and is not aimed at reproducing any particular detail of the known ring systems around Centaurs, trans-neptunian objects (TNOs), or dwarf planets.

2. Numerical simulations

2.1. Physical scenario and dynamical interactions

Following previous studies (Sicardy et al. 2019; Madeira et al. 2022; Giuliatti Winter et al. 2023), we approximated the gravitational field of the irregularly shaped central body of total mass M_c by a spherical component of mass M₀, physical radius R₀, and spin frequency Ω₀, plus a point-mass anomaly M₁, fixed on its surface². We note that M_c = M₀ + M₁. Given the oblateness parameters estimated for Chariklo (ϵ ∼ 0.2), Haumea (ϵ ∼ 0.61), and Quaoar (ϵ ≳ 0.1), we expected the ratio μ = M₁/M₀ to be considerable, even if quantitative values are beyond our possibilities. The barycentric positions of each mass component, M_k, (with k = 0, 1) is denoted by ρ_k, with modulus $|{\mathit{\rho }}_0|=R_0{M}_1/M_c$ and $|{\mathit{\rho }}_1|=R_0{M}_0/M_c$, both with angular frequency Ω₀.

We assumed a particle disk orbited the central mass and was the precursor of the observed ring structure. Its formation and origin lies beyond the scope of this work, but we presume it resulted from a partially disruptive collision of the central mass with an external body or from the fragmentation of a primordial satellite (see, e.g., Melita et al. 2017). We assumed the disk to be coplanar with the equator of the massive body and to be radially bounded between r_min and r_max. Additionally, we supposed a certain surface density profile Σ_disk(r), and a total mass m_disk. Any initial eccentricity distribution proved immaterial, rapidly reaching equilibrium values (denoted as e_eq, see next sub-section) in a few orbital periods. The size distribution of the particles is relevant only when modeling their collisional evolution.

The equations of motion of each disk particle may be written, in a barycentric reference frame, as:

$$\begin{array}{c}\hfill \ddot{\mathbf{r}}={\ddot{\mathbf{r}}}_{\mathrm{grav}}+{\ddot{\mathbf{r}}}_{\mathrm{col}}=-{\displaystyle \sum _{k=0}^1}\mathcal{G}{M}_k\frac{(\mathbf{r}-{\rho }_k)}{|\mathbf{r}-{\rho }_k{|}^3}-\eta n\dot{r}\widehat{\mathbf{r}}\end{array}$$(1)

where $\mathcal{G}$ is the gravitational constant. The dynamical system is thus described by a circular restricted three-body problem, ${\ddot{\mathbf{r}}}_{\mathrm{grav}}$, with non-Keplerian motion of the primaries, plus an external radial friction term, ${\ddot{\mathbf{r}}}_{\mathrm{col}}$. We disregarded mutual gravitational interactions between elements of the disk, treating them as massless particles. The external force ${\ddot{\mathbf{r}}}_{\mathrm{col}}$ is a simple analytical prescription that mimics the dynamical effects of mutual collisions between ring elements (Sicardy et al. 2019). In this expression, n is the mean motion of the particle, and η is a dimensionless coefficient related (among other things) to the particle size and density profile at the location. We explored values between 10⁻⁴ and 10⁻².

Since the external force is radial, it will preserve the angular momentum of each particle and thus of the disk as a whole. While this simplistic collisional model cannot compete with fully interacting N-body simulations, it does appear to reproduce many of its general features (Salo et al. 2021).

2.2. Equilibrium eccentricity and radial confinement

Conservation of angular momentum by ${\ddot{\mathbf{r}}}_{\mathrm{col}}$ implies that a(1 − e²) = const, and the orbits suffer a damping (proportional to η) while exhibiting almost no change in the semimajor axis. Were it not for the mass anomaly, an eccentric initial disk would end as a circular structure analogous to the Saturn rings. However, the non-spherical component of the central mass introduces an equilibrium eccentricity, e_eq, around which the particles oscillate (see Appendix A). Under the effects of an exterior non-conservative force such as our radial term, this stable solution morphs into an attractor, and any initial eccentricity distribution tends toward these values, hence the name “equilibrium eccentricity”.

The top frame of Fig. 1 shows the orbital evolution of two elements of an extended disk. The initial conditions are shown with large filled circles, while arrows indicate the flow. Regardless of the initial orbit chosen above (green) or below (gray) e_eq, the non-conservative force will rapidly lead the eccentricity to oscillate around the equilibrium eccentricity highlighted in orange. Once this occurs, the disk continues to evolve outward (see Appendix A), albeit at a much slower rate.

Fig. 1. Orbital evolution of two disk particles under the effects of gravitational perturbations and collisional drag (η = 10−3). The central body was defined by μ = 0.05, R0 = 115 km and Ω0/n1 = 0.46. Top: Eccentricity versus semimajor axis. Initial conditions are indicated in filled circles, while arrows show the direction of flow. The equilibrium eccentricity eeq, is highlighted in orange. Bottom: Semimajor axes as a function of time. Notice that the semimajor axes of the particles converge over time.

Perhaps the most interesting result of this section is shown in the bottom frame of Fig. 1, where we compare the semimajor axes of the two disk particles over time. As a grows, e_eq is strongly reduced. Even though the eccentricity gradient also decreases, it is insufficient to counteract the previous trend, and the rate of change of the semimajor axis decreases. The resulting effect is an accumulation of the disk particles into an increasingly narrow ring. The time necessary for this radial confinement depends on several factors, including (i) the initial width and location of the disk; (ii) the collision coefficient, η; (iii) the mass ratio, μ; and (iv) the spin rate of the central body, Ω₀. All other things being equal, this mechanism indicates that slow rotators are more efficient at converting disks into rings (Appendix A).

2.3. Rotational spin-down and resonance capture

Provided the confinement mechanism described above is more effective than any radial diffusion, the end product will be a stable narrow ring without the need of any shepherding satellite. As shown in the bottom frame of Fig. 1, the ring will continue to recede from the central body, albeit at an increasingly slowler rate. The eccentricity of the ring particles oscillate around e_eq, taking values of the order of e ∼ 10⁻⁴ − 10⁻², depending on the system parameters. The dynamics, however, are completely secular, and so far there is no reference to an SOR. The fact that all observed rings around minor planets favor SORs, particularly the 1/3, indicates that other dynamical phenomena must be at play.

For resonance capture to occur, orbital migration must be convergent (e.g., Beaugé & Cerioni 2022). In the case of exterior spin-orbit resonances, this implies that the frequency ratio n/Ω₀ must increase over time. Since the orbital dynamics described in the previous section leads to an increase of the semimajor axis of the disk particles and consequently a reduction in the orbital frequency n, convergent migration can only occur if the spin rate, Ω₀, were to also decrease at an even faster rate.

Although a secular variation in the body’s spin seems unlikely, it is nonetheless expected when assuming the conservation of the total angular momentum of the system. More importantly, its magnitude may be estimated from the back-reaction of the disk migration acting on the rotation of M_c.

Even if self-gravity was not considered in our numerical simulations, the disk is expected to have a certain mass m_disk, following a surface density profile, Σ_disk(r). By representing this population by a set of N particles, each of mass m_j, semimajor axis a_j, and eccentricity e_j, the orbital angular momentum of the disk may then be expressed as:

$$\begin{array}{c}\hfill L_{\mathrm{disk}}={\displaystyle \sum _{j=1}^N}{m}_j|{\mathbf{r}}_{\mathbf{j}}\times {\dot{\mathbf{r}}}_{\mathbf{j}}|={\displaystyle \sum _{j=1}^N}{m}_j\sqrt{\mathcal{G}{M}_c{a}_j(1-e_j^2)}.\end{array}$$(2)

Assuming that the disk retains radial symmetry throughout its orbital evolution, the center of mass of the complete system will coincide with that of the central mass, M_c. Consequently, the total angular momentum may be written in barycentric coordinates as:

$$\begin{array}{c}\hfill L_{\mathrm{tot}}=\mathrm{\Lambda }{\mathrm{\Omega }}_0+L_{\mathrm{disk}};\mathrm{where}\mathrm{\Lambda }=\frac{2}{5}{M}_0{R}_0^2+{\displaystyle \sum _{k=0}^1}{M}_k{\rho }_k^2,\end{array}$$(3)

with Λ = Λ(M₀, M₁, R₀), a function of the physical characteristics of the central body and independent of its rotation. The first term stems from the rotational motion of M₀, while the second is the contribution from the translational motion of the central dipole. Invariance of L_tot implies that any increase in the disk angular momentum due to outward migration will lead to a decrease in the rotational frequency of the central body.

Lacking an explicit expression for e_eq(a), our analytical insights could not go any further. However, we could perform numerical simulations introducing ad hoc spin-down rates ${\dot{\mathrm{\Omega }}}_0$, and study its effect on the disk dynamics. Fig. 2 presents three examples plotting the semimajor axes of disk particles as a function of time. We chose a Chariklo-like body for a central mass, with M_c = 6.3 × 10¹⁸ kg, a reference radius of R₀ = 115 km, and a rotational period P_rot ≡ 2π/Ω₀ = 7.004 hs (Leiva et al. 2017; Sicardy et al. 2019). The non-sphericity was modeled by a single mass anomaly with μ = 0.05. The initial disks were assumed to be circular with radial symmetry, between r_min = 2.5R₀ (slightly above the co-rotation radius) and r_max = 3.8R₀.

Fig. 2. Three sets of simulations of the dynamical evolution of particle disks considering different values for the collisional coefficient, η. In all cases, the mass anomaly was taken to be equal to μ = M1/M0 = 0.05, while a spin-down rate with e-folding time τs = 105 years was adopted for the rotational frequency. The blue curves show the semimajor axes associated with the 1/3 and 1/2 SORs. The vertical lines indicate the times at which the spin-down reached values of Ω0 equal to 90% and 80% of the initial rotation frequency.

Each disk, consisting of 35 massless particles under the effects of gravitational perturbations and radial drag, was integrated for 10⁵ years. We also included an ad hoc spin-down of the central mass with a characteristic timescale of ${\tau }_s=-{\mathrm{\Omega }}_0/{\dot{\mathrm{\Omega }}}_0$. The simulations shown in Fig. 2 used τ_s = 10⁵ years. Since the equilibrium eccentricity is inversely proportional to (Ω₀ − n), the spin-down helped fuel the outward migration of the disk and greatly reduced the evolution timescale. More importantly, variations in Ω₀ moved the position of the SORs outward.

Each frame shows results for a different value of η. The inner initial disk, corresponding to orbits with a semimajor axis of a ≲ 3.0R₀, has large values of equilibrium eccentricities e_eq, and thus suffers a very fast radial migration, reaching separations of the order of 3R₀ in a few decades. After this first stage, evolution continues more slowly. For η = 10⁻² (top panel), the initial disk converges to a narrow ring in timescales on the order of T ∼ 10³ years. However, the outward evolution of the ring proved faster than that of the SORs, leading to divergent migration and no resonance capture.

For η = 10⁻³ (center panel of Fig. 2) the outcome is very different. Even if the radial confinement takes longer (T ∼ 10⁴ years), resonance capture is observed early in the process for initial conditions above and below the location of the 1/3 SOR. The resonant dynamics is preserved throughout the evolution of the narrow ring, with eccentricities displaying short-period variations between 10⁻³ and 10⁻². The width of the ring continues to decrease as it recedes from the central body and M₀ spins down, reaching values analogous to those currently observed.

Finally, in the case of η = 10⁻⁴ (bottom panel of Fig. 2) we observed the formation of two rings. The first ring is narrow, captured again in the 1/3 SOR. It formed a bit earlier than in the previous example and consists of particles initially exterior to the resonance location. A second broader ring is formed closer to the central body and is ultimately trapped in the 1/2 SOR. However, this requires a spin-down of more than 50% of the initial value, and may thus be considered an unlikely result.

These simulations are illustrative and have been chosen to highlight the different possible outcomes. While the mechanism leading to radial confinement is extremely robust and was found in all simulations, resonance capture is more sensitive. The SORs require a minimum width for trapping to take place, which depends on the value of μ. We found no resonant ring for μ < 10⁻³, although the limit is a function of η. Similarly, convergent migration requires minimum values for the spin-down rates, and such a condition was not met for τ_s ≳ 10⁶ years. Within these ranges, outcomes analogous to the bottom frame of Fig. 2 were more frequent than that depicted in the center panel. However, our exploration of the parameter space is far from thorough.

Since our simulations are not self-consistent, we had to estimate the disk mass necessary to generate the assumed spin-down rate used for Fig. 2. Once again, we invoked the conservation of total angular momentum. Considering disk elements of equal mass, we rewrote expression (3) as:

$$\begin{array}{c}\hfill L_{\mathrm{tot}}=\mathrm{\Lambda }{\mathrm{\Omega }}_0+m_{\mathrm{disk}}{\widehat{L}}_{\mathrm{disk}}\end{array}$$(4)

where ${\widehat{L}}_{\mathrm{disk}}$ is the disk angular momentum per unit mass, a value that can be calculated from the numerical simulations at any given time. By calculating this quantity at two different times, say at t = 0 and when the spin was reduced by 20%, we could calculate the values of m_disk that preserves L_tot.

The results from Fig. 2 yield values for the total disk mass in the range of m_disk ∼ 8 × 10⁻³M₀, showing only a weak dependence with η. No significant changes were observed when varying the time at which the disk angular momentum is evaluated, indicating a fairly robust estimation. Although these values are much larger than the current estimations (Pan & Wu 2016), collisional grinding within the dense ring could have pulverized a large portion of the original mass later ejected by radiation pressure (Morgado et al. 2023). This process could also aid in reducing the spin-down and lead to a steady-state configuration.

3. Conclusions

In this work, we have proposed a dynamical model for the origin and evolution of dense narrow rings around minor bodies, starting from an extended debris disk and evolving toward dynamical structures similar to those observed in several Centaurs and TNOs. Based on first principles, the observed radial confinement is the natural outcome of the gravitational interactions with the central body and a radial external force. Although we assumed a radial Stokes-like acceleration, analogous results should be obtained for any non-conservative force resulting in strong eccentricity damping.

The back-reaction on the central body’s spin is also a crucial ingredient, leading to convergent migration and (at least for some disk particles) resonance capture. Coincidental with observed systems, the 1/3 SOR appears as the most preferable commensurability. The radial confinement is maintained by resonant dynamics and the eccentricity gradient, without the need of shepherding satellite companions. The eccentricities remain low (∼10⁻³ − 10⁻²) but are large enough to inhibit accretional collisions (e.g., Brilliantov et al. 2015), and thus perhaps explain the existence of rings beyond the Roche radius, as observed around Quaoar. Since no apsidal alignment was detected, the rings would be observed as circular.

These promising results require debris disks that are more massive than estimated for the observed rings. While we cannot rule out the limitations of our numerical method and models for particle interactions (collisions and self-gravity), this discrepancy could be partially solved by collisional grinding over long timescales, which are particularly active once the localized surface density increases as a consequence of radial confinement.

Acknowledgments

Simulations were carried out with the computing facilities hosted at IATE as well as in the High Performance Computing Center of the Universidad Nacional de Córdoba (CCAD-UNC). This research was funded by CONICET and Secyt/UNC.

References

Appendix A: The semi-secular and equilibrium eccentricities

In addition to the system parameters defined in 2.1, let θ₁ denote the azimuthal angle of the position of M₁ in a given inertial reference frame. In the two-body problem, where the motion of the perturber M₁ with respect to the center of mass is keplerian, i.e. ${\dot{\theta }}_1=n_1\equiv {(\mathcal{G}{M}_c/R_0^3)}^{1/2}$, the equilibrium eccentricity e_eq lies close to the semi-secular forced eccentricity e_ff, associated to stable fixed points in the restricted three-body problem, averaged over the orbital period of the perturber. For initial conditions exterior to the binary M₀ + M₁, an octupole-level expansion of the disturbing functions yields:

$$\begin{array}{c}\hfill e_{\mathrm{ff}}=\frac{3}{4}\frac{\mu }{{(1+\mu )}^2}(\frac{R_0}{a}{)}^2+\frac{45}{64}\frac{\mu (1+{\mu }^3)}{{(1+\mu )}^5}(\frac{R_0}{a}{)}^4.\end{array}$$(A.1)

A simplified version, limited to the quadrupole level, was found by Paardekooper et al. (2012). However, given the proximity of the future ring to M_c, an extension to octupole terms proves necessary. It is important to note that e_ff is equivalent to the classical secular forced eccentricity for exterior orbits (a > R₀) in the case where the perturber lies in a circular orbit.

The dashed black line in Fig. A.1 plots e_ff as function of a/R₀, assuming a mass ratio μ = 0.05. The exact values of the equilibrium eccentricity, in the case of M₁ in keplerian orbit, are shown in red, and were obtained from N-body simulations whose outputs were averaged over timescales much larger than orbital period. As expected, e_eq ≃ e_ff, and the discrepancy tends to zero for increasing semimajor axis. Although the difference could be due to higher order terms in the Legendre expansion of the disturbing function, its origin lies elsewhere and is much more interesting.

Fig. A.1. Schematic view of the different interactions acting on the initial particle disk. Curves of constant angular momentum (per unit mass) are shown in gray, with arrows indicating evolutionary routes due to a radial drag term. Red and orange curves highlight the equilibrium eccentricities of initial conditions undergoing both collisions and gravitational perturbations (assuming μ = 0.05 and R0 = 115 km). Both differ in the assumed angular frequency ${\dot{\theta }}_1$ of M1 around the center of mass of Mc. The red curve shows results assuming ${\dot{\theta }}_1$ equal to the keplerian mean-motion (n1). Conversely, the orange curve was calculated considering a sub-keplerian orbital frequency equal to the spin rate of M0, i.e. ${\dot{\theta }}_1={\mathrm{\Omega }}_0=0.46n_1$. For comparison, the black dashed curve shows the semi-secular equilibrium eccentricity eff, as given by equation (A.1).

We return to our minor body M_c = M₀ + M₁, where M₁ co-rotates with M₀. In such a case, the motion of the perturber with respect to the center of mass is sub-keplerian, with a corresponding angular ("orbital") frequency equal to Ω₀ < n₁. The orange curve in Fig. A.1 shows the equilibrium eccentricity in this case, considering Ω₀ = 0.46n₁, a value similar to the current spin of Chariklo. The values of this new equilibrium eccentricity are significantly higher than those obtained for ${\dot{\theta }}_1=n_1$ and the semi-secular limit e_ff.

Numerical simulations, backed by perturbation theory, indicate that the change in e_eq is related to the magnitude of short-period variations, which, by construction, are not considered in the secular or semi-secular disturbing functions. While a detailed model for the sub-keplerian equilibrium eccentricity is outside the scope of this Letter, its main features may be described by an expression of the type: e_eq² − e_ff² ∝ n/(k₁Ω₀ − k₂n)², with k₁ and k₂ integers. Thus, slow-rotating central bodies generate larger equilibrium eccentricities than fast rotators. As the orbital frequency of M₁ increases, the difference with e_ff decreases, as observed comparing the red and dashed black curves.

For low eccentricities (say ≲10⁻²), the curves of constant angular momentum are practically vertical en in (a, e) plane and the dominant effect is focused on the eccentricity. Thus, when the disk particle reaches the equilibrium eccentricity, the non-conservative term has a non-zero component along the gradient of e_eq(a). This eccentricity gradient then generates a secular change in the semimajor axis toward larger values. The rate of change in a may be written as:

$$\begin{array}{c}\hfill \frac{\mathit{da}}{\mathit{dt}}\simeq (\frac{\partial e_{\mathrm{eq}}}{\partial a}{)}^{-1}{\frac{\mathit{de}}{\mathit{dt}}|}_{\mathrm{col}}\simeq -\frac{1}{2}\eta n{e}_{\mathrm{eq}}(\frac{\partial e_{\mathrm{eq}}}{\partial a}{)}^{-1}.\end{array}$$(A.2)

Since the equilibrium eccentricity gradient is always negative, all disk particles recede from the central body over time.

All Figures

Fig. 1. Orbital evolution of two disk particles under the effects of gravitational perturbations and collisional drag (η = 10−3). The central body was defined by μ = 0.05, R0 = 115 km and Ω0/n1 = 0.46. Top: Eccentricity versus semimajor axis. Initial conditions are indicated in filled circles, while arrows show the direction of flow. The equilibrium eccentricity eeq, is highlighted in orange. Bottom: Semimajor axes as a function of time. Notice that the semimajor axes of the particles converge over time.

In the text

Fig. 2. Three sets of simulations of the dynamical evolution of particle disks considering different values for the collisional coefficient, η. In all cases, the mass anomaly was taken to be equal to μ = M1/M0 = 0.05, while a spin-down rate with e-folding time τs = 105 years was adopted for the rotational frequency. The blue curves show the semimajor axes associated with the 1/3 and 1/2 SORs. The vertical lines indicate the times at which the spin-down reached values of Ω0 equal to 90% and 80% of the initial rotation frequency.

In the text

Fig. A.1. Schematic view of the different interactions acting on the initial particle disk. Curves of constant angular momentum (per unit mass) are shown in gray, with arrows indicating evolutionary routes due to a radial drag term. Red and orange curves highlight the equilibrium eccentricities of initial conditions undergoing both collisions and gravitational perturbations (assuming μ = 0.05 and R0 = 115 km). Both differ in the assumed angular frequency ${\dot{\theta }}_1$ of M1 around the center of mass of Mc. The red curve shows results assuming ${\dot{\theta }}_1$ equal to the keplerian mean-motion (n1). Conversely, the orange curve was calculated considering a sub-keplerian orbital frequency equal to the spin rate of M0, i.e. ${\dot{\theta }}_1={\mathrm{\Omega }}_0=0.46n_1$. For comparison, the black dashed curve shows the semi-secular equilibrium eccentricity eff, as given by equation (A.1).

In the text