Michael Combi

Yinsi Shou

Nicolas Fougere

Valeriy Tenishev

Kathrin Altwegg

Martin Rubin

Dominique Bockelée-Morvan

Fabrizio Capaccioni

Yu-Chi Cheng

Uwe Fink

Tamas Gombosi

Kenneth C Hansen

Zhenguang Huang

David Marshall

Gabor Toth

The surface distributions of the

Keywords: Comets, Comets, nucleus, Comets, coma, Atmospheres, structure

come   L'archive ouverte pluridisciplinaire

Introduction

The Rosetta spacecraft orbited and maneuvered in the vicinity of the nucleus of comet 67P/Churyumov-Gerasimenko (67P hereafter) for 26 months from August 2014 through September 2016 recording various measurements from 11 instruments on the orbiter. There were 10 more instruments on the Philae Lander. The mission was able to monitor the comet carefully at heliocentric distances from 3.5 AU before perihelion, through perihelion at 1.24 AU and out to 3.5 AU after perihelion. As such the combination of temporal coverage, sensitivity and resolution of all types of measurements was unprecedented.

In previous work we (Fougere et al. 2016a(Fougere et al. , 2016b) ) used measurements by the Double Focused Mass Spectrometer (DFMS) and the Comet Pressure Sensor (COPS) subsystems of the Rosetta Ion and Neutral Analyzer (ROSINA) package to determine the potential surface distribution of the four major ions, H 2 O, CO 2 , CO and O 2 . Those results covered the period from first close encounter in early August 2014 through the end of February 2016, two months before the comet's second equinox. We had used the determined potential surface distribution of H 2 O and CO 2 to compare with images of these two species taken with the VIRTIS-M (Migliorini et al. 2016) moderate spectral resolution imaging spectrometer and with various point spectra taken with the VIRTIS-H high-resolution spectrometer [START_REF] Bockelée-Morvan | First observations of H2O and CO2 vapor in comet 67P/Churyumov-Gerasimenko made by VIRTIS onboard Rosetta[END_REF]2016). The spatial distributions, abundances and derived productions of H 2 O and CO 2 from the DFMS/COPS-determined model were in quite reasonable agreement with both VIRTIS-M images of H 2 O and CO 2 taken in March 2015 and taken with VIRTIS-H spectra taken from November 2014 through January 2015. However, the DFMS/COPS-determined model generally predicted higher abundances by up to a factor of four in the weeks around perihelion of H 2 O than determined by VIRTIS-H and somewhat higher abundances of CO 2 . In this paper we have improved our earlier approach and extended a similar analysis through the end of DFMS data, a month before the end of mission in September 2016.

Reduction and Calibration of the ROSINA DFMS Data

The results presented in this paper all stem from measurements by the DFMS and COPS subsystems of the ROSINA package. The ROSINA package of instruments was described by [START_REF] Balsiger | Rosina Rosetta Orbiter Spectrometer for Ion and Neutral Analysis[END_REF]. The DFMS instrument measures high resolution mass spectra with a resolution of up to m/Δm = 3000 and as such can determine molecular, elemental and isotopic abundances even for otherwise normally overlapping masses such as CO and N 2 [START_REF] Rubin | Molecular nitrogen in comet 67P/Churyumov-Gerasimenko indicates a low formation temperature[END_REF]. The reduction procedure for the data on H 2 O, CO 2 , CO and O 2 presented here has been described in great detail in a series of previous papers (Bieler et al. 2015a;Bieler et al. 2015b;Fougere et al. 2016a;Fougere et al. 2016b). Briefly, DFMS mass spectra were obtained in several individual mass channels for the four major species and were given an absolute calibration using the total gas densities derived from the nude gauge measurements by COPS.

The full set of DFMS data analyzed here is shown in Figure 1, which shows the measured densities times the square of the comet-spacecraft distance. It is roughly broken up into 3 time periods. In the first the north pole is nearly continually exposed to the Sun and the two dominant species, H 2 O and CO 2 , are generally anti-correlated with one another. The second period runs from the first equinox until the second equinox when enough or all of the south polar region is illuminated by the Sun and the two dominant species, H 2 O and CO 2 , are generally correlated with one another. The third period runs to the end of mission and returns to the geometry when and CO 2 , are generally anti-correlated again. We will discuss the relative behavior of all four species later in this paper. Uncertainties in the COPS measurements of the total gas density are ~10%. The DFMS densities of each major species are constructed from individual DFMS spectra at each mass channel and then scaled to the total measured by COPS. The final uncertainties in the individual species densities are ~ 20%.

Potential Surface Activity Distributions: Methods

In an analysis of ROSINA COPS data from the first 5 months of the Rosetta mission, August through December 2014, it was shown using a simple analytic model based on assuming that a uniform potential activity from the surface when illuminated could explain much of the density variation of the gas as the comet rotated and the spacecraft position changed (Bieler et al. 2015a). One exception was that it appeared that some additional activity was required from the so-called neck region of the "duck" shape of the nucleus that is not far from the north pole. This was tested and verified in the same paper using both a hydrodynamic model for the gas flow as well as a fully kinetic Direct Simulation Monte Carlo model. Figure 1. Heliocentric distance dependence of DFMS densities normalized by the inverse square of the spacecraft-comet distance. From the top down the plots are for the H 2 O (blue), CO 2 (orange), CO (green) and O 2 (purple). The lines in each give the power-law fits for the three separate time periods as given in Table 1, dotted for phase 1, solid for phase 2 and dashed for phase 3.

In order to examine this idea more fully and instead of assuming the potential surface activity was driven by illumination uniformly from the surface, Fougere et al. (2016a) analyzed the H 2 O and CO 2 gas densities measured by ROSINA DFMS for the period from August 2014 through May 2015 with a least squares inversion scheme based on expressing the potential surface activity as a spherical harmonics expansion with undetermined coefficients. This covered the time period from the near nucleus measurements soon after rendezvous until the comet reached its first equinox. The technique assumed the same analytic model connecting illumination-driven surface activity to densities at the spacecraft as introduced by Bieler et al. (2015a) and described and used by Fougere et al. (2016aFougere et al. ( , 2016b)). The final potential surface activity distributions for H 2 O and CO 2 were then used as the surface boundary conditions and a physically-based coma was generated using the fully-kinetic Direct Simulation Monte Carlo (DSMC) model.

The DSMC models based on the spherical harmonics inversions were also used by Fougere et al. (2016a) to synthesize a series of H 2 O and CO 2 images obtained in April 2015 by [START_REF] Migliorini | Water and carbon dioxide distribution in the 67P/Churyumov-Gerasimenko coma from VIRTIS-M infrared observations[END_REF] using the Rosetta Visible and Infrared Thermal Imaging Spectrometer (VIRTIS) medium resolution mapping channel and a set of H 2 O and CO 2 column abundances vector of length M elements containing the DFMS data and x is a vector of length N containing the undetermined coefficients of the spherical harmonics for which we are solving. There are 15000 to 25000 data points total for each species. The potential activity

𝑓 ! = 𝑥 ! 𝑌 ! (𝜗 ! , 𝜑 ! ) ! !!!
of any surface triangular facet is constrained to be positive definite. The gas flux of any triangle is given by 𝐹

! = 𝜆 !(! !"# )! ! ! !" !
where the function G give the variation of the flux with solar zenith angle defined by Bieler et al. (2015) as

𝐺(Θ !"# ) = 𝑎 !"!!! + (1 -𝑎 !"#!! )𝑐𝑜𝑠(Θ !"# ) which gives 𝐺(Θ !"# ) = 𝑎 !"#!! if
a surface triangle is not illuminated by the Sun. The factor λ is included here to enable correction of the approximate nature of the initial power-law assumption of the production rate as a function of heliocentric distance. Most of the values of λ range from 0.5 to 4 and serve mostly to correct densities for the initial assumed power laws in heliocentric distance.

Model Analysis

As in the papers by Fougere et al. (2016a) and (2016b) we have employed a 5th order spherical harmonics expansion with 25 undetermined coefficients to the DFMS measurements.

One change for this work is that having observed a noticeable difference in the activity between the pre-perihelion and the now extended post-perihelion variations of the four major volatiles, we have separated the data into three time groups: (1) from the beginning of the mission to the first equinox, (2) from the first to the second equinox, and (3) from the second equinox until the last data before the end of mission.

Figure 1 shows the evolution of the measured densities of the 4 species normalized by the square of the spacecraft-comet distance and fitted with power-laws for the three time periods shown as the overlain straight lines for each of the four species. Table 1 lists the values of the power-law exponents for the three time periods for the each of the four major volatile species.

Figure 2 shows the comparison of the measured DFMS number densities of H 2 O, CO 2 , CO and O 2 with those obtained by applying the potential surface activity distribution using the simple analytical model (Bieler et al. 2015a) that was used in the least squares fitting procedure.

The correlation coefficients for the model fits for the three phases (p1, p2, p3) are (0.971, 0.910, 0.895) for H 2 O, (0.924, 0.940, 0.952) for CO 2 , (0.832, 0.942, 0.947) for CO and (0.869, 0.559, 0.923) for O 2 , respectively. Tables 234give the spherical harmonics coefficients for the three phases, respectively.

.071106e-02 Y -2 sin(2ϕ) sin 2 (θ) cos(θ) 5.920391e-02 -5.579715e-01 -4.069729e-01 2.252411e-02 Y -1 sin(ϕ) sin(θ)5(cos2(θ) -1) -7.746035e-02 8.204630e-03 6.060292e-02 4.076148e-03 Y 0 5 cos 3 (θ) -3 cos(θ) -6.807962e-02 -4.696779e-01 -4.771903e-02 -1.344509e-01 Y 1 cos(ϕ) sin(θ)5(cos 2 (θ) -1) 1.293904e-01 -1.282020e-01 -2.840572e-02 1.707097e-02 Y 2 cos(2ϕ) sin 2 (θ) cos(θ) -4.824502e-02 1.196393e+00 4.563505e-01 -4.698604e-02 Y 3 cos(3ϕ) sin 3 (θ) -1.677468e-02 -1.740398e-01 -7.938387e-02 -9.831394e-03 Y -4 sin(4ϕ) sin 4 (θ) 4.762068e-02 -1.434847e-01 -2.638165e-02 6.271004e-03 Y -3 sin(3ϕ) sin 3 (θ) cos(θ) 7.163975e-02 4.958958e-01 2.879720e-01 -1.438314e-02 Y -2 sin(2ϕ) sin 2 (θ)(7 cos 2 (θ) -1) -1.857199e-02 2.723782e-01 1.180346e-01 -2.667406e-02 Y -1 sin(ϕ) sin(θ)(7 cos 3 (θ) -3 cos(θ)) -1.226328e-01 -2.080173e-01 -1.285519e-01 -5.684545e-02 Y 0 35 cos 4 (θ) -30 cos 2 (θ) + 3 9.768123e-02 8.176546e-02 4.075415e-02 6.118638e-02 Y 1 cos(ϕ) sin(θ)(7 cos 3 (θ) -3 cos(θ)) 7.062120e-02 3.730570e-02 1.036030e-02 1.607577e-02 Y 2 cos(2ϕ) sin 2 (θ)(7 cos 2 (θ) -1) 6.002072e-02 -3.559440e-01 -1.368811e-01 4.581718e-02 Y 3 cos(3ϕ) sin 3 (θ) cos(θ) 5.921993e-04 6.116564e-01 1.687368e-01 -3.823632e-03 Y 4 cos(4ϕ) sin 4 (θ) 1.702889e-02 -1.207067e-02 1.452672e-02 2.428081e-02
Table 4. Real spherical harmonics and the normalized coefficients for time period Phase 3.

Term Expression H 2 O coefficient CO 2 coefficient CO coefficient O 2 coefficient Y 0 0 1 6.457323e-01 3.002227e+00 8.823522e-01 3.938829e-01 Y -1 1 sin(ϕ) sin(θ) -7.648887e-01 2.001145e-01 2.379021e-01 -1.097899e+00 Y 0 1 cos(θ) -2.570970e-02 -6.338377e+00 -1.403207e+00 6.715108e-03 Y 1 1 cos(ϕ) sin(θ) 9.965626e-01 2.506014e+00 9.738838e-01 -1.383084e-02 Y -2 2 sin(2ϕ) sin2(θ) -1.193436e-01 -5.466062e-02 -5.966775e-02 -1.717304e-02 Y -1 2 sin(ϕ) sin(θ) cos(θ) 1.233427e-01 -4.396767e-02 -1.551214e-01 2.092213e-01 Y 0 2 3 cos 2 (θ) -1 -1.605573e-02 3.426135e+00 6.769170e-01 1.822769e-01 Y 1 2 cos(ϕ) sin(θ) cos(θ) -2.665689e-01 -1.688630e+00 -4.867386e-01 4.063866e-02 Y 2 2 cos(2ϕ) sin 2 (θ) -3.243013e-01 -1.322311e-01 4.061452e-02 -1.496820e-01 Y -3 3 sin(3ϕ) sin 3 (θ) 1.734523e-01 -2.713263e-01 -1.446861e-01 1.078008e-04 Y -2 3 sin(2ϕ) sin 2 (θ) cos(θ) 3.966970e-01 2.877859e-01 -1.811304e-01 9.942423e-02 Y -1 3 sin(ϕ) sin(θ)5(cos2(θ) -1) -1.705277e-01 4.682026e-02 5.778836e-02 -2.493769e-01 Y 0 3 5 cos 3 (θ) -3 cos(θ) 1.261865e-01 -2.076162e+00 -4.891634e-01 1.999027e-01 Y 1 3 cos(ϕ) sin(θ)5(cos 2 (θ) -1) 2.330877e-01 4.908625e-01 2.267818e-01 -1.099441e-02 Y 2 3 cos(2ϕ) sin 2 (θ) cos(θ) -5.384119e-02 6.938674e-01 2.201555e-01 1.672629e-01 Y 3 3 cos(3ϕ) sin 3 (θ) 1.190484e-01 -2.192789e-01 -1.541520e-01 -2.698796e-02 Y -4 4 sin(4ϕ) sin 4 (θ) 1.453714e-01 -7.812493e-02 7.363468e-02 2.176997e-02 Y -3 4 sin(3ϕ) sin 3 (θ) cos(θ) 6.825370e-02 5.849854e-01 5.239300e-01 -1.098071e-01 Y -2 4 sin(2ϕ) sin 2 (θ)(7 cos 2 (θ) -1) 3.076581e-03 -4.271857e-02 7.328999e-02 8.599415e-03 Y -1 4 sin(ϕ) sin(θ)(7 cos 3 (θ) -3 cos(θ)) -3.308553e-02 -8.651249e-02 -2.610426e-01 5.112065e-02 Y 0 4 35 cos 4 (θ) -30 cos 2 (θ) + 3 7.553597e-02 2.726439e-01 9.291173e-02 2.212457e-02 Y 1 4 cos(ϕ) sin(θ)(7 cos 3 (θ) -3 cos(θ)) -1.062888e-01 -2.617530e-01 -1.241888e-01 1.643481e-01 Y 2 4 cos(2ϕ) sin 2 (θ)(7 cos 2 (θ) -1) 1.755188e-01 -1.708951e-01 -1.398572e-01 -1.762693e-02 Y 3 4 cos(3ϕ) sin 3 (θ) cos(θ) -2.466417e-02 3.037333e-01 2.470421e-01 2.570353e-01 Y 4 4 cos(4ϕ) sin 4 (θ)
-1.072578e-01 -9.622150e-02 5.915493e-02 -1.031560e-02

There are likely changes in the near-surface composition throughout the mission that are responsible for the variations of the potential surface activity maps during the three mission phases. Such changes were seen in the color of the surface, which became more blue during the middle phase including perihelion [START_REF] Fornasier | Rosetta's comet 67P/Churyumov-Gerasimenko sheds its dusty mantle to reveal its icy nature[END_REF] indicating the presence of more water ice near the surface caused by the removal of surface dust in combination with deposition of water ice through internal recondensation (DeSanctis et al. 2015;Ciarniello et al. 2016).

Given the new surface distributions we have run new sets of DSMC models now spanning the entire Rosetta mission. Table 5 gives the times, heliocentric distances and solar latitudes for the 15 dates for which a set of 12 DSMC models were run, each set spanning one nucleus rotation. Figure 4 shows a comparison of the DFMS densities and those determined from the surface distributions and the DSMC model that follow power-law slopes discussed above.

Gas production rates

Figure 5 shows a plot of the four gas production rates throughout the mission plotted as functions of heliocentric distance and time. For the first two studied periods the variations of productions rates of the four species are similar to the previous results (Fougere et al. 2016b).

The most remarkable behavior of the post-second-equinox period (p3) is the leveling of the CO 2 and CO production rates and the precipitous drop of the H 2 O and O 2 production rates. Surprisingly, although the northern hemisphere of the comet is again primarily illuminated after the second equinox as it was early in the mission before the first equinox, the water production rate drops continuously to well below the CO 2 production rate and nearly equaling that of CO.

Throughout the early period before the first equinox when the northern hemisphere was primarily illuminated, the water production rate was factors of 2 to 3 larger than those of CO 2 and CO at the same large heliocentric distances (~3.5 AU) and with the subsolar latitude still in being the north but being +20 o instead of +45 o .

The result of calculating potential surface activity maps for three time periods compared with the previous results (Fougere et al. 2016b) which used a single set of maps for each species for all times up to approximately the second equinox is overall a lowering of the production rates by about 30% with the peak values shortly after perihelion coming down by over 35%. Fougere et al. (2016b) and [START_REF] Hansen | Evolution of water production of 67P/Churyumov-Gerasimenko: An empirical model and a multi-instrument study[END_REF] had a peak water production rate of 4.0 x 10 28 s -1 whereas the new map has a peak value of 2.8 x 10 28 s -1 and with values within a month or so of perihelion down to about 1 x 10 28 s -1 which is closer to previous determination of maximum water production rates in in 67P based on ground-based and space-based measurements [START_REF] Bertaux | The water production rate of Rosetta target Comet 67P/Churyumov-Gerasimenko near perihelion in 1996, 2002 and 2009 from Lyman α observations with SWAN/SOHO[END_REF][START_REF] Schleicher | Compositional and physical results for Rosetta's new target Comet 67P/Churyumov Gerasimenko from narrowband photometry and imaging[END_REF], Hansen et al. 2016). This narrows the gap in the previously reported discrepancy between the ROSINA DFMS water production rates (Fougere et al. 2016b;[START_REF] Kramer | Seasonal changes of the volatile density in the coma and on the surface of comet 67P/Churyumov-Gerasimenko[END_REF]) and some remote sensing values, but does not completely eliminate it. We return to this comparison in the following section. [START_REF] Kramer | Seasonal changes of the volatile density in the coma and on the surface of comet 67P/Churyumov-Gerasimenko[END_REF] performed a similar surface activity map exercise on ROSINA COPS data. Their peak value was 4.0 x 10 28 s -1 but only when using an initial temperature that resulted in an outflow velocity of 1.5 km s -1 , which is much larger than our maximum value of only 0.8 km s -1 that is more typical of previous measurements of outflow velocities in comets [START_REF] Tseng | Cometary water expansion velocity from OH line shapes[END_REF][START_REF] Shou | A New 3D Multi-fluid Dust Model: A Study of the Effects of Activity and Nucleus Rotation on Dust Grain Behavior at Comet 67P/Churyumov-Gerasimenko[END_REF]) and in reasonable agreement with those determined directly from observations of MIRO [START_REF] Lee | Spatial and diurnal variation of water outgassing on comet 67P/Churyumov-Gerasimenko observed from Rosetta/MIRO in August 2014[END_REF][START_REF] Biver | Long-term monitoring of the outgassing and composition of comet 67P/Churyumov-Gerasimenko with the Rosetta/MIRO instrument[END_REF]. For the [START_REF] Kramer | Seasonal changes of the volatile density in the coma and on the surface of comet 67P/Churyumov-Gerasimenko[END_REF] case of 0.8 km s -1 they find a peak production rate of about 2 x 10 28 s -1 , which is in better agreement with our new value. [START_REF] Läuter | Surface localization of gas sources on comet 67P/Churyumov-Gerasimenko based on DFMS/COPS data[END_REF] also used an inversion model like that of [START_REF] Kramer | Seasonal changes of the volatile density in the coma and on the surface of comet 67P/Churyumov-Gerasimenko[END_REF] in order to calculate the time-averaged surface production rates of the same major volatile species from the DFMS/COPS data set. There are several important differences between our method and theirs. They calculate essentially a rotation-averaged production fluxes from the surface, whereas we separate the effect of solar illumination explicitly and calculate the potential surface production distribution. This allows us to reproduce the coma at any given time for direct comparison with any given remote sensing line of sight and solar illumination direction. This is of critical importance for the kind of multi-instrument comparison we show in the following section of this paper. [START_REF] Läuter | Surface localization of gas sources on comet 67P/Churyumov-Gerasimenko based on DFMS/COPS data[END_REF] calculate the surface distribution using many surface triangles providing a potentially more detailed distribution of sources on the surface whereas, as discussed earlier, we describe the surface distributions with spherical harmonics. Our final results are achieved by using our inverted potential surface distribution as the input conditions for fully-kinetic DSMC models that account explicitly for the solar illumination variation of the rotation-induced activity changes as well as the resulting time-dependent surface temperature distribution given by the thermophysical nucleus model of Davidsson and Gutierrez (2004[START_REF] Davidsson | Nucleus properties of Comet 67P/Churyumov Gerasimenko estimated from non-gravitational force modeling[END_REF][START_REF] Schleicher | Compositional and physical results for Rosetta's new target Comet 67P/Churyumov Gerasimenko from narrowband photometry and imaging[END_REF] which yields outflow speeds that vary with local solar illumination, rotation and heliocentric distance. [START_REF] Läuter | Surface localization of gas sources on comet 67P/Churyumov-Gerasimenko based on DFMS/COPS data[END_REF] calculated production rates over 11 day intervals whereas we can determine production rates from individual measurements following the rotationally dependent solar illumination induced effects on surface heating and then can average the results over any desired time intervals when appropriate. Between the two approached there are potentially small differences in production rates owing to the different ways that outflow speeds are determined, but these differences are rather small. In any case there is reasonable agreement with the overall variation of production rates of all four major species between the two approaches. Synthetic fluorescence spectra were computed using the model of [START_REF] Crovisier | Fluorescence Cascades of Water and Carbon Dioxide in the Emission Spectrum of Comet 9P/Tempel 1. Deep Impact as a World Observatory Event: Synergies in Space, Time, and Wavelength[END_REF], assuming a Boltzmann distribution for the population of the rotational levels of the ground vibrational state, and an ortho-para ratio equal to 3. By spectral fitting, this time we derived the rotational temperature (T rot ), which best fits the relative intensities of the 2.775--2.910 µm hot-band rovibrational lines. The column density was then evaluated using the fluorescence efficiency (gfactor) in this wavelength range for a temperature equal to T rot (Villanueva, personal communication, and Villanueva et al. 2012). Column densities derived from the 2.9-µm hotbands are about 35% higher than those published in Bockelée-Morvan et al. ( 2016). The origin of this difference is still under study. Table A1 gives the results of the column densities from the 2.775--2.910 µm hot-band, those from the previous results from the 4.45--5.0 µm hot bands and the ratios between the two.

Comparison

We have calculated the average ratios of the column densities from the new 2.8 µm hot bands to the column densities derived from our surface inversion model determined from the ROSINA DFMS data. These are summarized in Table 6. As in the previous studies, the abundances determined from the model fit to the ROSINA DFMS data are generally larger than those from VIRTIS-H, but the differences vary between the major species, H 2 O and CO 2 , and with the time during the mission. Although we now have the potential surface inversions performed for three separate time intervals during the mission and new VIRTIS-H column densities, the large difference for H 2 O between the two instruments during the perihelion period has been reduced from a factor of 4 (Fougere et al. 2016b) to a factor of ~3, so the discrepancy is smaller but remains. For H 2 O during the early part of the mission the average ratio is 0.69, but there is a fair amount of scatter in the individual results. For CO 2 in the perihelion period the results are rather consistent with a ratio of 0. 90, but for the early mission the ratio for CO 2 is 0.52. Figure 6 shows the production rates of CO 2 over the whole missing calculated using the potential surface activity results scaled to the DFMS data in 2-day averages (or 4 rotations of the nucleus) compared with production rates derived from the VIRTIS-H column densities of CO 2 using the relative spatial coma distribution taken from our new DFMS determined model integrated along all the VIRTIS-H lines of sight. As expected the differences in production rates are consistent with the ratios of column density given in Table 6.

Our peak value of 2.8 x 10 28 s -1 derived from the ROSINA DFMS data for water production is 50% lower than that found in our previous work (Fougere et al. 2016b) by using a potential surface activity distribution determined from ROSINA DFMS measurements from the perihelion phase of the comet between the two equinoxes rather than one for the whole apparition. The comparison of production rates between DFMS and VIRTIS-H in the preperihelion period around January 2015 shows that while the scatter overlaps, the DFMSdetermined values are on the average larger than the VIRTIS-H values, as before. Similarly Here we used the water column densities obtained from the over 4000 individual observations, concentrated during the period between the two equinoxes to calculate water production rates using the coma distributions obtained from our model that was determined from the spherical harmonics inversion of the DFMS measurements. These data appear in a supplementary file associated with this paper (MIRO_Columns.txt). The MIRO column densities were derived from lookup tables, where measured MIRO spectral line area ratios for H 2 18 O and H 2 16 O were compared to simulated line area ratios. To get the simulated line area ratios, line-of-sight profiles for number density, velocity and kinetic temperature were made, assuming a [START_REF] Haser | Distribution d'intensité dans la tête d'une comète[END_REF] model comet with an expanding coma and decreasing temperature.

The profiles used in [START_REF] Marshall | Spatially resolved evolution of the local H 2 O production rates of coet 67P/Churyumov,Gerasimenko, from the MIRO instrument of Rosetta[END_REF] to describe the density, velocity and temperature are taken from parameterizations made by [START_REF] Lee | Spatial and diurnal variation of water outgassing on comet 67P/Churyumov-Gerasimenko observed from Rosetta/MIRO in August 2014[END_REF], which were inspired by DSMC simulations from [START_REF] Lee | Non-LTE radiative transfer for sub-millimeter water lines in Comet 67P/Churyumov-Gerasimenko[END_REF] and [START_REF] Davidsson | Gas kinetics and dust dynamics in low-density comet comae[END_REF]. By varying the parameters, simulated curves for the change in the line area ratio of H 2 16 O and H 2 18 O with column density were made and thus the measured line area ratios could be used to derive an estimate for the column density. The sensitivity of the MIRO column density extends mostly from 1-50 km from the surface. There is little sensitivity at the very surface and for long columns near perihelion (~400 km nucleus/spacecraft distance) there is not much contribution close to the spacecraft. The method used here for determining production rates from the MIRO nadirpointing column densities is the same as we have done for VIRTIS observations above and in the earlier papers (Fougere et al. 2016a and2016b). Just as shown by the original Haser model determinations by [START_REF] Marshall | Spatially resolved evolution of the local H 2 O production rates of coet 67P/Churyumov,Gerasimenko, from the MIRO instrument of Rosetta[END_REF], there is considerable scatter among the results for all 4000 values. This is similar to what we find if use all 25,000 individual DFMS measurements.

Here we show water production rates on 2-day averages over the whole apparition for both DFMS and MIRO measurements. This corresponds to about 4 rotations of the nucleus and thus yields rotation-averaged values. For the MIRO data we take 2-day averages, which include intervals that are determined by no less than 8 points each.

Figure 7 shows a comparison of all DFMS, MIRO and VIRTIS-H water production rates over the entire mission. The agreement between the DFMS and MIRO measurements is quite good and shows a difference between both of these and VIRTIS-H, especially in the perihelion period. One caveat to the calculation of the water production rates using our coma model and the MIRO data is that the MIRO observations were made with 30-minute integrations and our relative coma model column density was extracted for the mid-point time of the integration. [START_REF] Shinnaka | Imaging Observations of the Hydrogen Coma of Comet 67P/Churyumov-Gerasimenko in 2015 September by the PROCYON/LAICA[END_REF]. They found water production rates of (1.46±0.47) x 10 28 s -1 , (1.24±0.40)

x 10 28 s -1 and (1.30±0.42) x 10 28 s -1 on 7.40, 12.37 and 13.17 September 2015, respectively.

These values are in reasonable agreement with the contemporaneous measurements from ROSINA DFMS and MIRO but larger than those from VIRTIS-H.

Unusual composition variation from 1 August -31 October 2014

In the first 90 days of phase 1 while the water production rate drops somewhat, the CO 2 and CO production rates drop precipitously, and the O 2 production rate increases greatly. We have generated potential surface activity distributions for just the first 90 days to see if there was a change that could, in combination with solar illumination, be responsible for the unusual variations in production rates. Figure 8 shows a comparison of the potential surface activity distributions of all four species in period 1 (before the first equinox) with those for just the first 90 days. First, the water source distribution is not much different in the first 90 days compared to that determined for the whole of phase 1. The O 2 distribution for the first 90 days looks similar to that of water, but during the whole of phase 1 potential activity increases in the south. The distributions for CO 2 and CO do change.

During the first 90 days, CO has much more of its activity in the north than to the south compared with the north/south activity throughout phase 1 where the activity is more favored in the south. So there is more potential activity during the first 90 days when illumination heavily favors the north, but shifts to more to the south throughout the rest of phase 1. So it is likely there is some activation of the south polar region by the increasing solar exposure as the comet approached the first equinox. For CO 2 potential activity remains concentrated to the south, though the longitude distribution is more spread out in the remainder of phase 1 than in the first 90 days.

Mass and surface loss of comet 67P throughout the apparition

We have taken the volatile mass loss rate distribution over the surface of comet 67P using the potential surface activity and estimated the surface erosion loss as a function of location on the surface over the entire two-year mission. The results of the individual and total volatile mass loss rates at the end of the three studied time phases are given in Table 7. Not surprisingly nearly all of the volatile mass loss occurs in period 2 between the first and second equinoxes. We find a water mass loss of (4.9±1.5) x 10 9 kg and a total mass loss for the major volatiles of (6.1±1.9) x 10 9 kg. Consistent with our new lower production rates, this is below the total mass loss calculated by [START_REF] Hansen | Evolution of water production of 67P/Churyumov-Gerasimenko: An empirical model and a multi-instrument study[END_REF] using ROSINA COPS data and the previous coma models of Fougere et al. (2016b). 1.2 x 10 7 1.3 x 10 8 1.3 x 10 8 Total 2.7 x 10 8 6.0 x 10 9 6.1 x 10 9

Previously [START_REF] Hansen | Evolution of water production of 67P/Churyumov-Gerasimenko: An empirical model and a multi-instrument study[END_REF] found a water mass loss of 6.4 x 10 9 kg using ROSINA COPS data and the earlier coma model of Fougere et al. (2016b). Owing to our new modeling of the coma in three separate time intervals we find a total water mass loss of 31% lower than [START_REF] Hansen | Evolution of water production of 67P/Churyumov-Gerasimenko: An empirical model and a multi-instrument study[END_REF]. From their Haser model analysis of the nadir pointing MIRO observations, [START_REF] Marshall | Spatially resolved evolution of the local H 2 O production rates of coet 67P/Churyumov,Gerasimenko, from the MIRO instrument of Rosetta[END_REF] found a total water mass loss between the two equinoxes of (2.4±1.1) x 10 9 kg. However, using a 1D spherical Haser model on many nadir-pointing columns that came from much less active regions of the nucleus would likely underestimate the production rates and the resulting total water mass loss. Using a different MIRO data set that consisted of off-limb observations and a different modeling procedure that accounted for an aspherical coma, [START_REF] Biver | Long-term monitoring of the outgassing and composition of comet 67P/Churyumov-Gerasimenko with the Rosetta/MIRO instrument[END_REF] nonetheless found a similar total water mass loss of (2.4±0.1) x 10 9 kg. Our value of (4.9±1.5) x 10 9 kg is a factor of 2 larger but is consistent with the value of (4.8±1.5) x 10 9 kg found by [START_REF] Läuter | Surface localization of gas sources on comet 67P/Churyumov-Gerasimenko based on DFMS/COPS data[END_REF] from the same ROSINA DFMS data set but using a different model analysis procedure.

The total mass loss of the comet, of course, includes both the loss of the nonvolatile refractory component, or dust, as well as the volatiles. Both the dust production rate from the surface and the total net dust loss from the comet remain somewhat of a controversial topic (Fulle et al. 2016a(Fulle et al. , 2016b[START_REF] Fulle | The refractory-to-ice mass ratio in comets[END_REF][START_REF] Levasseur-Regourd | Cometary Dust[END_REF]. While several Rosetta instruments actually detected various sizes and numbers of dust particles and remote-sensing measurements of dust in coma were made by OSIRIS, VIRTIS and MIRO, there still remain some unknown parameters and incomplete sampling of dust particle sizes and masses. There is also a suggestion regarding gas activity lifting of either icy [START_REF] Keller | Insolation, erosion, and morphology of comet 67P/Churyumov-Gerasimenko[END_REF] or refractory [START_REF] Fulle | The refractory-to-ice mass ratio in comets[END_REF] particles, or both, that then fall back to the surface contributing to the surface mantle that may decrease or limit future activity [START_REF] Keller | Seasonal mass transfer on the nucleus of comet 67P/Chuyumov-Gerasimenko[END_REF], with most of the fall back likely transferring mass from the south to the north and consisting of larger and more granular particles.

For our simple calculation we calculate the total surface excavation depth assuming the total nucleus mass loss of 10.5±3.4 x 10 9 kg and total nucleus density 533 kg m -3 from the Radio Science Investigation (RSI) of [START_REF] Pätzold | A homogeneous nucleus for comet 67P/Churyumov-Gerasimenko from its gravity field[END_REF][START_REF] Pätzold | The Nucleus of comet 67P/Churyumov-Gerasimenko -Part I: The global view -nucleus mass, mass-loss, porosity, and implications[END_REF]. Over the whole 2-year period the volatile mass loss of the major species is (6.1±1.9) x 10 9 kg. Accounting approximately for the mass of all the minor volatile species taken together (LeRoy et al. 2015;[START_REF] Rubin | Elemental and Molecular Abundances in Comet 67P/Churyumov-Gerasimenko[END_REF] would raise the total volatile mass loss to about (6.3±2.0) x 10 9 kg. Figure 9 shows plots of the surface loss at the end of each of our three phases over the two-year mission, at the first equinox, at the second equinox and at the end of mission, respectively. Maximum erosion depths in the south polar region over the whole apparition of ~4 m were found for a refractory/ice loss ratio of 0.67, which is consistent with the average erosion depths in the south that were predicted by the model calculations of [START_REF] Keller | Insolation, erosion, and morphology of comet 67P/Churyumov-Gerasimenko[END_REF]. More than 95% of the mass loss occurs during the time between the first and second equinox.

Using the volatile loss from ROSINA would imply a net refractory mass loss of 4 x 10 9 kg or a net refractory/ice mass ratio of 1:1.5. If we take the full uncertainty ranges of the nucleus mass loss and our combined model plus ROSINA uncertainties we find a rather wide range for the refractory/ice mass ratio of 0.7 (-0.7 +1.0) or basically from 0 to 1.7. Depth results can be scaled simply enough for other dust/gas (refractory/ice) ratios, though clearly this does not account for variable dust/gas mass ratios over the surface. [START_REF] Biver | Long-term monitoring of the outgassing and composition of comet 67P/Churyumov-Gerasimenko with the Rosetta/MIRO instrument[END_REF] has recently made a similar estimate using different assumptions for the total volatile/ice loss from the comet starting with the water ice loss calculated from observations of MIRO over the whole mission. Although their water loss is about half of what we calculate from the ROSINA/DFMS results, their total volatile mass loss rate is about 2/3. Therefore, their refractory/ice mass loss ratio is 1.5:1.

Figure 9. Cumulative Erosion Depths at the end of each of three mission phases: 1st equinox, after 2nd equinox and end of mission. Note that by necessity the color contour look-up table is different for the first phase. A dust/gas (or refractory/ice) loss ratio of 1:1.5 is assumed for a nucleus of density 533 kg m -3 . Results can be scaled for other ratios.

A similar analysis by [START_REF] Biver | Long-term monitoring of the outgassing and composition of comet 67P/Churyumov-Gerasimenko with the Rosetta/MIRO instrument[END_REF] that includes also the nadir-pointing MIRO observations yields total volatile mass loss rate that is about 25% larger. These estimates are in contrast to the estimate by [START_REF] Fulle | The refractory-to-ice mass ratio in comets[END_REF] starting with the refractory mass loss from various dust measurements who find a refractory/ice mass loss ratio lower limit of 3.

Finally, these can also be compared with the dust-to-ice ratio calculated from the real part of the permittivity determined by the analysis of the propagation of the radio signal from the Comet Nucleus Sounding Experiment by Radiowave Transmission (CONSERT) instrument [START_REF] Kofman | Properties of the 67P/Churyumov-Gerasimenko interior revealed by CONSERT radar[END_REF][START_REF] Herique | Cosmochemical implications of CONSERT permittivity characterization of 67P/CG[END_REF]. The CONSERT measurements combined with their conclusion that the porosity of the nucleus is in the range of 75-85%, yields a dust-to-ice ratio in the nucleus in the range of 0.4 to 2.6 [START_REF] Kofman | Properties of the 67P/Churyumov-Gerasimenko interior revealed by CONSERT radar[END_REF]. This range includes the separate estimates from the ROSINA/DFMS and MIRO measurements at the lower end and go nearly to the lowest limit from the dust analysis by [START_REF] Fulle | The refractory-to-ice mass ratio in comets[END_REF]. [START_REF] Herique | Cosmochemical implications of CONSERT permittivity characterization of 67P/CG[END_REF] later assumed a priori a dust-to-ice ratio in the range of 4±2 from [START_REF] Rotundi | Dust measurements in the coma of comet 67P/Churyumov-Gerasimenko inbound to the Sun[END_REF] and from that could conclude that the dust-to-ice ratio is larger than 3. Most of the mass as measured in the coma by both mass spectrometry and remote sensing is lost during the period around southern summer, which includes perihelion. It is worth noting that [START_REF] Keller | Seasonal mass transfer on the nucleus of comet 67P/Chuyumov-Gerasimenko[END_REF] suggest that the composition of the coma over the southern hemisphere during southern summer is then likely to be much more representative of that of the nucleus itself.

Summary

We 1 with only some minor changes in the O 2 potential surface activity noted. This is not the case for CO 2 and CO. The potential surface activity distribution of CO has much more of its activity originating in the north than to the south during the first 90 days whereas throughout phase 1 the activity shifts to mostly from the south. The potential activity distribution of CO 2 remains concentrated to the south throughout the first 90 days and all of phase 1, but is spread out more in longitude in the remainder of phase 1 than during the first 90 days. It is also interesting to note that whereas water is the dominant species throughout the period before the second equinox, even at the beginning of the mission when the comet was at heliocentric distances larger than 3.5 AU, after the second equinox CO 2 becomes the dominant species in the coma and the water production rate falls to nearly as small as that of CO. Except for the unusual behavior during the first 90 days O 2 generally follows the variation of H 2 O even as the water production drops at large heliocentric distances after the second equinox.

The production rates from the ROSINA-determined model in the time interval around perihelion do still seem to vary with spacecraft latitude, although the variations are not nearly as large those seen in the density itself as the sub-spacecraft latitude moved north and south. This is worth future investigation with a more complicated model than the fifth order spherical harmonics expansion used.

Finally, we estimated the total mass loss distribution of the comet at stages throughout the mission. Not surprisingly most of the mass was lost during the few months around perihelion when the production rates were the largest. A total volatile mass loss of (6.3±2.0) x 10 9 kg was found. Also, not surprisingly most of the surface loss throughout the apparition occurs in the southern hemisphere because of its large solar exposure during southern summer that occurs just after perihelion. We also used the total mass loss of the nucleus reported by the RSI team together with the total volatiles mass loss that implied a refractory/ice ratio somewhat less than one.

Figure 3 .

 3 Figure 3. Potential surface activity distribution for the species H 2 O, CO 2 , CO, and O 2 for the three time periods: (1) before the first equinox, (2) from the first to the second equinox, and (3) after the second equinox. The color scale gives the relative distribution of potential surface activity in each plot on a scale from 0 to 1.

Figure 5 .

 5 Figure 5. Production rates of H 2 O (blue), CO 2 (orange), CO (green), and O 2 (purple) throughout the mission as a function of heliocentric distance, R H in AU (top), and as a function of time (bottom). The vertical lines represent the first equinox, perihelion, and second equinox, respectively.

  of H 2 O production rates from DFMS, VIRTIS-H, MIRO, and PROCYON/LAICA and CO 2 production rates from DFMS and VIRTIS-HIn previous work(Fougere et al. 2016a(Fougere et al. & 2016b) ) we compared H 2 O and CO 2 column densities measured by the VIRTIS-H high-resolution infrared mass spectrometer[START_REF] Bockelée-Morvan | First observations of H2O and CO2 vapor in comet 67P/Churyumov-Gerasimenko made by VIRTIS onboard Rosetta[END_REF] & 2016) with those from the model calculations resulting from the ROSINA DFMS data inversion. We also calculated production rates resulting from both sets of data using the spatial distribution predicted by the model inversion. As discussed already, because the VIRTIS remote sensing measurements were made through columns that sampled the coma much closer to the nucleus, those data can give a more detailed picture of the small scale activity distribution from the surface than do the ROSINA measurements at the spacecraft because the variable density structures seen by VIRTIS tend to smear out at larger distances from the nucleus[START_REF] Combi | Narrow dust jets in a diffuse gas coma: A natural product of small active regions on comets[END_REF] where ROSINA sampled the coma. In addition, the effective FoV of ROSINA is much larger than the apparent size of the nucleus. Uncertainties for the DFMS model-determined production rates can be estimated from the uncertainties in the densities and the correlation coefficients of the model inversions. The uncertainties in the model inversion are uncorrelated with those of the densities, so the final uncertainties in the production rates of H 2 O, CO 2 , CO and O 2 are 30%, 33%, 33% and 37%, respectively.The VIRTIS-H H 2 O column densities in[START_REF] Bockelée-Morvan | Evolution of CO2, CH4, and OCS abundances relative to H2O in the coma of comet 67P around perihelion from Rosetta/VIRTIS-H observations[END_REF] were derived from the analysis of faint, optically thin, hot-band fluorescence emissions falling in the 4.45--5.0 µm range (order 0 of the spectrograph) where the rotational temperature was assumed. For this paper, new determinations were performed this time using more intense hot-band transitions detected between 2.775 and 2.910 µm in order 4 (see Fig.1of[START_REF] Bockelée-Morvan | Evolution of CO2, CH4, and OCS abundances relative to H2O in the coma of comet 67P around perihelion from Rosetta/VIRTIS-H observations[END_REF]).

  around perihelion the difference between DFMS and VIRTIS-H remains, but is somewhat smaller than found byFougere et al. (2016b) as well as by[START_REF] Hansen | Evolution of water production of 67P/Churyumov-Gerasimenko: An empirical model and a multi-instrument study[END_REF] who used the coma distributions fromFougere et al. (2016b).

Figure 6 .

 6 Figure 6. Comparison of production rates of CO 2 determined from the ROSINA DFMS data (orange circles) with those from VIRTIS-H (black x's) from the observations of Bockelée-Morvan et al. (2015 and 2016) using the spatial distributions from the DSMC model. The vertical lines left to right represent the first equinox, perihelion, and second equinox, respectively.

Figure 7 .

 7 Figure 7. Comparison of Rosetta water production rates in comet 67P/Churyumov-Gerasimenko. The red +'s are the same water production rates in Figure 6 obtained by interpreting the VIRTIS-H measurements of Bockelée-Morvan et al. (2015 and 2016) with the ROSINA DFMS determined water potential surface distribution model. The blue circles are also the same two-day average values from the ROSINA DFMS measurements. The black circles are two-day averages determined from the column densities from the MIRO observations of Marshall et al. (2017) also using the coma distributions from our ROSINA DFMS determined water surface distribution. The blue X's are determined from the H Lyα observations by Shinnaka et al. (2017), the sizes of which are indicative of the given uncertainties. The vertical lines left to right represent the first equinox, perihelion, and second equinox, respectively.

  [START_REF] Marshall | Spatially resolved evolution of the local H 2 O production rates of coet 67P/Churyumov,Gerasimenko, from the MIRO instrument of Rosetta[END_REF] and[START_REF] Shinnaka | Imaging Observations of the Hydrogen Coma of Comet 67P/Churyumov-Gerasimenko in 2015 September by the PROCYON/LAICA[END_REF] show comparisons of various water production rates including those determined from ground-based and near Earth satellite observations from previous apparitions going back to 1982. These represent a wide range of types of observations including H Lyα and OH in the radio and UV. They cover a considerable range in values but generally spanning the range of the values shown here from the 2014-2016 apparition.

Figure 8 .

 8 Figure 8. Comparison of potential surface activity distributions during the first 90 days with period 1, August 2014 -May 2015. The color scale gives the relative distribution in each plot on a scale from 0 to 1 of potential surface activity.

  present a model inversion of ROSINA DFMS density measurements of H 2 O, CO 2 , CO and O 2 , in the Rosetta target comet 67P/Churyumov-Gerasimenko throughout the mission which calculates the potential surface activity of these species in response to solar illumination as the nucleus rotates and moves in its orbit over the 25 months of the near nucleus phase of the mission. The potential surface activity for the four species was calculated during three time intervals of the orbit: (1) from the start of the near nucleus phase, 1 August 2014 to the first equinox; (2) from the first equinox through perihelion to the second equinox; and (3) from the second equinox until August 2016 when the DFMS measurements ended. Using the potential surface activity distributions determined, 15 sets of DSMC models were computed covering the two-year period, with each set consisting of 12 models spanning the ~12-hour rotation period of the nucleus. Interpolation and scaling of these 180 model calculations then provides a realistic description of the full variation of the major volatile species at any time during the mission. The results of the DSMC models were then shown to compare very favorably with the ROSINA DFMS densities, consisting of more than 50,000 individual measurements.The coma models were used to calculate the production rates of the four species throughout the mission. The coma model distributions were then used to analyze remote sensing measurements of H 2 O and CO 2 made by the VIRTIS-H infrared spectrometer and of H 2 O by the MIRO microwave instrument, and production rates of those species were compared with those determined from the model fit to the ROSINA DFMS densities. The discrepancy between the ROSINA DFMS production rates and those from VIRTIS-H obtained using different hot bands during the month around perihelion was improved somewhat, but remains and is still under study. However, the H 2 O production rates determined from nadir-pointing MIRO data and the DFMS-determined potential surface activity model were in much better agreement with those from DFMS throughout the mission, though the scatter in those results was larger.The MIRO observations used were viewing nadir columns directly to the nucleus often from near terminator locations. ROSINA was sampling the densities from similar spacecraft locations as MIRO. PROCYON was observing the entire large-scale hydrogen coma and should provide a reasonable global estimate of the production rate independent of pointing. VIRTIS was mostly viewing the active sector nearer the sub-solar point. However, ROSINA was sampling more directly the same regions of the nucleus but at other times when the illumination conditions were different. The difference in production rates around perihelion between the ROSINAdetermined model, the nadir-pointing MIRO and the PROCYON values determined from H-Lyα on the one hand from the lower values from our analysis using the same model on the VIRTIS columns still needs further investigation.The variation of production rates was shown to be different during the first 90 days of the mission than during the remainder of the pre-first-equinox time period. The water production rate decreases somewhat during this period but those of CO 2 and CO drop precipitously, and the O 2 production rate increases greatly. Using the same spherical harmonics inversion technique we looked at the potential surface activities of the four species during just the first 90 days (August 1 to October 31, 2014) and compared them with those throughout phase 1 from the beginning of the mission through the first equinox. The potential surface activity distributions of O 2 and water look similar to one another and look similar during the first 90 days as they do throughout phase

Table 1 .

 1 Parameters used in the analytical model to de-trend the input DFMS data for computing the surface activity distribution inversions. Best fit of potential surface activity distribution (black) compared with the ROSINA DFMS measurements of H 2 O (blue), CO 2 (orange), CO (green), and O 2 (purple) using just the simple analytical model ofBieler et al. (2015a). The left and right vertical lines give the times of the first and second equinox, respectively.

	Species	β-period 1	β-period 2	β-period 3	a night
	H 2 O	4.96	5.17	8.86	2 per cent
	CO 2	0.85	2.37	-0.83	2 per cent
	CO	-0.15	4.44	0.69	10 per cent
	O 2	5.71	4.94	7.17	2 per cent
	Figure 3 shows the potential surface activity distributions for the four major volatiles as a
	function of latitude and longitude on the nucleus of 67P. A number of interesting things can be
	seen by inspection. First the potential surface source distributions of H 2 O and O 2 are very similar
	for each of the three time periods. This further enhances the conclusions of Bieler et al. (2015b)
	that the O 2 is predominantly associated with H 2 O ice and despite the high volatility of frozen O 2 ,
	the release of O 2 in comets is mostly controlled by H 2 O.
	In contrast that of CO 2 has most of its dominant activity in the south with little in the
	north compared with H 2 O and O 2 . This is consistent with the first reported DFMS results by
	Hässig et al. (2015) and with VIRTIS results first reported by Bockelée-Morvan et al. (2015). In
	addition the potential surface activity distributions are quite similar throughout the whole
	apparition, unlike those for the other species that change noticeably over the three studied
	periods. The distributions of CO and CO 2 are very similar for the second and third time periods,
	that is, after the first equinox, but before the first equinox the distribution of the CO surface
	source is quite different. CO has a major source in the north coinciding with H 2 O as well as
	activity in the south similar to but not quite as strong as CO 2 .

Table 3 .

 3 Real spherical harmonics and the normalized coefficients for time period Phase 2.

	Term Y 0	Expression 1	H 2 O coefficient 4.895834e-01	CO 2 coefficient 1.375611e+00	CO coefficient 4.609846e-01	O 2 coefficient 2.407584e-01
	Y -1	sin(ϕ) sin(θ)	-3.628282e-01	-4.683337e-02	2.602201e-01	-1.642184e-02
	Y 0	cos(θ)	-2.644014e-01	-2.335596e+00	-5.886954e-01	-1.625394e-01
	Y 1	cos(ϕ) sin(θ)	5.399610e-01	-4.373107e-01	-1.177043e-01	5.914185e-02
	Y -2	sin(2ϕ) sin2(θ)	-8.288829e-02	9.081417e-02	6.886639e-02	2.291460e-02
	Y -1	sin(ϕ) sin(θ) cos(θ)	-1.459630e-02	-3.300306e-01	-2.445945e-01	-4.651805e-03
	Y 0	3 cos 2 (θ) -1	3.385707e-01	1.239087e+00	2.656085e-01	1.044153e-01
	Y 1	cos(ϕ) sin(θ) cos(θ)	-3.197988e-02	-7.304284e-02	3.411257e-02	-8.059736e-04
	Y 2	cos(2ϕ) sin 2 (θ)	-1.070801e-01	-3.457324e-01	-9.687787e-02	-1.179830e-01
	Y -3	sin(3ϕ) sin 3 (θ)	-7.809985e-04	-2.156722e-01	-1.349597e-01	3

Table 5 .

 5 Geometries for the DSMC performed simulations. For each of the 15 dates in the table 12 hourly simulations were run covering one nucleus rotation each.

	Date	R au (au)	Solar Latitude ( o )
	2014-08-23	3.50	43.6
	2014-12-23	2.71	32.6
	2015-03-04	2.18	20.3
	2015-05-06	1.70	1.7
	2015-06-08	1.48	-12.9
	2015-07-01	1.35	-25.3
	2015-07-22	1.27	-37.3
	2015-08-16	1.24	-49.0
	2015-09-25	1.35	-49.3
	2015-10-28	1.54	-37.8
	2015-12-01	1.78	-25.8
	2016-01-14	2.12	-13.4
	2016-03-30	2.69	1.30
	2016-06-06	3.15	10.3
	2016-09-25	3.80	20.0

Table 6 .

 6 Ratios of H 2 O and CO 2 column densities from VIRTIS-H spectra and the DFMS fitted model

	Species and Time Period	VIRTIS-H/DFMS Model Ratio
	H 2 O Early Mission	0.69
	CO 2 Early Mission	0.52
	H 2 O Perihelion Period	0.34
	CO 2 Perihelion Period	0.90

Table 7 .

 7 Mass loss in kg from gas sublimation of H 2 O, CO 2 , CO and O 2 .

	Species	Before 1st Equinox	Before 2nd Equinox	End of Mission
	H 2 O CO 2 CO	2.2 x 10 8 2.3 x 10 7 1.2 x 10 7	4.8 x 10 9 8.6 x 10 8 2.0 x 10 8	4.9 x 10 9 8.8 x 10 8 2.1 x 10 8
	O 2

Table A1 .

 A1 Revised VIRTIS-H column densities based on the 2.775-2.910 µm hot bands Heliocentric distance [au] Column #4 LOS mean distance to the nuclei center [km] Column #5 mean LOS position angle relative to the Sun (PA(Sun)=27 o ) [deg] Column #6 H 2 O column density from 4.450-5.030 micron [10 20 m -2 ]; uncertainty in parentheses. Column #7 H 2 O column density from 2.775-2.910 micron [10 20 m -2 ]; uncertainty in parentheses. Column #8 Rotational temperature from H 2 O hot band [K] Column #9 Ratio of column densities from 2.775-2.910 micron hot band to from 4.450-5.030 micron hot band *** More detail on observing/geometric conditions please read Bockelée-Morvan et al. (2016)

	Obs ID	Start Time	r h	ρ	P.A.	N 4µm (H 2 O)	N 2µm (H 2 O)	T (rot)	ratio
	T1_00395011055	2015-07-08T21:11:00.3	1.316 2.82	259	9.25 (2.88) 10.96 (0.18) 128.7 (2.6) 1.12
	T1_00395322550	2015-07-12T11:46:52.6	1.302 2.70	274	10.10 (0.77)	9.80 (0.18)	125.7 (3.2) 0.95
	T1_00395742154	2015-07-17T08:15:55.9	1.286 3.10	272	6.33 (0.36)	7.48 (0.13)	120.7 (2.4) 1.14
	T1_00396199623	2015-07-22T15:23:57.9	1.271 2.69	266	10.00 (0.41) 12.25 (0.20) 141.8 (2.9) 1.19
	T1_00396220410	2015-07-22T21:10:24.7	1.270 2.98	266	9.79 (0.28) 12.38 (0.17) 126.4 (2.2) 1.21
	T1_00396659230	2015-07-27T23:06:40.5	1.259 4.01	247	6.48 (0.17)	6.14 (0.12)	111.4 (3.1) 0.92
	T1_00396826054	2015-07-29T21:20:59.8	1.255 3.56	265	6.62 (0.44)	9.53 (0.16)	125.9 (2.4) 1.38
	T1_00396842044	2015-07-30T01:47:25.8	1.255 3.59	265	7.27 (0.48)	9.82 (0.16)	110.7 (2.4) 1.31
	T1_00397038715	2015-08-01T08:25:16.8	1.252 3.08	270	7.10 (0.39) 10.28 (0.19) 131.5 (2.7) 1.41
	T1_00397871165	2015-08-10T23:32:29.0	1.244 6.23	286	4.19 (0.40)	7.35 (0.12)	93.8 (3.3)	1.72
	T1_00398347798	2015-08-16T11:56:25.1	1.244 3.58	267	9.21 (0.63) 12.70 (0.22) 121.0 (2.6) 1.36
	T1_00398603680	2015-08-19T11:08:06.5	1.246 2.96	263	10.20 (0.50) 12.82 (0.23) 136.3 (2.6) 1.25
	T1_00398619671	2015-08-19T15:34:33.4	1.246 2.97	262	11.60 (0.64) 12.42 (0.18) 129.3 (2.1) 1.06
	T1_00398640452	2015-08-19T21:20:58.5	1.246 3.47	262	8.90 (0.61) 11.27 (0.20) 138.9 (2.6) 1.24
	T1_00398729980	2015-08-20T22:13:06.6	1.247 4.58	258	8.32 (0.38) 11.00 (0.20) 113.6 (2.6) 1.31
	T1_00398970868	2015-08-23T17:07:54.5	1.250 3.78	256	8.27 (0.33) 11.13 (0.20) 118.7 (2.8) 1.31
	T1_00398986563	2015-08-23T21:29:29.5	1.250 3.82	255	7.44 (0.36) 12.55 (0.18) 119.9 (2.2) 1.65
	T1_00399198695	2015-08-26T08:17:52.6	1.254 4.25	264	6.87 (0.36)	9.47 (0.20)	116.3 (3.1) 1.32
	T1_00399208316	2015-08-26T10:58:13.4	1.254 3.95	263	9.30 (0.99) 11.43 (0.17) 112.4 (2.5) 1.21
	T1_00400195914	2015-09-06T21:25:16.9	1.280 4.92	266	5.25 (0.49)	7.03 (0.11)	111.6 (2.5) 1.31
	T1_00400216374	2015-09-07T03:06:16.9	1.281 4.87	268	5.64 (0.43)	7.82 (0.12)	112.1 (2.6) 1.34
	T1_00400232094	2015-09-07T07:28:16.8	1.282 4.83	268	3.64 (0.60)	7.94 (0.12)	106.6 (2.5) 2.12
	T1_00400252549	2015-09-07T13:09:16.0	1.282 4.83	270	5.35 (0.35)	8.42 (0.11)	107.1 (2.5) 1.54
	T1_00400407713	2015-09-09T08:15:16.0	1.288 3.14	263	8.22 (0.49) 12.50 (0.19) 125.9 (2.3) 1.46
	T1_00400433767	2015-09-09T15:29:34.0	1.289 3.72	262	7.42 (0.30) 10.83 (0.16) 109.8 (2.0) 1.40
	T1_00400470243	2015-09-10T01:37:26.0	1.290 5.78	263	4.93 (0.18)	6.99 (0.13)	89.2 (3.1)	1.39
	T1_00400748481	2015-09-13T06:54:44.2	1.301 3.58	265	7.13 (0.29) 10.07 (0.12) 119.8 (2.1) 1.35
	T1_00400765976	2015-09-13T11:46:19.3	1.302 3.61	265	10.00 (1.00) 12.22 (0.20) 124.5 (2.7) 1.20
	T1_00401717336	2015-09-24T11:55:21.0	1.347 5.02	275	3.68 (0.19)	7.12 (0.17)	98.5 (4.0)	1.88
	T1_00401994267	2015-09-27T16:57:54.5	1.362 6.88	282	2.41 (0.18)	3.31 (0.09)	77.1 (6.9)	1.34
	Notes to Table A1
	Column #1 Observation ID
	Column #2 Start time of observation
	Column #3

Acknowledgements

This work was supported by contracts JPL1266313 and JPL 1266314 from the US Rosetta Project and NASA grants NNX14AG84G from the Planetary Atmospheres Program and 80NSSC18K1280 from the Rosetta Data Analysis Program. Work at UoB was funded by the State of Bern, the Swiss National Science Foundation (200021_165869, 200020_182418), and the European Space Agency PRODEX Program. The authors would like to thank ASI, Italy;


F Raponi

G Capaccioni

T Filacchione

C Suhasaria

D Popa

L Kappel

V Moroz

Vinogradoff

Vito Mennella
email: vito.mennella@inaf.it

M Ciarniello

A Raponi

F Capaccioni

G Filacchione

T Suhasaria

C Popa

D Kappel

L Moroz

V Vinogradoff

A Pommerol

B Rousseau

I Istiqomah

E Quirico

P Beck

D Bockelee-Morvan

R W Carlson

C Pilorget

Hydroxylated Mg-rich Amorphous Silicates: A New

Keywords: astrochemistry -comets, general -methods, laboratory, atomic processestechniques, spectroscopic

published or not. The documents may come

INTRODUCTION

The payload instruments on-board ESA's Rosetta mission have observed the evolution of the Jupiterfamily comet 67P/Churymov-Gerasimenko (67P/CG) in a period from August 2014 to September 2016 and encompassing the perihelion passage occurring in August 2015 at a heliocentric distance of 1.24 AU. The data collected have provided unprecedented spatial, spectral and temporal coverage of the changes occurring on a cometary surface and surrounding inner coma before, during and after perihelion passage (Taylor et al. 2017). Determining the nucleus surface composition was one of the primary scientific objectives of the Visible, Infrared and Thermal Imaging Spectrometer (VIRTIS) (Coradini et al. 2007).

When observed at global scale, the 67P/CG surface appears dominated by a dust layer, remarkably uniform in composition. The overall spectral reflectance is characterized by a low albedo (0.06) and red slope at visible wavelengths while a wide absorption band centered at 3.2 µm is a distinctive feature in the infrared range (Capaccioni et al. 2015;Ciarniello et al. 2015). Temporal and spatial spectral variability at global scale has been observed on the 67P/CG surface (Ciarniello et al. 2016) in addition to spatial variations on meter-sized scales where color changes have been correlated with surface composition inhomogeneities and active sources (Filacchione et al. 2019). In particular, local water ice enrichments show higher albedo and flat-to-blue spectral response in contrast to the more dark and red reflectance typical of the other areas. The lack of the 1.5 and 2.0 µm ice absorption bands in the typical dark areas indicates a very dehydrated surface with the water ice abundance <1% on average (Raponi et al. 2016). Although the surface appears dehydrated at global scales, crystalline water ice has been spectroscopically identified at infrared wavelengths on specific areas of the nucleus during the pre-perihelion phase (De Sanctis et al. 2015;Filacchione et al. 2016a;Barucci et al. 2016). Also CO 2 ice has been detected for the first time on a cometary nucleus surface by VIRTIS in a small area located in the meridional Anhur region (Filacchione et al. 2016b).

The most striking spectral feature is the broad band between 2.8 and 3.6 µm centered at 3.2 µm, whose attribution is challenging. The band falls in a spectral region where the contribution to the absorption of X-H stretching vibrations in different functional groups is expected. C-H vibrations in aliphatic and aromatic C bonding have recently been established as sub-structures of the observed absorption band in newly calibrated VIRTIS spectra (Raponi et al. 2020). O-H vibrational modes in water ice, structural water or chemisorbed hydroxyl groups absorb from 2.7 to 3.1 µm but cannot explain the entire feature. Quirico et al. (2016) suggested several possibilities, OH groups in carboxylic acids inserted in a macromolecular organic solid and NH 4 + (ammonium) ions were proposed as plausible candidates to encompass the broad width of the cometary feature. More recently on the basis of laboratory experiments, the 3.2 µm band and some of its sub-features have been interpreted to be ammonium (NH 4 + ) salts (such as ammonium formate, sulfate, citrate) mixed with dark cometary surface material (Poch et al. 2020).

Silicates, together with powders of coal and iron sulphides, have been used to model the red visible and near infrared slope and low albedo of 67P/CG (Rousseau et al. 2018). Detection in the coma of a mineral phase, predominantly composed of silicates and carbonaceous components, was carried out by the COSIMA experiment aboard Rosetta (Bardyn et al. 2017). The presence of silicates in comets was first reported by Bregman et al. (1987) and later confirmed by Campins & Ryan (1989) in comet 1/Halley and by Crovisier et al. (1997) in Hale-Bopp (C/1995 O1). During their journey from the birth sites around evolved stars through interstellar medium to comets, silicates interact with ions (mainly protons), which induce many effects such as amorphization, sputtering and implantation with the formation of OH groups (Demyk et al. 2001(Demyk et al. , 2004;;Carrez et al. 2002;Jäger et al. 2003a;Brucato et al. 2004;Djouadi et al. 2005). The abundance of amorphous silicates in ISM is 98 % with respect to the crystalline counterpart, whereas the abundance varies between 70-98 % for Jupiter family comets (Sugita et al. 2005). Furthermore, the silicates are primarily Mg-rich e.g. Kelley et al. (2017). Direct evidence of silicate processing by ions has come from analysis of Interplanetary Dust Particles (IDPs) (Bradley 1994;Bradley et al. 1999Bradley et al. , 2005)). During their journey, silicates can also interact with hydrogen atoms both in diffuse and dense regions of the ISM.

In this paper, we experimentally investigate the consequence of the interaction of H atoms with analogues of interstellar silicate grains. Hydrogen atoms, with an energy of a few tens of meV, form hydroxyl groups in Mg-rich amorphous silicates with different Mg/Si ratios. The resulting IR spectrum is characterized by a broad absorption band centered at 3.2 µm similar to the one observed in 67P/CG spectrum.

EXPERIMENTS AND RESULTS

Hydroxylation of silicates has been reproduced in controlled laboratory conditions. The amorphous silicates were produced using a Nd-YAG solid state pulsed laser at 1024 nm equipped with II and IV harmonics at 532 and 266 nm, respectively. By using a focusing lens, the power density on the target was measured by LaserMate-P power meter (Coherent) to be 1×10 8 W cm -2 . The target to be ablated was mounted inside a chamber, filled with a 10 mbar Ar atmosphere, on a rotating sample holder placed in front of the laser beam. Silicate grains were collected at a distance of 5 cm from the target on ZnSe windows and carbon stubs, respectively, for IR spectroscopy and electron microscopy. We have used pellets made out of oxide mixtures in two different proportions (MgO:SiO 2 in a 1:1 and 2:1 mixture) and a natural forsterite olivine (dunite, Webster-Addie, N. Carolina, acquired from Wards Science) as targets. The oxide mixtures were prepared by weighing the exact amounts of MgO and SiO 2 in order to have the starting stoichiometric composition of enstatite (MgSiO 3 ) for the first sample (SIL1) and forsterite (Mg 2 SiO 4 ) for the second sample (SIL2). The mixtures were produced in the 1. The chemical composition resulting from the process of vaporization and condensation was analysed with an Electron Dispersive X-ray (EDX) system (Oxford Inca Energy 350) linked to a Field Emission Scanning Electron Microscope (FESEM) (ZEISS SUPRA 25). The FESEM images show that mainly a fluffy structure of amorphous grain aggregates characterizes the morphology of the silicate samples. For a detailed characterization of texture, grain dimension and size distribution of the samples see Rotundi et al. (2000). The atomic ratio of the samples obtained by EDX is also reported in Table 1. Taking into account the variations due to non-homogeneity of the samples, the Mg/Si ratio is in agreement with the corresponding value of the starting oxide mixtures for SIL1, while it is slightly lower in the case of SIL2. The values, Fe/Mg = 0.11 and (Mg+Fe)/Si = 1.9, indicate that the composition of the sample obtained from natural olivine (SIL3) is Mg-rich.

Figure 1 shows the optical depth spectra of the silicates acquired with a Fourier Transform Infrared (FTIR) spectrometer (Bruker Equinox 55). Spectra acquisition was performed in air at room temperature, in the range 2-16 µm at a spectral resolution of 2 cm -1 . The spectra are characterized by a strong band around 10 µm due to the Si-O stretching mode, whose profile slightly varies with the sample composition. The band is centered at 9.8 µm for SIL1 and it shifts to 10.2 µm as the Mg/Si ratio increases (for SIL2 and SIL3 samples). In the latter cases, an additional band at around 10.8 µm is also evident. Jäger et al. (2003b) also reported a similar redshift in the band position with increasing MgO/SiO 2 ratio in their amorphous magnesium silicate samples, prepared by the sol-gel method. The absorption increase for wavelengths longer than 13 µm is due 1. The spectra of SIL2 and SIL3 are offset in ordinate by 0.3 and 0.5, respectively.

to the low energy wing of the Si-O bending feature at around 20 µm. The small double band at 6.5 and 7 µm, more intense in the case of SIL2, is caused by traces of magnesium carbonate, formed during sample preparation (Jäger et al. 2003b). Figure 2 shows the silicate infrared spectra in the 3 µm spectral range after continuum removal. A broad band centered at 2.88 µm, produced by adsorbed water and chemisorbed hydroxyl groups, due to ambient contamination during the sample handling at the time of pellet preparation and also at the time of laser ablation, characterizes this region. Two narrow bands at 2.67 µm and 2.69 µm are respectively attributed to isolated Si-OH and Mg-OH groups [START_REF] Iler | The chemistry of silica (John Wiley and Sons Inc) Jäger[END_REF]Shinoda et al. 2002). Their relative intensity varies according to the Mg/Si ratio. In the case of SIL3, produced from natural olivine, only a weak band, resulting from the blend of the two features, is observed. Similar bands have been detected in silicates produced by the sol-gel method (Jäger et al. 2003b), suggesting that laser ablation and sol-gel techniques provide materials with chemical properties expected for cosmic silicates. In fact, the formation of silicates in space with a small content of isolated Si-OH and Mg-OH bonds can be considered highly probable, since H 2 O, the most abundant oxygen-bearing molecule, participates in the formation of olivine and pyroxene-type silicates around evolved stars [START_REF] Gail | Dust formation in M stars[END_REF]Jäger et al. 2003b). To study the interaction of hydrogen with silicates, our samples were placed inside a vacuum chamber, with a base pressure of 2×10 -8 mbar, on a rotating sample holder at the right angle intersection between an atomic hydrogen beam and the IR beam of an FTIR spectrophotometer (Bruker Vertex 80). H atoms were produced by microwave dissociation of molecular hydrogen; they have a Maxwellian velocity distribution at 300 K (Mennella 2006). After H atom exposure, the sample was rotated by 90 • to face the IR beam, and we measured the IR spectrum with a resolution of 2 cm -1 in the range 2-5 µm (5000-2000 cm -1 ). For any spectral measurement of the silicate before and after irradiation, a background was acquired. Figure 3 shows the spectral changes, obtained by subtracting the pre-irradiated spectrum from the irradiated one, induced by H atom irradiation in amorphous silicates. The activation of a broad band centered at about 3.2 µm takes place in all the silicates. At shorter wavelengths, an absorbance decrease in the irradiated sample (resulting in a negative optical depth difference in Figure 3) centered at 2.76 µm (broad band) and bands at 2.67 and 2.69 µm (sharp doublet) are observed. This behavior indicates that the H atoms induced a depletion of the O-H bonds already present in the starting silicate samples.

In addition, a weak feature at 3.36 µm with a shoulder at 3.39 µm is possibly due to asymmetric C-H stretching modes of aliphatic CH 3 and CH 2 groups. They originated from processing of carbonate, also present as substructure of the broad 3.2 µm band. The faint feature at around 4.25 µm is associated to uncompensated ambient CO 2 . Irradiation of a blank ZnSe substrate with hydrogen (see Figure 3) does not induce any variation in the spectrum, demonstrating that the observed spectral changes are due to the interaction of hydrogen with silicate grains.

To investigate whether the presence of OH groups in the pristine silicates affects their interaction with hydrogen, we tried to eliminate/reduce them by annealing one of the samples (SIL3) at 300 • C for 1 hour before exposure to H atoms. In the case of the annealed SIL3, only a small absorbance decrease centered at 2.76 µm is present after hydrogen irradiation (see Figure 3), testifying that the starting sample was almost fully depleted of the chemisorbed hydroxyl groups. Despite that, the activation of the 3.2 µm band does not vary. Hydroxylation of Mg-rich amorphous silicates by H atoms is therefore independent of the OH groups already present in silicates and of the specific composition of the samples.

Figure 3. Spectral changes of silicates upon irradiation (shown after subtraction of the initial spectrum) of 1×10 18 H atoms cm -2 . The spectrum of an irradiated ZnSe window (blank) after an irradiation is also reported for comparison. The feature at around 4.25 µm is associated to uncompensated ambient CO2. Note that vertical is inverted and therefore absorbance increases downwards.

DISCUSSION AND CONCLUSIONS

The spectral variations induced in Mg-rich silicates by atomic hydrogen are clear evidence for chemical reac-tions producing hydroxyl groups on silicates. According to density functional quantum mechanical calculations, the formation of OH groups in Mg-rich silicate nanoclusters by impinging H atoms takes place without any significant energy barrier (Kerkeni et al. 2017). Moreover, the observed spectral evolution with H atom irradiation suggests that exchange reactions (abstraction and addition of H atoms in OH) are active. The equilibrium value for the number of OH groups depends on the cross section of the abstraction and addition processes. Consequently, an absorbance decrease as in the region of the bands at 2.76 µm (broad band) and 2.67 and 2.69 µm (sharp doublet) indicates a net reduction in the number of available OH bonds due to abstraction of H atoms already bonded to these groups. At longer wavelengths the exchange reactions induced by H atoms determines, as a net result, the formation of an absorption band; this is the case for the observed broad 3.2 µm band.

Formation of OH bonds with a peak at around 2.8 µm has been reported after implantation of keV to MeV energy range protons in silicates with different composition and structure (e.g. Zeller et al. (1966); Djouadi et al. (2011); Schaible & Baragiola (2014)). Exchange of H/D in hydroxyl groups has also been observed after 1.1 keV D + implantation of Apollo lunar soils [START_REF] Ichimura | [END_REF]. Hydroxylated silicates produced by proton implantation have also been considered as a potential reservoir of water during the formation of terrestrial planets (Djouadi et al. 2011).

There are two substantial differences between the present work and previous implantation experiments. The first one concerns the interaction energy. Unlike proton implantation in silicates, where ions with high energy have been considered, the present experiments show the possibility to hydroxylize amorphous silicates via reactions with thermal (a few tens of meV) hydrogen atoms. This result extends the possibility to hydroxylate silicates in different environments: from the interstellar medium (ISM), where atomic hydrogen (10-100 K) is abundant, to the surface of airless planetary bodies (solar protons) such as the Moon.

The second difference is in the width of the band around 3.2 µm. In the present experiments, it is much wider than that obtained by proton implantation due to multiple and independent overlapping absorption of hydroxyl groups. These groups are formed by hydrogen atom interaction on amorphous silicate defects, which determine a wider range of binding sites enhancing both physisorption and chemisorption of hydrogen.

The profile of the band formed by H atoms in the considered silicates (see Figure 3) is reminiscent of the 3.2 µm absorption band in the reflectance spectrum of 67P/CG observed by VIRTIS onboard Rosetta (Raponi et al. 2020). Figure 4 shows the comparison of the cometary profile with that of annealed SIL3. The similarity in profile suggests that OH groups of hydroxylated amorphous silicates can be a constituent of the 67P/CG nucleus surface, together with aliphatic carbons, ammonium salts, water ice and possibly PAHs, Polyphenols (Quirico et al. 2016;Carlson et al. 2020;Poch et al. 2020;Raponi et al. 2020).

To quantitatively evaluate the spectral contribution of hydroxylated amorphous silicate to the 3.2 µm band, we have modeled the average reflectance (I/F) spectrum of the comet 67P/CG surface with a similar approach already used to estimate the spectral contribution of water ice (Raponi et al. 2020). The Hapke radiative transfer model [START_REF] Hapke | Theory of reflectance and emittance spectroscopy[END_REF]) was adopted to calculate the reflectance of an intimate mixture of hydroxylated amorphous silicate with a synthetic continuum representing a dark-red terrain. Model fits to the average 67P/CG reflectance spectrum indicate that hydroxylated silicate can reproduce on average 65% of the 3.2 µm cometary band intensity with an abundance of hydroxylated amorphous silicate of 15%. An example of the resulting modeling is shown in Figure 5. For more detail see Table 2 and Figure 9 in the Appendix. The estimated contribution should be considered as an upper limit due to the involved assumptions (see Appendix). It is worth stressing that for an accurate quantitative estimation, one should also consider the contributions of the other proposed carriers of the 3.2 µm band as mentioned above.

The presence of a hydroxylated fraction in silicates on the surface of a comet can represent an evolutionary link between a primitive object of the solar system and dust in the ISM. The spectral homogeneity of the nucleus surface observed by VIRTIS indicates that the composition of the refractory surface materials does not change despite the rejuvenating processes active during the perihelion passages of 67P/CG (Ciarniello et al. 2016). The refractory material seen at the surface must be then representative of the bulk pristine material of the nucleus.

As for all cosmic dust, silicates participate in the life cycle of matter that brings dust from their circumstellar birth sites to the diffuse interstellar medium and dense clouds where stars form. The formation of planetary systems, like our solar system, can also take place during this cycle. Since silicates are involved in the global evolution of the interstellar matter, from diffuse medium to comets, they interact with hydrogen atoms in the ISM. The laboratory H atom fluence of 1×10 18 atoms cm -2 used in the present experiments as reported in Figure 3 corresponds to an H atom exposure time of 4×10 3 and 4×10 5 years, respectively, in interstellar diffuse and dense regions (considering the respective H atom flux of 8×10 6 and 9×10 4 atoms s -1 cm -2 for these regions, Mennella (2010)). In addition, we observed the formation of the 3.2 µm band already at a lower fluence of ca. 10 16 atoms cm -2 . Therefore, hydroxylation can take place on a much shorter time scale compared to the residence time of 3×10 7 years of dust grains in both dense and diffuse regions of the ISM. However, for a complete understanding of the degree of hydroxylation of the silicates in different regions of the ISM, we also need to know the effects of the dehydrogenation processes induced by energetic processing such as by UV photons and cosmic rays.

The evolutionary link suggested by the present results is consistent with the detection of the aliphatic matter in 67P/CG (Raponi et al. 2020). We thus can trace for both aliphatic carbon and hydroxylated amorphous silicates a common origin in the ISM during the journey of the dust from its birth-sites around evolved stars to our solar system. phous silicate material, synthetic optical constants have been computed. Upon H atom irradiation, the optical properties of the amorphous silicates are modified by formation of the absorption feature centered at 3.2 µm. This effect can be interpreted as an increase of the value of the imaginary part of the refractive index k, ∆k, in correspondence to the 3.2 µm absorption feature. To estimate ∆k, we represent the contribution of OH groups formed on the fluffy and amorphous grain aggregates in consideration by assuming a homogenous film of hydroxylated amorphous silicate of equivalent thickness d. Then the measured optical depth difference (irradiated -pristine) can be expressed as:

∆τ (λ) = 4π λ d ∆k. (A1)
The equivalent thickness d can be derived from the integrated optical depth through the following relation:

band ∆τ (ν) dν = N A = ρ d N A m A, ( A2
)
where N is the column density, A is the integrated band strength and N A is Avogadro's number. ρ and m are the density and the molecular weight of the film, respectively. We evaluated the integral in the left hand term of Eq. A2 from the annealed SIL3 reported in Figure 6 and adopt A = 2×10 -17 cm molecule -1 , as in the case of the OH groups in natural olivine (Schaible & Baragiola 2014). Moreover, we assumed ρ = 1 g cm -3 and m = 101 g mol -1 for hydroxylated silicate film with a starting silicate stoichiometry as [MgO.SiO 2 ]. We obtain a starting value of d = 0.35 µm. It is then possible to derive the corresponding increase of the imaginary part of the refractive index from Eq. A1 as seen in Figure 7. To see how the final results are affected by our assumptions, we also varied d by a factor of 10 higher and lower, reflecting possible uncertainties in the assumptions on A, ρ and m.

Since the H atom processing of Mg-rich amorphous silicates does not lead to variations of the absorption spectrum outside the 3.2 µm band, it is reasonable to assume the imaginary part of the refractive index k b of the unprocessed silicate as a baseline. However, optical constant values of the unprocessed silicates analyzed in the present work are missing. Therefore, we adopted the optical constants of similar silicates. The choice is further justified by the fact that H atom irradiation modifies different Mg-rich amorphous silicates in a similar way, as also evidenced by the reported measurements. We considered real and imaginary parts of optical constants of Mg-rich glass from Dorschner et al. (1995) For each synthetic refractive index values, a fit to the 67P/CG spectrum is performed in a similar way as Raponi et al. (2020) did for water ice. We forced the fit to match the continuum level of the cometary spectrum in the range between 2.7-2.9 µm and 3.5-3.6 µm to ensure that the modeled absorption is not more intense than the one shown by 67P/CG. The best fits to the 67P/CG spectrum are reported in Figure 9 while the best fit parameters i.e. the amount of hydroxylated silicates with respect to dark-red terrain, R, and the corresponding grain size, S, are listed in Table 2. Note that R values reported in Table 2 represent the corresponding fractional cross section seen by the light propagating within the medium. These quantities can be interpreted as volumetric abundances under the assumption that the featureless dark-red terrain has the same grain size as the hydroxylated silicate. The derived R values are constant at 0.15, and F is very similar for the three cases with an average value of around 0.65. Conversely, the obtained S values increase with d. This can be explained by the fact that ∆k as obtained from Eq. A1 is lower for larger d, thus producing a smaller absorption for a given grain size. Given this, in order to reproduce the observed spectral contrast, a larger value of S is required in the modeling. One should note that a rigorous application of the Hapke's model would limit it's application to the geometric optics regime, where particles are larger than wavelength. This is the case for the derived value of S = 0.5 µm which implies that small particle sizes (likely sub-µm) are required to match the observed spectrum with the adopted set of optical constants. . From top to bottom: best fits (blue dashed lines) to average 67P/CG reflectance spectrum (black lines) obtained with synthetic optical constants computed using Mgrich glass (Dorschner et al. 1995) as baseline for different values of d.

Figure 1 .

 1 Figure1. Infrared spectra of the silicates reported in Table1. The spectra of SIL2 and SIL3 are offset in ordinate by 0.3 and 0.5, respectively.

Figure 2 .

 2 Figure2. Infrared spectra, after continuum removal, of the silicates reported in Table1in the 2.5 -3.3 µm range. The spectra of SIL2 and SIL1 are offset in ordinate by 0.006 and 0.012, respectively.

Figure 4 .

 4 Figure 4. Qualitative comparison of 3.2 µm absorption band in SIL3 (blue line) first annealed at 300 • C and then irradiated with 1×10 18 H atoms cm -2 to the observed absorption band in the 67P/CG reflectance spectrum (black line). The left vertical axis is for the SIL3 spectrum whereas the right vertical axis refers to that of the cometary band.

Figure 5 .

 5 Figure 5. Fit to the average 67P/CG reflectance spectrum (black line, Raponi et al. (2020)) by the modeled hydroxylated amorphous silicate mixed with the dark-red endmember (blue dashed line) using Hapke's model for 5 µm grain size and 15% of hydroxylated amorphous silicate abundance. See the Appendix for the meaning of d and other parameters.

  optical constants, a synthetic imaginary part of the refractive index k s = k b +∆k is derived as shown in Figure8, for the Mg-rich glass for different values of d.

Figure 7 .

 7 Figure 7. ∆k as a function of wavelength for d = 0.35 µm.

Figure 8 .

 8 Figure 8. Imaginary part of the synthetic refractive index derived by using Mg-rich glass as baseline. The three curves represent the different values of ks obtained for d values of 0.035, 0.35 and 3.5 µm from top to bottom, respectively.

Figure 9

 9 Figure9. From top to bottom: best fits (blue dashed lines) to average 67P/CG reflectance spectrum (black lines) obtained with synthetic optical constants computed using Mgrich glass(Dorschner et al. 1995) as baseline for different values of d.

Table 1 .

 1 Atomic ratios of silicates produced by laser ablation period

	Sample name Starting material Measured ratios
			Mg/Si	Fe/Si
	SIL1	MgO:SiO2 (1:1)	1.1	-
	SIL2	MgO:SiO2 (2:1)	1.7	-
	SIL3	Natural Olivine	1.7	0.2
	form of pellets (13 mm in diameter and a few millime-
	ters in thickness) by pressing finely ground oxides at 10
	tons. The third sample (SIL3) was prepared from nat-
	ural olivine target. The three silicates (SIL1, SIL2 and
	SIL3) prepared in the present work are listed in Table

  Table 1 in the 2.5 -3.3 µm range. The spectra of SIL2 and SIL1 are offset in ordinate by 0.006 and 0.012, respectively.

  . Further, to incorporate the effects of H atom processing on the Optical depth difference (after continuum removal) for annealed SIL3. As mentioned in the text, the H atom irradiation process abstracts H atoms in OH groups already present in the sample, introducing a spurious negative signature in the ∆τ spectrum between 2.65 and 2.85 µm. To eliminate this spurious signature, we discard the measurement in this range and replace it by a polynomial (blue segment) that connects the part of the optical depth difference spectrum shortward of 2.65 µm with the region longward of 2.85 µm.

	0.008
	0.006
	0.004
	Dt
	0.002
	0.000
	2.0 -0.002	2.5	3.0	3.5	4.0	4.5
			wavelength (µm)
	Figure 6.

Table 2 .

 2 Best fit parameters (R, S) of the model to the average reflectance spectrum of 67P/CG. Equivalent film Hydroxylated silicate Grain size Fraction of the 3.2 µm band intensity

	thickness	abundance		due to hydroxylated silicates
	d	R	S	F
	µm		µm
	0.035	0.15	0.5	0.62
	0.35	0.15	5	0.66
	3.5	0.15	50	0.67

ACKNOWLEDGMENTS

This work takes advantage of the collaboration of the ISSI international team "Comet 67P/Churyumov-Gerasimenko Surface Composition as a Playground for Radiative Transfer Modelling and Laboratory Measurements", number 397.

APPENDIX

A. MODELING OF THE SPECTRAL CONTRIBUTION OF HYDROXYLATED SILICATES TO THE 3.2 µm FEATURE OF 67P/CG

We apply the Hapke model [START_REF] Hapke | Theory of reflectance and emittance spectroscopy[END_REF] to derive an upper limit of the amount of hydroxylated Mg-rich amorphous silicates contributing to the 3.2 µm band in-tensity of 67P/CG (Raponi et al. 2020). To this aim, 67P/CG's surface is modeled as an intimate mixture of a featureless dark-red synthetic component, representing predominantly the optically opaque phases in the comet regolith, and a hydroxylated Mg-rich amorphous silicate. The spectral output of the featureless end-member is modeled as a continuum by means of a 2nd degree polynomial, while for the hydroxylated Mg-rich amor-