ORIGIN OF MOLECULAR OXYGEN IN COMET 67P/CHURYUMOV–GERASIMENKO

The Rosetta Orbiter Spectrometer for Ion and Neutral Analysis (ROSINA) Double Focusing Mass Spectrometer (DFMS) on board the Rosetta spacecraft (Balsiger et al. 2007) enabled the detection of O₂ in the coma of comet 67P/Churyumov–Gerasimenko (67P/C–G) with local abundances in the 1%–10% range and a mean value of 3.80 ± 0.85% (Bieler et al. 2015). A subsequent reinvestigation of the 1P/Halley data from the Giotto Neutral Mass Spectrometer also indicates that the coma of comet 1P/Halley should contain O₂ with an abundance of 3.7 ± 1.7% with respect to water, suggesting that this molecule may be a rather common parent species in comets (Rubin et al. 2015b).

To investigate the origin of O₂ in 67P/C–G, Bieler et al. (2015) considered the possibility of O₂ production via the radiolysis of water ice incorporated within the nucleus. Based on 67P/C–G's known orbital history, they estimated that any O₂ produced during the residence time of 67P/C–G in the Kuiper Belt was quickly lost during the first pass or two around the Sun. The authors further found that radiolysis on closer orbit to the Sun would most likely only affect the top few micrometers of the nucleus' active surface. In this case, the O₂/H₂O ratio produced in these conditions would decrease with depth. Because they did not observe any variation of the O₂/H₂O ratio during the sampling period, Bieler et al. (2015) ruled out the hypothesis of O₂ production via the radiolysis and determined that O₂ must have been incorporated into 67P/C–G at the time of its formation in the protosolar nebula (PSN).

In order to explain how O₂ could have been incorporated into the ices of 67P/C–G, we investigate here the radiolysis hypothesis at epochs prior to the formation of comets, when icy grains were the dominant solid phase of the outer PSN. Furthermore, we examine the different trapping scenarios of O₂ that could explain its presence. Because some recent works suggest that this comet may have been accreted from a mixture of clathrates and pure crystalline ices formed in the PSN (Luspay-Kuti et al. 2016; Mousis et al. 2016), we study the propensity for O₂ trapping in clathrates, and also evaluate if its condensation as pure crystalline ice is consistent with the comet's inferred composition. Among all these investigated mechanisms, we find that clathration of O₂ is efficient in the PSN and that radiolysis can explain the formation of O₂ and its stabilization in icy grains. However, to produce enough O₂ molecules, the radiolysis of icy grains must have happened in a low-density environment such as the presolar cloud.

We first investigate the possibility of radiolytic production of O₂ in icy grains present in the outer PSN prior to their agglomeration by 67P/C–G. The energy available for radiolysis is provided by the galactic cosmic ray flux (CRF) impacting icy grains. In the following, since galactic CRF can penetrate into water ice down to depths of a few tens of meters (Cooper et al. 1998), we only consider icy grains with sizes below this limit, implying that no H₂O ice can be out of reach of radiolysis. In our calculations, we use the energy range and CRF distribution from Yeghikyan (2011) and Cooper et al. (2003), respectively. The CRF energy dose absorbed by icy grains located at 30 au from the Sun is within the ∼(5–60) × 10¹⁶ eV kg⁻¹ yr⁻¹ range, depending on the disk's surface density (between 10 and 10³ g cm⁻²; see Hersant et al. 2001).

O₂ is produced by radiolysis of water ice through the chemical reaction 2 H₂O $\longrightarrow$ 2 H₂ + O₂, with an amount of energy needed to alter one H₂O molecule W_r = 235 eV (Johnson 1991). H₂ is then rapidly lost from the water ice due to its fast diffusion. Further reactions with O₂ are precluded because the diffusion of these molecules is slowed down by the disk's low temperatures (Johnson 1990). We have thus assumed that all the energy absorbed by water ice is used to form O₂. To reach the molecular ratio O₂/H₂O (1%–10%) measured by Bieler et al. (2015) in 67P/C–G, cosmic rays must alter twice as many H₂O molecules in icy grains. The time τ needed to reach this ratio is then given by

$$\begin{eqnarray}&&\tau =\displaystyle \frac{{W}_{r}\cdot {N}_{A}}{{E}_{\mathrm{CR}}\cdot {M}_{{{\rm{H}}}_{2}{\rm{O}}}}\times {f}_{{{\rm{H}}}_{2}{\rm{O}}}\end{eqnarray} \tag{ 1 }$$

where N_A (mol⁻¹) is the Avogadro constant, ${M}_{{{\rm{H}}}_{2}{\rm{O}}}$ (kg mol⁻¹) is the molar mass of water, E_CR (eV kg⁻¹ yr⁻¹) is the CRF energy dose received by water ice and ${f}_{{{\rm{H}}}_{2}{\rm{O}}}$ is the fraction of altered H₂O molecules, which corresponds to two times the fraction of O₂ produced.

Figure 1 shows the results of our calculations. An O₂ fraction in the 1%–10% range is reached in ∼0.25–2.5 Gyr at the aforementioned nominal CRF value (Case 1). These extremely long time periods are incompatible with the lifetime of icy grains in the PSN (a few 10⁴ year; Weidenschilling & Cuzzi 1993). If icy grains have grown to sizes larger than tens of meters in the PSN, then the deepest layers should remain unaltered. In this case, even longer timescales would be needed for O₂ formation. However, the CRF may have undergone significant enhancements throughout the history of the solar system, by a factor of ∼3 during its passages through the Milky Way's spiral arms (a few tens of Myr every 400–500 Myr; Effenberger et al. 2012; Werner et al. 2015; Alexeev 2016), or even by a factor of ∼100 during a few kyr because of a close supernova explosion (<30 pc; Fields & Ellis 1999). Such enhancements can decrease the time needed to form O₂ by up to a factor of 100, which is still too long for our consideration.

Zoom In Zoom Out Reset image size

We also consider the possibility of an icy grain receiving the maximum CRF energy dose estimated by Yeghikyan (2011), namely ∼1.20 × 10²⁰ eV kg⁻¹ yr⁻¹. This value leads to a time τ in the ∼1–10 Myr range (see case 2 of Figure 1), or ∼10–100 kyr with a CRF enhanced by a factor of 100. However, such a high value of E_CR corresponds to a surface density of 10⁻³ g cm⁻², which can only be reached in molecular clouds. In such environments, the column densities would be low enough to form 1%–10% of O₂ in the icy grains even on very short timescales. Therefore, to incorporate significant amounts of O₂ produced via radiolysis of icy grains, cometary grains must have formed in the presolar cloud prior to disk formation.

An important question is whether O₂ molecules produced via radiolysis of ice grains can remain stabilized within the water icy matrix of 67P/C–G. The stabilization energy is defined as the difference between the energy of the system of O₂ interacting with the ice and the sum of the energies of the pure ice and O₂ at infinite separation. To investigate this problem, a sampling of the representative structures of O₂ in solid water ice has been obtained using a strategy based on first principle periodic density functional theory quantum calculations, that has been proven to be appropriate for modeling bulk and surface ice structures (Lattelais et al. 2011, 2015; Ellinger et al. 2015). Among the different forms, we considered the apolar variety of hexagonal ice Ih because these structures have a balanced distribution of alternating hydrogen and oxygen avoiding computational artifacts for surface optimizations and at the same time reproduce the bulk properties (Casassa et al. 2005). How O₂ behaves as a function of the number of H₂O molecules removed is illustrative of the storage capability of the ice as a function of porosity. The results of our calculations, performed using the Vienna ab initio simulation package (Kresse & Hafner 1993, 1994), are presented below.

• 1.  
  Starting with no H₂O removed, i.e., the pure cristalline ice, we found no stabilization for the inclusion of O₂ in the hexagonal lattice. It is in fact an endothermic process.
• 2.  
  With one H₂O removed, and replaced by one O₂, we have a substitution structure whose stabilization, in the order of 10⁻³ eV, is meaningless.
• 3.  
  With 2, 3, and 4 adjacent H₂O molecules removed from the hexagonal lattice we obtained well defined cavities that, after reconstruction, show different shapes according to the positions of the entities removed. The stabilization energies were found to be on the order of 0.2–0.3 eV, going to 0.4–0.5 eV for an embedded O₂ dimer. A typical structure of embedding is illustrated in Figure 2 where O₂ is stabilized with an energy of ∼0.23 eV. This energy is on the order of that of a water dimer, which means that the presence of O₂ should not perturb the ice structure until it is ejected into the coma via sublimation with the surrounding H₂O molecules.

Zoom In Zoom Out Reset image size

It should be stressed that the formation of one O₂ requires at least the destruction of two H₂O. The present simulation is fully consistent with the aforementioned radiolysis hypothesis, where the irradiation process is at the origin of both the formation of O₂ and the development of the cavity in which it remains sequestered. Similar results are obtained in the case of O₂ stabilization in amorphous ice.

One possible source of O₂ in the nucleus of 67P/C–G is the trapping of O₂ in clathrates that formed in the PSN prior to having been agglomerated by the comet as it formed. This is supported by recent works showing that the Ar/CO and N₂/CO ratios and the time variation of other volatile species measured in 67P/C–G's coma are found to be consistent with the presence of clathrates in its nucleus (Luspay-Kuti et al. 2016; Mousis et al. 2016). To investigate the amount of O₂ that could have been trapped in clathrates and now be present in 67P/C–G, we use the same statistical thermodynamic model as the one described in Mousis et al. (2010, 2016), which is used to estimate the composition of these crystalline structures formed in the PSN. To evaluate the trapping efficiency of O₂, we consider a gas constituted of O₂ and CO. After H₂O, CO is one of the dominant gases found in 67P/C–G (Le Roy et al. 2015) and in most of comets (Bockelée-Morvan et al. 2004, p. 391). The Kihara parameters for the molecule-water interactions employed in our calculations are derived from Mohammadi et al. (2003) for O₂ and from Mohammadi et al. (2005) for CO. These represent the most recent sets of data found in the literature for the two species. We refer the reader to the model description provided in Mousis et al. (2010) for further details.

When clathrates destabilize in the nucleus, the trapped volatiles are released prior to water sublimation, implying that the water vapor measured at the time of the O₂ sampling by ROSINA should be derived from the vaporization of crystalline ice layers located closer to the surface. Hence, the O₂ depletion is better quantified by comparing the O₂/CO ratio in clathrates and the coma value since these two species are expected to be released simultaneously from destabilized clathrates. Figure 3 represents the value of the O₂/CO ratio in structure I clathrates⁷ as a function of the O₂/CO ratio in the coexisting gas phase at a chosen disk's temperature of ∼45 K. This value is within the temperature range needed for clathrates to form in the PSN from a gaseous mixture of protosolar composition that reproduces the Ar/CO and N₂/CO ratios measured in 67P/C–G's coma (Mousis et al. 2016). We find that, whatever the O₂/CO ratio considered for the initial PSN gas phase, it is enriched by a factor of ∼1.4–1.8 in the formed clathrate. Figure 3 also shows that the O₂/CO ratio must be in the 0.026–0.24 range in the PSN gas phase for the clathrate trapping mechanism to agree with the measured range of O₂/H₂O in the coma (∼1%), assuming that all cavities are filled by guest molecules and that the CO/H₂O abundance ratio in the coma corresponds to the sampled value (∼2.7%–21%; Le Roy et al. 2015). This range of O₂/CO ratios is consistent with values obtained at distances beyond ∼5 au in a T Tauri disk (Walsh et al. 2015). Therefore, our calculations show that the clathration of O₂ in the PSN is a realistic mechanism to account for the O₂/H₂O ratio observed by ROSINA in 67P/C–G's coma.

Zoom In Zoom Out Reset image size

An alternative possibility for the observed presence of O₂ in the coma of 67P/C–G is that the O₂ could have been agglomerated as pure crystalline ice by the nucleus forming at cooler PSN temperatures than those required for clathration. To investigate this scenario, we calculated the temperature dependence of the equilibrium curves of O₂, CO, N₂, and Ar pure crystalline ices via the use of the polynomial relations reported by Fray & Schmitt (2009). To derive the partial pressures for each gas, we assumed that O, C, N, and Ar exist in protosolar abundances in the PSN (Lodders et al. 2009), with all C and all N in the forms of CO and N₂, respectively. The partial pressure of O₂ is derived from the O₂/CO gas phase ratio (∼33%) predicted beyond the snowline of a T Tauri disk via an extensive chemical model (Walsh et al. 2015). The equilibrium curves of O₂, CO, N₂, and Ar pure crystalline ices are represented along with the equilibrium curve of the CO–N₂–Ar multiple guest clathrate proposed by Mousis et al. (2016) to explain 67P/C–G's composition, as a function of the total PSN pressure in Figure 4. Because the CO–N₂–Ar multiple guest clathrate is by far dominated by CO (see Figure 1 of Mousis et al. 2016), we assume that its partial pressure is the same as for CO crystalline ice. The equilibrium curve of the clathrate is taken from Lectez et al. (2015).

Zoom In Zoom Out Reset image size

From the examination of the condensation sequence presented in Figure 4, we find that the hypothesis of O₂ agglomeration as pure crystalline ice is inconsistent with 67P/C–G's current composition. The fact that Ar/CO and N₂/CO ratios are found to be significantly depleted by factors of ∼90 and 10 in 67P/C–G's coma, respectively, compared to the protosolar values (Balsiger et al. 2015; Rubin et al. 2015a; Mousis et al. 2016), implies that Ar and N₂ cannot form substantial amounts of pure crystalline ices at the formation location of the comet in the PSN (Mousis et al. 2016). Instead, it has been proposed that these volatiles were mostly trapped in CO-dominated clathrates (Mousis et al. 2016). Under these circumstances, because the equilibrium curve of O₂ ice is in the vicinity of those of Ar and N₂ ices, the incorporation of O₂ in this form would require the trapping of larger amounts of Ar and N₂, incidentally leading to quasi protosolar Ar/CO and N₂/CO ratios. This does not agree with the depleted ratios observed in 67P/C–G.

In this study, we have investigated several scenarios that may explain the presence of molecular oxygen in the nucleus of 67P/C–G. Our results are the following:

• 1.  
  Even with a strong CRF enhancement due to the presence of a nearby supernova, we find that the radiolysis of icy grains is not fast enough in the PSN to create amounts of O₂ comparable with those observed in 67P/C–G. Instead, icy grains must be placed in low-density environments such as molecular clouds to allow radiolysis to work efficiently. The irradiation process also favors the stabilization of O₂ molecules in the icy matrix via the development of cavities and is compatible with both amorphous and crystalline ice structures.
• 2.  
  O₂ can be efficiently trapped in clathrates formed in the PSN. The O₂/CO ratio in the clathrate phase is up to ∼2 times the O₂/CO ratio in the coexisting PSN gas phase.
• 3.  
  The incorporation of O₂ as pure crystalline ice is unlikely in 67P/C–G because the condensation of this species in the PSN would imply much larger abundances of Ar and N₂ than those observed in the coma.

Based on these results, and assuming that radiolysis has been the only mechanism for producing O₂, we find that the formation of 67P/C–G is possible in a dense and early PSN in the framework of two extreme scenarios: (1) agglomeration from pristine amorphous icy grains/particles formed in the ISM and (2) agglomeration from multiple guest clathrates including O₂ that formed during the cooling of the disk subsequent to the vaporization of the amorphous icy grains entering the PSN. However, scenario 1 was found inconsistent with ROSINA pre-perihelion observations of volatile abundances in the coma. In contrast, Mousis et al. (2016) and Luspay-Kuti et al. (2016) have shown that scenario 2 could match these data if 67P/C–G agglomerated from a mixture of clathrates and crystalline ices that condensed in the PSN. Also, scenario 2 is compatible with a possible chemical production of O₂ in the PSN gas phase (Walsh et al. 2015). In this picture, whatever the considered source, i.e., radiolysis of ISM grains or/and PSN gas phase chemistry, O₂ is efficiently entrapped in clathrates prior to their agglomeration by 67P/C–G.

On the other hand, with the incorporation of O₂ in the cavities created by CRF in the icy matrix, scenario 1 naturally provides an explanation for the strong correlation found between the O₂ and H₂O production rates observed in 67P/C–G's coma (Bieler et al. 2015). If this scenario is correct, this would make implausible the accretion of 67P/C–G from clathrates and crystalline ices originating from the PSN. Meanwhile, a way to reconcile scenario 2 with the strong O₂–H₂O correlation would be to assume that the icy grains initially formed as in scenario 1. These icy grains/particles would have then subsequently experienced an amorphous-to-crystalline phase transition in the 130–150 K temperature range when entering the disk (Kouchi et al. 1994; Maldoni et al. 2003; Ciesla 2014). In this alternative scenario, all volatiles initially adsorbed by ISM amorphous ice would be released in the PSN gas phase during phase transition. With the cooling of the disk, these volatiles would have been later trapped in the clathrates formed with the crystallized icy grains. The case of O₂ is unique because, due to its formation process, this molecule is inserted into the icy matrix. In spite of the phase transition, O₂ would remain stable within the icy matrix because the strength of the interaction between O₂ and the surrounding H₂O molecules is expected not to decrease (eventually increase) upon crystallization. In this scenario, CO, Ar, and N₂ would be trapped in clathrates with O₂ remaining embedded in water, in a way consistent with the observed correlation.

To conclude, further post-perihelion ROSINA data, in particular the precise measurements of the relative abundances of the different volatiles as a function of geography and time, are needed to disentangle between the existing formation scenarios. It is also possible that only the in situ sampling of a nucleus by a future lander will provide a definitive answer to the question of the formation conditions of 67P/C–G and other Jupiter Family Comets in the PSN.

O.M. acknowledges support from CNES. This work has been partly carried out thanks to the support of the A*MIDEX project (n^o ANR-11-IDEX-0001-02) funded by the "Investissements d'Avenir" French Government program, managed by the French National Research Agency (ANR). This work also benefited from the support of CNRS-INSU national program for planetology (PNP). R.M. was supported by the Belgian Science Policy Office through the Solar-Terrestrial Centre of Excellence and by PRODEX/ROSETTA/ROSINA PEA 4000107705. J.I.L. acknowledges support from JWST. K.E.M. acknowledges support from JPL Subcontract 1345493.