# Primordial Origin of Methane on Eris and Makemake
## Supported by D/H Ratios

Olivier Mousis, Aaron Werlen, Tom Benest Couzinou, Antoine Schneeberger

► To cite this version:

Olivier Mousis, Aaron Werlen, Tom Benest Couzinou, Antoine Schneeberger. Primordial Origin of Methane on Eris and Makemake Supported by D/H Ratios. The Astrophysical Journal Letters, 2025, 983 (1), pp.L12. 10.3847/2041-8213/adc134. hal-05486747

HAL Id: hal-05486747

https://hal.science/hal-05486747v1

Submitted on 23 Feb 2026

HAL is a multi-disciplinary open access
archive for the deposit and dissemination of sci-
entific research documents, whether they are pub-
lished or not. The documents may come from
teaching and research institutions in France or
abroad, or from public or private research centers.

L'archive ouverte pluridisciplinaire HAL, est
destinée au dépôt et à la diffusion de documents
scientifiques de niveau recherche, publiés ou non,
émanant des établissements d'enseignement et de
recherche français ou étrangers, des laboratoires
publics ou privés.

<figure><img src="image_1.png" /></figure>

Distributed under a Creative Commons CC BY 4.0 - Attribution - International License
OPEN ACCESS

<figure><img src="image_2.png" /></figure>

CrossMark

# Primordial Origin of Methane on Eris and Makemake Supported by D/H Ratios

Olivier Mousis<fnref n="1" />,<fnref n="2" /> ID, Aaron Werlen<fnref n="1" />,<fnref n="3" /> ID, Tom Benest Couzinou<fnref n="1" /> ID, and Antoine Schneeberger<fnref n="1" /> ID

Received 2024 December 30; revised 2025 March 16; accepted 2025 March 16; published 2025 April 4

## Abstract

Deuterium, a heavy isotope of hydrogen, is a key tracer of the formation of the solar system. Recent James Webb Space Telescope observations have expanded the data set of deuterium-to-hydrogen (D/H) ratios in methane on the KBOs Eris and Makemake, providing new insights into their origins. This study examines the elevated D/H ratios in methane on these KBOs in the context of protosolar nebula (PSN) dynamics and chemistry, proposing a primordial origin for the methane, in contrast to previous hypotheses suggesting abiotic production by internal heating. A time-dependent disk model coupled with a deuterium chemistry module was used to simulate the isotopic exchange between methane and hydrogen. Observational constraints, including the D/H ratio measured in methane in comet 67P/Churyumov-Gerasimenko, were used to refine the primordial D/H abundance. The simulations show that the observed D/H ratios in methane on Eris and Makemake are consistent with a primordial origin. The results suggest that the methane on these Kuiper Belt Objects likely originated from the PSN, similar to cometary methane, and was sequestered in solid form—either as pure condensates or clathrates—within its building blocks prior to accretion. These results provide a simple explanation for the high D/H ratios in methane on Eris and Makemake, without the need to invoke internal production mechanisms.

*Unified Astronomy Thesaurus concepts: Trans-Neptunian objects (1705); Cosmochemistry (331); Astrobiology (74); Protoplanetary disks (1300); Solar system formation (1530)*

## 1. Introduction

Deuterium, a heavy isotope of hydrogen, was produced during the birth of the Universe with a primordial deuterium-to-hydrogen (D/H) abundance of ~2.5 × 10⁻⁵ (P. A. Kislitsyn et al. 2024) and is subsequently destroyed in stellar interiors. As a key tracer of the formation and evolution of the solar system, the D/H ratio varies with temperature and chemical conditions, providing critical insights into the thermodynamic state of the early protosolar nebula (PSN) and the origin of planetary bodies. In the cold outer regions of the PSN, deuterium fractionation produces deuterium-enriched ices and organics, typically close to values observed in the interstellar medium (ISM), while the warmer inner regions acquire lower D/H ratios. This creates a spatial and thermal record of the protoplanetary disk, preserved in the D/H ratios measured across different reservoirs of the solar system (A. Drouart et al. 1999; O. Mousis et al. 2000; C. Ceccarelli et al. 2014).

Recent observations with the James Webb Space Telescope (JWST) have revealed the presence of monodeuterated methane (CH₃D) on the surfaces of the Kuiper Belt Objects (KBOs) Eris and Makemake. The measured D/H ratios in methane were (2.5 ± 0.5) × 10⁻⁴ for Eris and (2.9 ± 0.6) ± 10⁻⁴ for Makemake (W. M. Grundy et al. 2024). Previously, D/H measurements in CH₄ were limited to Titan's atmosphere, where the Cassini Composite Infrared Spectrometer determined a value of (1.59 ± 0.33) × 10⁻⁴ (C. A. Nixon et al. 2012) and a small number of comets. For comets C/2004 Q2 (Machholz), C/2001 Q4 (NEAT), and C/2007 N3 (Lulin), only upper limits were reported, with D/H values of <0.005, <0.01, and <0.0075,

respectively (H. Kawakita et al. 2005; B. P. Bonev et al. 2009; E. L. Gibb et al. 2012). In contrast, a precise D/H measurement of (2.41 ± 0.29) × 10⁻<fnref n="3" /> was obtained for comet 67P/Churyumov-Gerasimenko (67P/C-G), using the ROSINA instrument on board the Rosetta mission (D. R. Müller et al. 2022).

The D/H measurements in CH₄ on Makemake and Eris, as interpreted by C. R. Glein et al. (2024), provide evidence for abiotic methane production. According to this interpretation, CH₄ production may have been triggered by significant radiogenic heating, potentially driving hydrothermal circulation at the base of an ice-covered ocean and facilitating geochemical transformations of volatile species. C. R. Glein et al. (2024) also suggest that metamorphic reactions involving accreted organic matter may have occurred due to heating in the deeper interior, resulting in the production of thermogenic methane.

Here, we propose an alternative interpretation, suggesting that the D/H ratios measured in methane on Makemake and Eris align with a primordial origin, similar to the CH₄ observed in numerous comets. In this scenario, the methane on these KBOs would have originated in the PSN and been accreted during their formation, similar to the process hypothesized for methane in comets. To investigate this hypothesis, we use a time-dependent protoplanetary disk model coupled with a module that simulates the deuterium exchange between CH₄ and H₂ in the gas phase of the PSN. A key uncertainty in this analysis is the primordial deuterium abundance acquired by CH₄ from the ISM during its incorporation into the PSN, since no reliable measurements of the D/H ratio in interstellar methane exist. To date, the only tentative ISM detection potentially associated with CH₃D does not precisely match any cataloged spectral lines (N. Sakai et al. 2012). Consequently, we adopt the D/H ratio measured in comet 67P/C-G as a conservative estimate for our calculations.

<figure><img src="image_3.png" /></figure>

Original content from this work may be used under the terms of the [Creative Commons Attribution 4.0 licence](https://creativecommons.org/licenses/by/4.0/). Any further distribution of this work must maintain attribution to the author(s) and the title of the work, journal citation and DOI.

<footnote n="1">¹Aix-Marseille Université, CNRS, CNES, Institut Origines, LAM, Marseille, France; olivier.mousis@lam.fr</footnote>

<footnote n="2">²Institut Universitaire de France (IUF), France</footnote>

<footnote n="3">³ETH Zurich, Institute for Particle Physics and Astrophysics, Wolfgang-Pauli-Strasse 27, 8093, Zurich, Switzerland</footnote>
Section 2 provides a detailed description of the PSN model and its accompanying deuterium chemistry module. Section 3 presents the simulation results, while Section 4 outlines the conclusions and discusses their implications.

## 2. Model

### 2.1. Model of PSN

To compute the thermodynamic evolution of the PSN, we use a one-dimensional $\alpha$-viscous accretion-disk model based on the work of A. Aguichine et al. (2020). The viscosity of the disk is given by $\nu = \alpha c_s^2 / \Omega_K$, following the prescription of N. I. Shakura & R. A. Sunyaev (1973), where $\alpha$ is the viscosity parameter. Here, $c_s = \sqrt{R_g T / \mu}$ is the isothermal sound velocity, where $T$ is the midplane temperature of the disk, $R_g$ is the ideal-gas constant, and $\mu = 2.31 \text{ g mol}^{-1}$ is the mean molar mass of the gas (K. Lodders 2019). The Keplerian frequency $\Omega_K = \sqrt{GM_* / r^3}$ depends on the gravitational constant $G$, the mass of the star $M_*$, and the distance from the star $r$. The evolution of the disk is governed by changes in the surface density of the gas $\Sigma_g$ and its midplane temperature $T$. The surface density of the gas is calculated using the equation from D. Lynden-Bell & J. E. Pringle (1974):

$$ \frac{\partial \Sigma_g}{\partial t} = \frac{3}{r} \frac{\partial}{\partial r} \left[ r^{\frac{1}{2}} \frac{\partial}{\partial r} (r^{\frac{1}{2}} \Sigma_g \nu) \right], \quad (1) $$

which we rewrite as

$$ \begin{gathered} \frac{\partial \Sigma_g}{\partial t} = \frac{1}{2\pi r} \frac{\partial \dot{M}}{\partial r} \\ \dot{M} = 3\pi \Sigma_g \nu \left( 1 + 2 \frac{d \ln(\Sigma_g \nu)}{d \ln r} \right), \end{gathered} \quad (2) $$

where $\dot{M}$ is the mass accretion rate.

We use the following self-similar solutions as initial conditions to solve these equations (D. Lynden-Bell & J. E. Pringle 1974; A. Aguichine et al. 2020):

$$ \left\{ \begin{aligned} \Sigma_{g,0} &= \frac{M_{\text{acc},0}}{3\pi\nu} e^{-(r/c)^{0.5}} \\ \dot{M}_0 &= \dot{M}_{\text{acc},0} \left(1 - \left(\frac{r}{r_c}\right)^{0.5}\right) e^{-(r/c)^{0.5}}, \end{aligned} \right. \quad (3) $$

with $r_c$ being the centrifugal radius and $\dot{M}_{\text{acc},0}$ being the initial mass accretion rate.

The midplane temperature is computed by summing all the heating sources (R. Hueso & T. Guillot 2005):

$$ T^4 = \frac{1}{2\sigma_{\text{sb}}} \left( \frac{3}{8}\tau_R + \frac{1}{2\tau_P} \right) \Sigma_g \nu \Omega_K^2 + T_{\text{amb}}^4, \quad (4) $$

The first term represents the viscous heating contribution, where $\sigma_{\text{sb}}$ is the Boltzmann constant and $\tau_R$ and $\tau_P$ are the Rosseland and Planck optical depths, respectively. These are related by $\tau_P = 2.4\tau_R$ and $\tau_R = \frac{\kappa_R \Sigma_g}{2}$, where $\kappa_R = \kappa_0 \rho^a T^b$ is the mean Rosseland opacity, with $\rho$ being the midplane density and $\kappa_0$, $a$, and $b$ being determined from empirical data obtained from Table 3 of K. R. Bell & D. N. C. Lin (1994). The last term in Equation (4) accounts for the contribution of the interstellar ambient temperature $T_{\text{amb}} = 10 \text{ K}$. In the temperature calculation, irradiation from the central star is neglected, due to dust accumulation at the ice-line positions (A. Aguichine et al. 2020; K. Ohno & T. Ueda 2021; A. Schneeberger et al. 2023).

<figure><img src="image_4.png" /><figcaption>**Figure 1.** Temperature ($T$), pressure ($P$), and surface density ($\Sigma_g$) profiles of the gas in the PSN midplane, shown for various stages of its evolution.</figcaption></figure>

The $\alpha$ parameter is typically set within the range of $10^{-4}-10^{-2}$ (L. Hartmann et al. 1998; R. Hueso & T. Guillot 2005; S. J. Desch et al. 2017), while the disk mass accretion rate, $\dot{M}_{\text{acc}}$, is constrained by observations to values between $10^{-9}$ and $10^{-6} M_{\odot} \text{ yr}^{-1}$ (L. Hartmann et al. 1998; E. Gullbring et al. 1998). Based on these studies, we define our nominal PSN model using $\alpha = 10^{-3}$ and an initial accretion rate of $\dot{M}_{\text{acc},0} = 5 \times 10^{-7} M_{\odot} \text{ yr}^{-1}$ (L. Hartmann et al. 1998). Additionally, we assume an initial disk mass of $0.1 M_{\odot}$, primarily distributed within a radius of approximately 150 au (A. Aguichine et al. 2020; A. Schneeberger et al. 2023; O. Mousis et al. 2024). The midplane temperature, midplane pressure, and surface density profiles of the disk at various evolutionary stages are shown in Figure 1.

### 2.2. Deuterium Chemistry Module

The enrichment factor, $f$, is defined as the ratio of the D/H ratio in the deuterated species of interest to the D/H ratio in molecular hydrogen. In the case of methane, this ratio can be
expressed as

$$f = \frac{(\mathrm{D}/\mathrm{H})_{\mathrm{CH}_4}}{(\mathrm{D}/\mathrm{H})_{\mathrm{H}_2}} = \frac{\frac{\mathrm{CH}_3\mathrm{D}}{\mathrm{CH}_4}}{\frac{\mathrm{HD}}{\mathrm{H}_2}}. \qquad (5)$$

In the following, the D/H ratio in molecular hydrogen is set to the protosolar value, defined as $2.1 \times 10^{-5}$ (L. J. Hallis 2017). The evolution of $f$ in the PSN is governed by the following differential equation (A. Drouart et al. 1999; O. Mousis et al. 2000, 2002b):

$$\frac{\partial f}{\partial r} = k(T)P(A(T) - f) + \frac{1}{\Sigma r} \frac{\partial}{\partial r} \left( \kappa r \Sigma \frac{\partial f}{\partial r} \right). \qquad (6)$$

The first term on the right-hand side of Equation (6) represents the isotopic exchange between HD and the deuterated molecule of interest, CH₃D. The function $k(T)$ denotes the rate of isotopic exchange, $P$ is the total pressure of the disk, and $A(T)$ represents the equilibrium fractionation factor, which depends solely on the temperature $T$. The definition of $f$ is more general than that of $A(T)$, as it does not assume that the system under consideration is in equilibrium. The value of $k(T)$ is obtained from laboratory experiments (C. Lécluse & F. Robert 1994; C. Lécluse et al. 1996), while $A(T)$ is derived from quantum mechanical calculations by P. Richet et al. (1977) and subsequently fitted using the analytical expression proposed by C. Lécluse et al. (1996). Equation (6) holds under the assumption that isotopic exchange takes place exclusively between neutral molecules. As shown in O. Mousis et al. (2002b), $k(T)$ decreases by 7 orders of magnitude as $T$ decreases from 1000 K to 200 K. This means that at temperatures below 200 K, isotope exchange is practically inhibited. At high temperatures, it is also negligible, because equilibrium is quickly reached, so $f$ becomes equal to $A(T)$ and the term $(f - A(T))$ converges to 0. The second term on the right-hand side of Equation (6) represents the turbulent diffusion within the PSN. This process depends on the local surface density, $\Sigma(R, t)$, and the diffusivity, $\kappa$, defined as the ratio of the turbulent viscosity, $\nu$, to the Prandtl number, $Pr$. The Prandtl number is expressed following the relation provided by A. Drouart et al. (1999):

$$Pr = -\frac{4}{3}r\frac{\partial_r\Sigma_g}{\Sigma_g}. \qquad (7)$$

The value of $Pr$ varies with the opacity of the PSN but remains below unity, which serves as its upper limit (B. Dubrulle & U. Frisch 1991).

The integration of Equation (6) is governed by the spatial and temporal boundary conditions. For the PSN, the spatial boundary conditions are defined as $\partial f/\partial r = 0$ at $r = 1$ au and at the radial distance $r$ corresponding to the condensation or clathration point of methane for any $t$.<fnref n="4" /> The temporal boundary condition is $f(r) = \text{constant}$ at $t = 0$, applicable to all $r$. This

assumption reflects the hypothesis that the interstellar ices that accreted onto the protoplanetary disk were uniformly enriched throughout the presolar cloud, irrespective of heliocentric distance. Deviations from this assumption would have a negligible impact on the evolution of $f$ across the PSN (O. Mousis et al. 2002b). Equation (6) applies to the isotopic exchange between CH₃D and H₂ in the vapor phase. In particular, deuterium exchange is inhibited at temperatures below 200 K. As a result, the evolution of the D/H ratio in the cooling PSN, following the period when the temperature was 200 K, was driven solely by turbulent mixing of the gas.

The initial value $f_0$, which corresponds to $f$ at $t = 0$, is poorly constrained. Observations reveal substantial deuterium enrichment in ISM ices (C. Ceccarelli et al. 2014; L. I. Cleeves et al. 2016; L. Piani et al. 2021). This enrichment is mainly attributed to significant deuterium fractionation in prestellar cores, driven by their extremely low temperatures and the high propensity of molecules to accrete onto icy mantles on dust grains. These processes provide a robust explanation for the elevated deuterium fractionation observed in various molecules within prestellar cores (E. Roueff et al. 2005; C. Ceccarelli et al. 2014). Several studies have proposed that interstellar ices vaporized during accretion onto the protoplanetary disk, at heliocentric distances that potentially extend to 30 au (P. Cassen & K. M. Chick 1997). Upon vaporization, CH₃D exchanged deuterium with H₂ in the PSN, reducing the deuterium fractionation factor $f$, because the equilibrium values at nebula temperatures are lower than those in the ISM. To date, the only available determination of a purported “pristine” D/H ratio in CH₄ does not stem from ISM ices but is instead based on measurements obtained by the ROSINA instrument on board the Rosetta spacecraft. These measurements, conducted on Comet 67P/C-G, provide a value of $(2.41 \pm 0.29) \times 10^{-3}$ (D. R. Müller et al. 2022). The central value of this measurement will be used as the default initial condition for our simulations, yielding $f_0 \sim 115$ relative to the protosolar value.

We consider that CH₄ either condensed at approximately 28 K in the PSN or was trapped as clathrate at around 55 K (O. Mousis et al. 2009a) during the cooling of the PSN. Consequently, the D/H ratio in CH₄ within microscopic grains embedded in the nebula, and later in planetesimals, reflects the ratio acquired during the condensation or clathration of methane, depending on the crystallization scenario considered. In the following, our calculations consistently represent the D-enrichment value acquired by methane at the point of condensation or clathration in the PSN.

### 3. Results

Figures 2 and 3 show the evolution of the $f$-profile for CH₄ in the PSN, as determined from our nominal model. The profiles are calculated along the trajectories of the methane ice line and the methane clathrate line, corresponding to condensation and entrapment temperatures of 28 K and 55 K, respectively. These calculations are compared with the D/H measurements in CH₄ obtained by the Cassini spacecraft in the Titan atmosphere and by JWST on the surfaces of the KBOs Makemake and Eris (see Table 1). We consider three initial conditions, corresponding to 0.5, one, and two times the central enrichment measured on Comet 67P/C-G. The deuterium enrichment profiles show a sharp gradient between the 11–12 au and 8–9 au regions, before leveling off in the inner regions of the PSN. This behavior is attributed to the slow temporal evolution of the disk’s

<footnote n="4">⁴ $\partial f/\partial r = 0$ at $r = 1$ au and the condensation/clathration radius of CH₄ implies that the mass flux is neglected at these boundaries due to the poorly constrained physical conditions. This assumption is justified by the rapid decline of the mass flux, which is significant only during the earliest stages of disk evolution. At these early times, high temperatures drive chemical fractionation rather than diffusion, rapidly bringing isotopic enrichment to unity in the inner disk ($T > 1000$ K), while maintaining its maximum value in the colder outer regions ($T < 300$ K). At later stages, diffusion becomes the dominant process, but by then the mass flux is negligible. More details of the mathematical model can be found in A. Drouart et al. (1999).</footnote>
<figure><img src="image_5.png" /><figcaption>**Figure 2.** Deuterium enrichment profile *f* calculated along the trajectory of the methane ice line in the PSN, assuming a condensation temperature of 28 K. The thick horizontal lines indicate the nominal enrichment values observed for Makemake (gray), Eris (green), and Titan (red). The dashed lines in matching colors represent the $1\sigma$ uncertainties associated with these measurements. From left to right, the initial enrichment values were set at 0.5, one, and two times the central value of 115 measured in Comet 67P/C-G. The epochs where *f* intersects the central values observed for the three bodies are highlighted, corresponding to a very narrow variation in the location of the methane ice line within the PSN.</figcaption></figure>

<figure><img src="image_6.png" /><figcaption>**Figure 3.** Similar to Figure 2 but showing the deuterium enrichment profile *f* calculated along the trajectory of the methane clathration in the PSN, with an assumed entrapment temperature of 55 K.</figcaption></figure>

<table><caption>**Table 1**
D/H Measurements in CH$_4$ in 67P/C-G, Titan, Eris, and Makemake</caption><thead><tr><th>Object</th><th>Value</th><th>Reference</th></tr></thead><tbody><tr><td>67P/C-G</td><td>(2.41 ± 0.29) × 10<sup>-3</sup></td><td>D. R. Müller et al. (2022)</td></tr><tr><td>Titan</td><td>(1.59 ± 0.33) × 10<sup>-4</sup></td><td>C. A. Nixon et al. (2012)</td></tr><tr><td>Eris</td><td>(2.5 ± 0.5) × 10<sup>-4</sup></td><td>W. M. Grundy et al. (2024)</td></tr><tr><td>Makemake</td><td>(2.9 ± 0.6) × 10<sup>-4</sup></td><td>W. M. Grundy et al. (2024)</td></tr></tbody></table>

thermodynamic profiles, which causes a gradual inward progression of the two considered ice lines. The slow progression allows sufficient time for efficient isotopic exchange to take place in the gas phase before methane either condenses or becomes trapped in clathrates.

The calculated PSN enrichment profiles agree consistently with all observed measurements. The central values measured in Makemake, Eris, and Titan intersect the *f*-profile at different time intervals. For Makemake, the intersections occur at $5.8 \times 10^4$, $9.1 \times 10^4$, and $1.3 \times 10^5$ yr. For Eris, they occur at $6.5 \times 10^4$, $1.0 \times 10^5$, and $1.5 \times 10^5$ yr. For Titan, the intersections occur at $9.0 \times 10^4$, $1.3 \times 10^5$, and $2.0 \times 10^5$ yr. These values at the methane ice line correspond to $f_0$ values of 0.5,

one, and two times the central value of 115 measured in Comet 67P/C-G. Our calculations assume that the three bodies that accreted icy solids crystallized in a narrow region of the PSN, typically spanning a few tenths of an astronomical unit. This region lies between 11–12 au and 8–9 au, depending on whether solid methane forms from pure condensate or clathrate. Our results support the hypothesis that the building blocks of giant planet systems and KBOs originated within a compact region of the PSN, in agreement with predictions from dynamical models (K. Tsiganis et al. 2005; D. Nesvorný 2018). This interpretation aligns with the D/H ratio measured in H$_2$O on Enceladus, a moon of Saturn located closer to the planet than Titan. Data from the Cassini–Huygens mission revealed a cometary D/H ratio ($\sim 2.9 \times 10^{-4}$), suggesting that the building blocks of Enceladus formed in a region of the PSN similar to that of many comets (O. Mousis et al. 2009b; J. H. Waite et al. 2009; R. N. Clark et al. 2019). This conclusion also extends to Titan, as it accreted in a colder region of Saturn’s circumplanetary disk from building blocks that formed earlier in the PSN (O. Mousis et al. 2002a).

Additional tests assuming up to five times the central enrichment measured in Comet 67P/C-G were also performed and yielded consistent results. This suggests that our model
effectively reproduces the D/H ratios observed in Titan's atmosphere and on the surfaces of Eris and Makemake, without requiring an internal source of CH₄ from hydrothermal or metamorphic processes (S. K. Atreya et al. 2006; O. Mousis et al. 2009a, 2009c; C. R. Glein et al. 2024).

## 4. Discussion and Conclusions

The deuterium enrichments in CH₄ observed on Titan, Makemake, and Eris can be matched by our model for any value of the initial accretion rate, $\dot{M}_{acc,0}$, in the range $5 \times 10^{-7}$ to $10^{-6} M_{\odot} yr^{-1}$, assuming $\alpha = 10^{-3}$. At constant $\alpha$, higher accretion rates result in a warmer protoplanetary disk, causing the deuterium enrichment profile to shift outward to higher heliocentric distances while maintaining a similar overall shape. For example, assuming $\alpha = 10^{-3}$ and an initial accretion rate of $10^{-6} M_{\odot} yr^{-1}$ instead of $5 \times 10^{-7} M_{\odot} yr^{-1}$, the *f*-value crossings for the central measurements of Titan, Makemake, and Eris would shift about 2 au farther from the Sun. On the other hand, an initial accretion rate lower than $5 \times 10^{-7}$ would result in colder temperatures in the outer regions of the PSN, preventing the *f*-profile from decreasing to the values observed in the three bodies.

The calculated *f*-profiles also show significant sensitivity to variations in the $\alpha$-parameter. Our results indicate that for any $\alpha$-value greater than $10^{-3}$, the *f*-profile does not decrease sufficiently to simultaneously intersect all three measurements. Specifically, with $\alpha = 2 \times 10^{-3}$ and the same initial accretion rate as in our nominal model, the *f*-profile only reaches the upper limit of Makemake's D/H measurement. This is because the *f*-value depends on turbulent viscosity (see Equation (6)), and higher values of $\alpha$ lead to a higher enrichment of deuterium for the same set of other parameters.

Our model allows for some tuning through the viscosity and accretion-rate parameters, enabling a range of D/H ratios between the protosolar value and that measured in 67P. However, this flexibility is intrinsic to disk evolution models, where temporal evolution naturally governs isotopic fractionation. The model remains constrained by physical processes such as vapor diffusion, isotopic exchange, and condensation, which shape the range of possible outcomes.

A key test of our approach would be independent measurements of D/H ratios in methane across different cometary bodies. If comets were found to exhibit D/H ratios in methane comparable to those measured on Makemake and Eris, or intermediate between the protosolar value and that observed in 67P, this would provide a strong consistency check for the D/H pattern predicted by our model and a foundation for further observational validation. Additionally, a precise determination of the D/H ratio in interstellar methane would be helpful for refining the initial conditions of our simulations.

Our model also predicts the formation of CH₄-rich ices after approximately 60–150 kyr of PSN evolution, with deuterium enrichments matching the values observed in Makemake and Eris. As the CH₄ ice line moves inward, CH₄-rich microscopic grains form, agglomerate with other grains, and grow into pebbles. Once the abundance of icy pebbles exceeds a critical threshold near the ice line, planetesimal formation is triggered via streaming instability, with timescales reaching several megayears for planetesimals of the order of 100 km in size (A. Johansen et al. 2015; J. Drążkowska & Y. Alibert 2017; C. T. Lenz et al. 2020).

The high D/H ratio measured in CH₄ in 67P/C-G aligns with the previously reported D/H value in H₂O ((5.3 ± 0.7) × 10⁻⁴) by K. Altwegg et al. (2015), suggesting that the comet was assembled from nearly pristine materials that experienced minimal reprocessing in the PSN. However, a recent analysis by K. E. Mandt et al. (2024) revises the D/H ratio in H₂O downward to ~ (2.2–2.9) × 10⁻⁴—namely, 1.2–1.6 times the terrestrial value—after accounting for dust contamination that may have artificially elevated the original measurement from the ROSINA instrument. If confirmed, this revision implies that the comet's water underwent significant isotopic exchange with hydrogen before crystallizing in the PSN, whereas CH₄ remained linked to a more primitive reservoir. This coexistence could result from large-scale mixing of solid materials across the PSN, a process already invoked to explain isotopic heterogeneity in meteorites (E. Deloule 1997; E. Deloule et al. 1998).

An important marker of primordial ice is CO, which is not detected on Eris and Makemake, according to the W. M. Grundy et al. (2024) analysis of the JWST data. To reconcile a primordial origin of CH₄ with an apparent absence of CO, several chemical and physical processes must be considered. One possibility, as suggested by C. R. Glein & J. H. Waite (2018), is that CO was destroyed by reactions in a subsurface liquid-water ocean. Another possibility is that CH₄ was preferentially enriched at its ice line in the PSN relative to CO and other volatiles, contributing to the observed depletion (O. Mousis et al. 2021). This hypothesis supports the near absence of detected volatiles other than CH₄ and N₂ on Makemake and Eris (W. M. Grundy et al. 2024), with the observed nitrogen possibly coming from organic matter (K. E. Miller et al. 2019).

Our model suggests that the methane observed on Eris and Makemake could originate from the PSN, while not excluding the possibility of endogenous production within these bodies (C. R. Glein et al. 2024). Although this does not confirm endogenous methane production, the differentiation of Eris (F. Nimmo & M. E. Brown 2023) and the potential cryovolcanic activity on Makemake (C. Kiss et al. 2024) are in agreement with this hypothesis. However, our study indicates the potential coexistence of two methane reservoirs: one primordial and another produced in situ. Furthermore, our findings support the idea that these KBOs accreted from solids formed in an environment with a high diversity of ices.

Finally, our model focuses only on gas diffusion and does not consider the transport of grains and pebbles with different D/H ratios across the disk. This contrasts with the detection of crystalline silicates in comets, suggesting their formation near the protosun followed by efficient radial transport (D. Bockelée-Morvan et al. 2002; O. Mousis et al. 2007; D. Brownlee 2014). Addressing this limitation is an important direction for future work, as the resulting deuterium distribution in solids formed from dust and pebbles may differ from the idealized patterns presented here. A similar consideration applies to C. R. Glein et al. (2024), who proposed that the D/H ratio in H₂O and CH₄ remains constant in bodies that accreted these volatiles from the PSN. This assumption seems questionable, as it assumes identical chemical pathways for both molecules from the ISM to the PSN and neglects the potential effects of radial mixing on their pure condensates.

## Acknowledgments

The project leading to this publication has received funding from the Excellence Initiative of Aix-Marseille Université—
A\*Midex, a French “Investissements d'Avenir program,” AMX-21-IET-018. This research is part of the project FACOM (ANR-22-CE49-0005-01_ACT) and has benefited from funding provided by l'Agence Nationale de la Recherche (ANR) under the Generic Call for Proposals 2022. O.M. acknowledges funding from CNES.

## ORCID iDs

Olivier Mousis https://orcid.org/0000-0001-5323-6453

Aaron Werlen https://orcid.org/0009-0005-1133-7586

Tom Benest Couzinou https://orcid.org/0000-0002-8719-7867

Antoine Schneeberger https://orcid.org/0000-0002-3289-2432

## References

Aguichine, A., Mousis, O., Devouard, B., & Ronnet, T. 2020, ApJ, 901, 97

Altwegg, K., Balsiger, H., Bar-Nun, A., et al. 2015, Sci, 347, 1261952

Atreya, S. K., Adams, E. Y., Niemann, H. B., et al. 2006, P&SS, 54, 1177

Bell, K. R., & Lin, D. N. C. 1994, ApJ, 427, 987

Bockelée-Morvan, D., Gautier, D., Hersant, F., Huré, J. M., & Robert, F. 2002, A&A, 384, 1107

Bonev, B. P., Mumma, M. J., Gibb, E. L., et al. 2009, ApJ, 699, 1563

Brownlee, D. 2014, AREPS, 42, 179

Cassen, P., & Chick, K. M. 1997, in AIP Conf. Ser. 402, Astrophysical Implications of the Laboratory Study of Presolar Materials, ed. T. J. Bernatowicz & E. Zinner (Melville, NY: AIP), 697

Ceccarelli, C., Caselli, P., Bockelée-Morvan, D., et al. 2014, in Protostars and Planets VI, ed. H. Beuther et al. (Tucson, AZ: Univ. Arizona Press), 859

Clark, R. N., Brown, R. H., Cruikshank, D. P., & Swayze, G. A. 2019, Icar, 321, 791

Cleeves, L. I., Bergin, E. A., O'D. Alexander, C. M., et al. 2016, ApJ, 819, 13

Deloule, E. 1997, M&PSA, 32, A34

Deloule, E., Robert, F., & Doukhan, J. C. 1998, GeCoA, 62, 3367

Desch, S. J., Estrada, P. R., Kalyaan, A., & Cuzzi, J. N. 2017, ApJ, 840, 86

Drażkowska, J., & Alibert, Y. 2017, A&A, 608, A92

Drouart, A., Dubrulle, B., Gautier, D., & Robert, F. 1999, Icar, 140, 129

Dubrulle, B., & Frisch, U. 1991, PhRvA, 43, 5355

Gibb, E. L., Bonev, B. P., Villanueva, G., et al. 2012, ApJ, 750, 102

Glein, C. R., Grundy, W. M., Lunine, J. I., et al. 2024, Icar, 412, 115999

Glein, C. R., & Waite, J. H. 2018, Icar, 313, 79

Grundy, W. M., Wong, I., Glein, C. R., et al. 2024, Icar, 411, 115923

Gullbring, E., Hartmann, L., Briceno, C., & Calvet, N. 1998, ApJ, 492, 323

Hallis, L. J. 2017, RSPTA, 375, 20150390

Hartmann, L., Calvet, N., Gullbring, E., & D'Alessio, P. 1998, ApJ, 495, 385

Hueso, R., & Guillot, T. 2005, A&A, 442, 703

Johansen, A., Mac Low, M.-M., Lacerda, P., & Bizzarro, M. 2015, SciA, 1, 1500109

Kawakita, H., Watanabe, J.-i., Furusho, R., Fuse, T., & Boice, D. C. 2005, ApJL, 623, L49

Kislitsyn, P. A., Balashev, S. A., Murphy, M. T., et al. 2024, MNRAS, 528, 4068

Kiss, C., Müller, T. G., Farkas-Takács, A., et al. 2024, ApJL, 976, L9

Lécluse, C., & Robert, F. 1994, GeCoA, 58, 2927

Lécluse, C., Robert, F., Gautier, D., & Guiraud, M. 1996, P&SS, 44, 1579

Lenz, C. T., Klahr, H., Birnstiel, T., Kretke, K., & Stammler, S. 2020, A&A, 640, A61

Lodders, K. 2019, arXiv:1912.00844

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

Mandt, K. E., Lustig-Yaeger, J., Luspay-Kuti, A., et al. 2024, SciA, 10, eadp2191

Miller, K. E., Glein, C. R., & Waite, J. H., Jr. 2019, ApJ, 871, 59

Mousis, O., Aguichine, A., Bouquet, A., et al. 2021, PSJ, 2, 72

Mousis, O., Gautier, D., & Bockelée-Morvan, D. 2002a, Icar, 156, 162

Mousis, O., Gautier, D., Bockelée-Morvan, D., et al. 2000, Icar, 148, 513

Mousis, O., Gautier, D., & Coustenis, A. 2002b, Icar, 159, 156

Mousis, O., Lunine, J. I., Pasek, M., et al. 2009a, Icar, 204, 749

Mousis, O., Lunine, J. I., Thomas, C., et al. 2009c, ApJ, 691, 1780

Mousis, O., Lunine, J. I., Waite, J. H., et al. 2009b, ApJL, 701, L39

Mousis, O., Petit, J. M., Wurm, G., et al. 2007, A&A, 466, L9

Mousis, O., Schneeberger, A., Cavalié, T., et al. 2024, PSJ, 5, 173

Müller, D. R., Altwegg, K., Berthelier, J. J., et al. 2022, A&A, 662, A69

Nesvorný, D. 2018, ARA&A, 56, 137

Nimmo, F., & Brown, M. E. 2023, SciA, 9, eadi9201

Nixon, C. A., Temelso, B., Vinatier, S., et al. 2012, ApJ, 749, 159

Ohno, K., & Ueda, T. 2021, A&A, 651, L2

Piani, L., Marrocchi, Y., Vacher, L. G., Yurimoto, H., & Bizzarro, M. 2021, E&PSL, 567, 117008

Richet, P., Bottinga, Y., & Janoy, M. 1977, AREPS, 5, 65

Roueff, E., Lis, D. C., van der Tak, F. F. S., Gerin, M., & Goldsmith, P. F. 2005, A&A, 438, 585

Sakai, N., Shirley, Y. L., Sakai, T., et al. 2012, ApJL, 758, L4

Schneeberger, A., Mousis, O., Aguichine, A., & Lunine, J. I. 2023, A&A, 670, A28

Shakura, N. I., & Sunyaev, R. A. 1973, A&A, 24, 337

Tsiganis, K., Gomes, R., Morbidelli, A., & Levison, H. F. 2005, Natur, 435, 459

Waite, J. H., Jr., Lewis, W. S., Magee, B. A., et al. 2009, Natur, 460, 487