Estimate of the erosion rate from H2O mass-loss measurements from SWAN/SOHO in previous perihelions of comet 67P/Churyumov-Gerasimenko and connection with observed rotation rate variations

Context. The SWAN Lyman α photometer onboard SOHO monitored the hydrogen cloud around comet 67P / Churyumov-Gerasimenko (67P) postperihelion at the three last perihelions in 1996, 2002, and 2009. Aims. Combining the SWAN results with some new Rosetta data allows estimating the erosion rate of the comet, quantiﬁed by the thickness of a layer that is disposed of at each orbit. Methods. By integrating the production rates measured with SWAN in time and adding some estimates for periods that are not covered by SWAN, we estimate the total H 2 O mass loss per orbit to be 2 . 7 ± 0 . 4 × 10 9 kg. We also tried to explain the observed change in the rotation rate during the 2009 orbit (period decrease of 1285 s) and the change observed by Rosetta from June 2014 to February 2015 (period increase of 32 s and 98 s up to 17 May 2015) with three di ﬀ erent mechanisms: sublimation-induced torque, thermal dilatation, and separation between the two lobes. Results. The total ejected mass depends on dust-to-gas mass ratios (4 ± 2) determined from Rosetta. This means that a layer of 1 . 0 ± 0 . 5 m thickness is lost at each orbit. The outgassing-induced torque may explain the observed changes in the rotation rate around perihelion in 2009 and recent changes. The torque decelerated the rotation from August 2014 to 17 May 2015, at which time it changed sign and began to accelerate the rotation, consistent with the average behavior observed for the 2009 apparition. Conclusions. The thickness of lost material needs to be kept in mind when interpreting all surface features. At 1 m ± 0.5 m, the erosion rate per orbit is high and supports the idea that the composition of the material that is measured in the coma (gas and solid) is indeed representative of the bulk material of the nucleus. We also argued that monitoring the rotation rate yields a very accurate and precious indicator of the global activity of the comet with which other activity measurements can be compared.


Introduction
During this first and historical rendezvous with the comet 67P/Churyumov-Gerasimenko (67P), the scientific payload of the Rosetta ESA spacecraft collected an unprecedented amount of new data.The Optical, Spectroscopic, and Infrared Remote Imaging System (OSIRIS) revealed breathtaking landscapes of surface features (Sierks et al. 2015;Thomas et al. 2015), the formation of which is hard to understand, to say the least (pits, terraces, basins, sand dunes, fractures, layers, cliffs, boulders, etc.).It is our opinion that any interpretation must take into account the erosion rate of the comet, quantified by the thickness of the material removed at each orbit.While more will be known on this subject at the end of the Rosetta mission (well after perihelion in mid-2015), it was felt important to publish estimates based on previous perihelion passages, as an early guide for the interpretation of the observed features.This is possible by combining data collected from Solar Wind Anisotropies (SWAN) onboard the Solar and Heliospheric Observatory (SOHO) at the last three perihelions (Bertaux et al. 2014) with some data already collected from Rosetta, such as the bulk density, the surface of the nucleus of 67P, the gas production rate, and a new mass dust-to-gas ratio, which is somewhat higher than previously derived.
Comet 103P/Hartley 2 was visited by the EPOXI mission in 2010, and Thomas et al. (2013a) used the Lα measurements of observations made with SWAN/SOHO (Combi et al. 2011a,b) as we do here, to constrain the overall loss of H 2 O over one orbit.Assuming a dust-to-gas ratio of 1 and 30% in mass of CO 2 , and a density of 520 kg/m3, they found a total massloss of 1.2 × 10 10 kg for the 1997 apparition and 4 × 10 9 kg for the 2010 apparition, which translates into a layer of 4.4 m lost in 1997, that is, 3% of the total mass of the nucleus, and 1% in 2010.For comet 9P/Tempel 1 (visited in 2005 by the Deep Impact mission and in 2011 by the Stardust spacecraft), the estimate of the layer thickness was much lower (≈1/3 m per orbit or ≈3 × 10 −4 of total mass per orbit, Thomas et al. 2013b).
Here we used actual SWAN/SOHO measurements of Q(H 2 O) as reported in Bertaux et al. (2014) from comet 67P combined with several parameters determined from the current encounter of Rosetta with the comet.
In the following section, we briefly describe the SWAN data (with a production rate of Q(H 2 O) as a function of time).Then we introduce the method we used to integrate the SWAN data The error bars do not include a 30% systematic uncertainty.The change in the level of the production rate near perihelion from 1996 to 2009 is remarkable.However, the activity levels are more similar by 25 days after perihelion.These data are taken from Bertaux et al. (2014).The dashed line is an interpolation between a 2002 point at 18 days after perihelion and a 2009 point at 45 days after perihelion to complement the 2002 SWAN measurements after the last point at 28 days.The solid line is an interpolation between a point at 1.4 AU suggested by SWAN values in 2009 and the Rosina/COPS estimate at 2.6 AU.
over time, as well as some extrapolations that together yielded an estimated total H 2 O loss rate for the three last orbits of the perihelion passages in 1996, 2002, and 2009.Finally, we take into account the contribution of some other gases (mainly CO and CO 2 ), the gas density measurements of Rosetta (Bieler et al. 2015), and the Rosetta-observed dust-to-gas mass ratio to estimate the thickness of material lost at each orbit.In Sect. 3 we discuss the kinematics of the nucleus rotation during two epochs to verify whether the observed outgassing rate is responsible for the observed changes in the rotation period.The first epoch refers to the orbit centered on the perihelion in 2009, where a decrease of the rotation period of 1285 s was measured (Mottola et al. 2014).The other epoch refers to the changes observed since the approach in 2014 until May 2015, which showed a reverse trend: an increase of the rotation period of up to 98 s.For the sake of completion, we also discuss the effect of a possible thermal dilatation on the rotation as well as of the progressive separation of the two lobes of the nucleus of comet 67P.

Mass-loss estimates
When comets are within 3-4 AU from the Sun, solar UV photolysis of H 2 O and then of the OH radical is the main source of atomic hydrogen in the coma, producing a huge envelope that is illuminated through resonance scattering by solar hydrogen Lyman-α radiation at 121.566 nm that becomes detectable by space-borne Lyman-α photometers, as was first reported for comet Bennett (Bertaux & Blamont 1970).When these observations are compared to appropriate models, it provides a reliable measure of the water production rate and its variation in time in comets, as was done for comet Bennett (Bertaux et al. 1973) and many more comets since.The primary scientific mission of SWAN/SOHO is to determine the latitude distribution of the solar wind and its changing pattern during the solar cycle (Bertaux et al. 1995).
An important secondary objective of SWAN is the monitoring of outgassing of comets.From its viewpoint in a halo orbit around the L1 Lagrange point between the Earth and Sun, at 1.5 × 10 6 km, SOHO is beyond the geocorona, and the cometary H Lyman α is not impeded by the exospheric Lyman α emission of the terrestrials H atoms, as is the case for observatories that orbit at low altitude, like the Hubble Space Telescope (HST).The usual observation mode of SWAN for monitoring the solar-wind distribution is the mode designed to collect full-sky Lyman-α maps.This means that active comets may be detected in these maps, and this was the case for the 67P apparitions of 1996 and 2002.For interesting comets, a dedicated observation campaign was organized in the SWAN program.This was the case for 46P/Wirtanen near its 1997 perihelion (when the comet was still the target of the Rosetta mission) and of 67P around the perihelion in 2009, when it was decided to change the target of Rosetta because of some concerns about the reliability of the Ariane launch rocket before the scheduled launch to 46P/Wirtanen.The SWAN data for 46P/Wirtanen were analyzed with a simple Haser model (Bertaux et al. 1999), while subsequent analyses (including the current data for 67P) to retrieve Q(H 2 O) from Lyman-α maps were made with a more sophisticated model developed by Michael Combi at the University of Michigan that is described in detail in Mäkinen & Combi (2005).
The measurements of Q(H 2 O) from SWAN are displayed in Fig. 1 as a function of time for the three apparitions in 1996, 2002, and 2009 with 4, 10, and 28 measurements, respectively.All observations were collected after perihelion except for two in 2002.The variable error bars are due to interferences from hot stars in the FOV of SWAN, which vary with the position in the sky.In addition, there is an overall uncertainty of about ±30% resulting from the typical uncertainties in instrument calibration, model, and model parameters.To extract an estimate for the total H 2 O sublimation along one full orbit by integrating Q(H 2 O) versus time, we assumed that the production is equal before and after perihelion.Then, the post-perihelion period was divided into five time sections.
1. From perihelion to the end of the SWAN measurements.2. From the end of the SWAN measurements to the time when a heliocentric distance of 1.4 AU is reached.The production rate decreases seriously at this point.3. From the time when 1.4 AU is reached to 2.6 AU. 4. From the time when 2.6 AU is reached to 3.6 AU, we rely on the in situ measurements of the total gas density and production rate of the Rosetta COPS sensor in 2014 (Bieler et al. 2015). 5. Beyond 3.6 AU, we have no estimate of the gas production, but there is an estimate of the dust production rate from the analysis of many comet or coma ground-based images taken in 2007 and 2008 (Snodgrass et al. 2013).
Ejected masses of H 2 O, dust, and gas are displayed in Table 1, where we indicate whether the quantity is estimated over half an orbit or a full orbit.The details of estimates are described below.
For period 1, we ignored the two points before perihelion for 2002 and numerically integrated during the period of actual measurements (line 1 of Table 1).
For period 2 and the perihelion in 1996, there are only four SWAN measurement points from 5.1 to 34 days after perihelion, at which point the distance to the Sun was 1.362 AU, and a distance of 1.4 AU was reached 45 days after perihelion.We linearly extrapolated the two last points of SWAN between 34 and 45 days after perihelion (line 2 of Table 1).
For period 2 and the perihelion in 2002, there are ten SWAN measurement points, eight of which were taken 2.3 to 28 days after perihelion, at which time the distance to the Sun was 1.33 AU, and a distance of 1.4 AU was reached 45 days after perihelion.To estimate Q(H 2 O) from 1.33 to 1.4 AU, we interpolated and extrapolated logarithmically between the point in 2002 at 1.31 AU and the Q(H 2 O) value for 2009 at 1.361 AU.The resulting curve is the dashed line in Fig. 1, obtained somewhat arbitrarily with mixing the passages in 2002 and 2009 (line 2 of Table 1).
For period 2 and the perihelion in 2009, there are 28 SWAN measurement points taken from 2.2 to 50.3 days after perihelion, at which time the distance to the Sun was 1.39 AU, and a distance of 1.4 AU was reached 53.5 days after perihelion.We ignored the last SWAN point, which is possibly an outlier on the low side, and added a production rate of 3× 10 27 mol/s for 3.5 days (line 2 of Table 1).
Line 3 of Table 1 is the H 2 O production from perihelion to 1.4 AU, the sum of lines 1 and 2 (half orbit).
For period 3, we interpolated (linearly in a log-log plot of Q(H 2 O) versus the heliocentric distance Rh) between a value of Q(H 2 O) = 2.1 × 10 27 mol/s at a distance of 1.4 AU (as suggested by 2009 SWAN measurements) and the actual COPS measurements at 2.6 AU of Q(H 2 O) = 1 × 10 26 /1.2, at the end of 2014, as derived in Bieler et al. (2015).The division by 1.2 was made to transform COPS measurement of the total density into a production rate of H 2 O, taking into account 10% of CO 2 and 10% of CO.This represents a rough estimate of the average ratios CO 2 /H 2 O and CO/H 2 O as measured by the Rosina mass spectrometer, although these ratios are found to be variable, depending on the location of Rosetta w.r.t. the nucleus (Hässig et al. 2015).The corresponding Q(H 2 O) is represented by a solid black line in Fig. 1.For the time integration away from SWAN measurements, we used the actual ephemeris of the comet for all three perihelions (Rh distance versus time) that are found on the Institut de Mécanique Céleste et de Calcul des Ephémérides (IMCCE) server for the ephemeris of solar system objects1 .
The total mass production of H 2 O from perihelion up to 2.6 AU (line 5 of Table 1) was computed by summing the productions for the three first periods.The total gas mass ejection (line 6) is obtained by adding 45% to account for 10% CO 2 and 10% CO outgassing (as said above).The associated error is estimated to be about 40%.
For period 4 from 2.6 to 3.6 AU, the COPS in situ total gas measurements were transformed by Bieler et al. (2015) into a production rate by scaling to sublimation models, taking into account the actual shape of the nucleus, the actual solar illumination, and the expansion of the coma (Bieler et al. 2015).Taken together, COPS data indicate that the production rate was constant from early August 2014 (at 3.6 AU) to 1 November 2014 at a value Qtotal = 5 × 10 25 mol/s, in agreement with the estimate of 4 × 10 25 H 2 O mol/s at 3.5 AU derived from the MIRO submillimeter instrument (Gulkis et al. 2015).After 1 November, COPS data indicate that the production rate doubled and remained constant at 1 × 10 26 mol/s up to 1 January, 2015, around 2.6 AU.Therefore, we integrated these two total gas production rates according to these two constant values, yielding a gas mass of 3.08 × 10 7 kg (line 7 of Table 1).
For period 5 beyond 3.6 AU, the analysis of Snodgrass et al. (2013) indicates a dust content in the coma varying as a power law of the heliocentric distance Rh.The retrieved quantity is Afρ, (where A is the albedo, f is the filling factor, and ρ is the radius of integration in the coma, in cm).Measurements were fitted with single power-law fits pre-and post-perihelion, given by Af ρ = 958 × Rh −3.18 (cm) and Af ρ = 1552 × Rh −3.35 (cm), respectively.
The quantity Af ρ is (with simple assumptions) proportional to the dust production rate, with a convenient relationship where Af ρ in cm ≈ Q d in kg s −1 (A' Hearn et al. 1995).We used this relationship and integrated in time along the orbit.Snodgrass et al. (2013) indicated that there was no detectable activity at 4.4 AU in September 2007.The integration from 3.6 to 4.3 AU yields 1.52 × 10 8 kg and 1.95 × 10 8 kg, respectively, before and after perihelion, with an average of 1.73 × 10 8 kg for a half-orbit (line 8 of Table 1).If the integration is pursued from 4.3 to aphelion at 5.68 AU, it would add another 3 × 10 8 kg.But since no activity was recorded beyond 4.3 AU, we preferred not to add this somewhat dubious extrapolation outside of the range where it was established.This dust production (3.6 to 4.4 AU) also corresponds to some gas production that must be added to the whole mass estimate.Taking a dust-to-gas ratio of 2 or 6 (extreme values considered here, see discussion below), this yields a small extra mass of gas of 1.73 × 10 8 kg (line 12) or 5.8 × 10 7 kg (line 13) for a full orbit beyond 3.6 AU.
During period 4 from 2.6 to 3.6 AU, the gas mass as retrieved from COPS is 3.0 × 10 7 kg (line 7 of Table 1).Therefore, the total gas (up to 3.6 AU, still on half-orbit) is obtained by summing lines 6 and 7, to achieve the result presented in line 9.
We now rely on a dust-to-gas ratio to estimate the total ejected mass.From a combined analysis of OSIRIS images of individual grains and of the GIADA dust collector (which collects different grain sizes), the dust-to-gas ratio was estimated to be 4 ± 2 (Rotundi et al. 2015).This wide range reflects the fact that the dust mass is dominated by the largest grains, which cutoff size is difficult to determine from photometric measurements of dust (e.g., Fulle et al. 2004).Since we need an average estimate over the full orbit and since we are not sure that the measured ratio at 3.5 AU is valid over the full orbit, we considered two extreme hypotheses to estimate the total mass loss: 2 and 6 for the dust-to-gas ratio, called hypothesis 1 and hypothesis 2, respectively.
Assuming that the loss is symmetrical w.r.t. the perihelion, we may now estimate the total losses of dust for the two dustto-gas ratio hypotheses, lines 10 and 11, respectively; adding the gas of line 9 and lines 12 or 13 yields the total mass ejected of gas+dust on lines 14 and 15, from 1 to 3 × 10 10 kg, depending on assumed dust-to-gas ratio.
This means that a significant fraction of the comet total mass (10 13 kg, Pätzold et al. 2015) is lost at each perihelion passage at 1.3 AU: between 0.1 and 0.3 percent.The volume of the lost material (lines 16 and 17 of Table 1) was computed from the measured average density of 470 ± 45 kg/m 3 derived from the mass and the total volume of 21.4 ± 2.0 km 3 estimated from OSIRIS images (Sierks et al. 2015).Following this, we derived the thickness of the lost material from the surface of the nucleus, which area was estimated to be 47.4 km 2 from the shape model of the nucleus (Preusker et al. 2015), as analyzed by Jorda et al. (2015).We found values of between 0.50 and 1.37 m (lines 18 and 19 of Table 1).The total H 2 O mass per orbit (line 20) was 2.3 to 3.1 × 10 9 kg, depending on the apparition.
Considering the results of Table 1, one may compute the relative productions rates of H 2 O within the four orbit sections from perihelion to 1.4 AU, 1.4-2.6AU, 2.6-3.6 AU, and 3.6-4.6AU.We found that H 2 O production is dominant within 1.4 AU, representing 73-76% of the total, the other sections representing 20%, 1.8%, and 1.7 to 5% of the total (according to the assumed dust-to-gas ratio in the outer part of the orbit), respectively.Our SWAN/SOHO measurements of 2009 cover most of the inner section from perihelion to 1.4 AU, and the absolute error is about 30%.The main factors of uncertainty in our estimate of the total mass loss come from the hypothesis of symmetry about perihelion (loss only depends on the distance), and the assumed average dust-to-gas ratio.
Taking the two extreme hypotheses of 4 ± 2 for the dust-togas ratio from Rotundi et al. (2015), this wide range probably also includes the other sources of uncertainties that affect the result, which are − the assumption that the water and dust production rates are symmetric versus the perihelion; − the assumption that the dust-to-gas ratio is constant over the whole orbit, while we relied on the Rosetta numbers measured at 3.5 AU pre-perihelion (Rotundi et al. 2015).
Neither of these two assumptions is consistent with the analysis of Fulle et al. (2004), who found (from the analysis of dust brightness mapping at 2002 perihelion passage of 67 P) that the dust production rate was approximately constant within 150 days from perihelion, but at 100 kg/s before perihelion, and 10 kg/s after perihelion, and found the opposite behavior of the gaseous coma observations of OH, CN, C 3 , C 2 , and NH during the passage in 1982, with more gas found after perihelion than before (A'Hearn et al. 1995;Cochran et al. 1992).This would imply a large decrease of the dust-to-gas ratio in the time before and after perihelion.By assuming a symmetry about perihelion, we might overestimate the gas production (since we mainly considered post-perihelion SWAN H 2 O data), but we might underestimate the dust production before perihelion with a constant dust-to-gas ratio.One effect might more or less compensate for the other in estimating the total ejected mass.Notwithstanding some contradictions, and hoping that all systematic effects will not conspire to go in the same direction away from our central estimate, we summarize our results by describing the thickness of the lost layer as 1.0 ± 0.5 m per orbit, which encompasses a range of a factor of 3. Clearly, the situation will be greatly clarified after the passage in 2015, when post-perihelion data will be analyzed, and better estimates of the thickness of the lost layer will be possible, hopefully directly by OSIRIS NAC (Narrow Angle Camera) images of terrain changes.This value of 1.0 ± 0.5 m per orbit is an average over the whole nucleus.It is certainly lower in some areas, and it is certainly more in other areas.Moreover, it is not clear at this point which fraction of the solid material is detached from the nucleus and escapes for ever, and which fraction returns to the ground after each perihelion passage.It is clear that the returned material constitutes the blanket layer of "smooth" material seen in many areas.There is some evidence, however, that this layer is not very thick in many areas.For instance, an image of the socalled pits (Fig. 8 in Sierks et al. 2015) shows the inside wall of the pit and its "goosebump features", while the external part of the pit is covered by a blanket.The thickness of the blanket seen from the side is much thinner than the characteristic scale, ∼3 m, of the goosebumps.
In most parts of the comet, the whole smooth blanket is therefore probably completely blown out at each perihelion, with a only small fraction falling back on the ground.This would prevent the formation of a solid crust, which has sometimes been invoked to describe the surface of the nucleus.One exception might be the smooth terrain in the region of the neck, denominated the Hapi region (Thomas et al. 2015), which is not illuminated closer to the Sun, as pointed out by the referee.Keller et al. (2015a) have also addressed the question of erosion of comet 67P in a totally different way.They used the actual shape of the nucleus described by more than 10 5 facets, and a "simple" ice sublimation model, considering three cases: dirty ice directly exposed to the Sun (model A), or covered by a porous dust layer of 50 µm (model B) or 1 mm (model C).They concluded that only a small fraction of the surface nucleus is active (but spread everywhere like on a chess board).There is also a large difference of erosion between the north part (erosion around 1 m or slightly below), while the southern hemisphere, which is fully exposed near perihelion, would experience more than 4 m, reaching 10 m at some places.They found water mass losses for models A, B, C of 6.5 × 10 10 kg, 4.9 × 10 10 kg, and 1.5 × 10 10 kg, respectively, and an ejected material layer thickness of 14.5, 11, and 3.4 m.This is more than our estimates, which is based on actual H 2 O outgassing measurements.The authors checked that their model predicts a higher H 2 O outgassing rate than observed by SWAN at perihelion (factor of 16 for their favorite model B), meaning probably that only 6% of the southern hemisphere surface is active, as described by their model.
We may compare our result for 67P to erosion studies performed on two other comets that have been encountered by space missions, for which the size of the nucleus and its area were also determined.For comet 9P/Tempel 1, a loss of ≈1/3 m per orbit or ≈3 × 10 −4 of the total mass per orbit was already mentioned (Thomas et al. 2013b).For comet 103P/Hartley 2, this was 3% of the total mass of the nucleus, and 1% in 2010 (Thomas et al. 2013a).Therefore, at a loss rate of 0.1-0.3%mass per orbit, comet 67P is in between these two cases on a log scale.We note that for the three comets that have been visited by space missions, the estimate of the lost layer thickness varies over a smaller range (from 1/3 m to 4.4 m per orbit, about one order of magnitude range), while the loss of relative mass per orbit ranges over two orders of magnitude.
On the other hand, Groussin & Lamy (2003) have estimated the erosion rate for comet 46P/Wirtanen to be 0.5 m per orbit.They used a radius of the comet of 0.6 km derived from HST observations, and their surface activity model was also adjusted to the observed perihelion production rates of H 2 O.This value is quite comparable to our value found for comet 67P, but it also implies a much greater relative mass loss because the comet is smaller.
Since the encounter with Jupiter in 1959, comet 67 P passed only eight times through its perihelion.Given our erosion rate of 1 m per orbit, this would correspond to a cumulated erosion of only 8 m thickness.
Orbital simulations of comet 67P backward in time (Groussin et al. 2007) showed that before 1959, the perihelion distance was larger than 2.7 AU for a duration of about 400 years, implying little erosion and little geological evolution of the nucleus.At this time (around year 1597), there was a close planetary encounter with a broad variety of possible scenarios before that encounter (chaotic behavior).Some of them include periods of time when the comet had a perihelion around 1 AU, while others have an orbit with large perihelions, and therefore little erosion.Therefore, if some features observed now imply an erosion much stronger than 8 ± 4 m, it would mean that the former scenarios (with small perihelions) are favored.We present an example of such an interpretation, admittedly highly speculative at this time.There are a number of circular features (Fig. 9 of Auger et al. 2015) that we could interpret as being the remnants of ancient pits, similar to the numerous circular pits detected (Sierks et al. 2015;Vincent et al. 2015) on the nucleus.The material infalling in the bottom of the pit would remain there and could form some consolidated material.When the erosion would peel off the nucleus around the pit down to the bottom of the pit (or even below), it would remain a circular feature similar to those observed.Vincent et al. (2015) showed that some pits have a depth of up to 20-200 m, which is much larger than the erosion since 1959.Therefore, if the observed circular features are shown to be remnant bottoms of ancient pits, this would be a strong indication of substantial erosion before 1959, favoring the scenarios before 1597 in which comet 67P experienced a period of its orbital lifetime in the inner solar system, with a perihelion distance to the Sun much smaller than 2.7 AU, rather than other scenarios where the comet perihelion was always larger than 2-3 AU.Another important result of the simulations of Groussin et al. (2007) is the dynamical lifetime of comet 67P, which is estimated to be about 105 000 yr.
We now investigate the role of the material ejection in the observed changes in the spin rate of the nucleus.

Observations
The rotation of the nucleus of 67P is mainly a pure rotation around the axis of maximum inertia (Sierks et al. 2015).In Table 2 we list the epochs, periods, rotation rates and changes in rotation rates that are discussed for comet 67P.Epoch 5 for 26 May, 2015 was added during the revision process of this paper with new available data, and the discussion was updated accordingly.Here we consider the measured changes in rotation rate both during approach and escort phase (August 2014 to May 2015, epochs 3 to 5 of Table 2, Fig. 2) and the change between 2004-2007 (before 2009 perihelion, epoch 1) and well after the 2009 perihelion, during the 2014 approach phase (Mottola et al. 2014).During this approach phase (in the period 23 March to 24 June, 2014, epoch 2), while the nucleus was still unresolved in the OSIRIS NAC, a time series of photometric measurements allowed accurately pinpointing a sidereal rotation period of 12.4043 h (Mottola et al. 2014) or 12 h 24 mn 15.48 s corresponding to a rotation rate ω = 0.000140704 rad/s (Table 2).This period was significantly shorter (by 1285 s) than the one determined from ESO NTT ground-based observations at large distance from the Sun in the period from February 2004 to September 2007 (epoch 1 of Table 2), well before the perihelion in 2009 (Lowry et al. 2012): P sid = 12.76137 ± 0.00006 h.This decrease in period most certainly occurred around the perihelion passage in 2009, most probably forced by sublimation-induced torques.
On the other hand, the combination of NAVCAM images, tracking of Rosetta, and some OSIRIS images has allowed the ESOC Rosetta navigation team to measure the rotation period of the nucleus of 67P and its change since the beginning A38, page 5 of 10  This recent behavior of decreasing rotation rate is the inverse of what occurred during the perihelion passage in 2009.Before discussing the implications, we recall a few kinematics equations governing the laws of the change in the angular momentum of a rotating body (boldface indicates vectors), as schematically illustrated in Fig. 3: in which I is the moment of inertia of the nucleus, ω is the rotation rate (in radians/s) = 2π/P (P period in s), r is one vector joining the center of mass of the nucleus to one point of its surface at distance r, d f is the elementary force exerted by outgassing from the surface, X is the vectorial product, and the integral is made over the whole surface of the nucleus.According to Mottola et al. (2014), there are very few observed changes in the rotation period of comets, and (to the best of our knowledge) they were all interpreted with the moment of inertia being kept fixed.In this case, Eq. ( 1) becomes Here we first keep the possibility to change the moment of inertia.This is motivated by the two-lobe structure of the nucleus of 67P, the presence of a crevice in the neck (Hapi region), and the possible increase of the distance between the two lobes as a sign of separation.There is also the case of thermal dilatation, which would increase the moment of inertia.For a spherical nucleus, there is a rotation rate ω L for which the centrifugal force is equal to the gravity at equator, which is independent of the radius, where ρ is the density and G is the gravitational constant.Taking ρ = 470 kg/m 3 gives ω L = 3.62 × 10 −4 rad/s, corresponding to a period slightly shorter than five hours.It is likely that a rotation faster than the limit ω L will result in a disruption of the comet nucleus.Indeed, in his review of known periods of comets (determined at that time by regular patterns in the coma assigned to the A38, page 6 of 10 rotation of the nucleus), Whipple (1982) did not find any comet with a faster rotation among 47 comets for which a rotation period could be measured.More recently, Weissman et al. (2004) derived a more sophisticated formula than our spherical case of Eq. ( 3) for the case of an oblate nucleus.The oblateness may be obtained from carefully analyzing the light curve of the bare nucleus.They found that the shortest period of the 14-comet sample available (oblateness and rotation) was 5.6 h, just above the limiting value ω L = 5 h for a spherical nucleus at 470 kg/m 3 .If the same decrease of rotation rate (as in 2009) is experienced at each perihelion of comet 67P, the limit ω L will be reached in 56 orbits, around year 2375.
Going back in time, the rotation rate would have been 0 only 35 orbits in the past, but nothing would preclude that the rotation would reverse its sense with the same exerted torques.However, the extrapolation of the orbit back in time cannot be made accurately before 1959, when a close encounter with Jupiter was experienced, and before this, the perihelion distance was much larger and the corresponding torque smaller (Groussin et al. 2007).

Case without torque
In a first approach, we ignore the effect of torques and examine the observed changes of the rotation rate between epochs 1 and 2, then between epochs 3 and 4 or epoch 5 (Table 2).Between two epochs, the change in rotation rate ∆ω and the change in momentum of inertia ∆ I are related by We consider two phenomena that could induce a change in the moment of inertia I.These are thermal dilatation of the nucleus associated with the thermal wave penetrating the nucleus while the comet is approaching its perihelion, and the other is an increase of the distance between the two lobes of 67P.Form epoch 1 to epoch 2 (before and after perihelion in 2009, respectively) the angular velocity increased by 3.94 × 10 −6 rad/s (the period increased by 21.4 mn).Since the dilatation is connected to the thermal wave, it should be roughly symmetric around perihelion, with a retraction after perihelion.Therefore we do not expect the velocity increase to be due to this phenomenon.And since an increase of the rotation rate ω is seen, it should induce (Eq.( 3)) a decrease of I, while the separation of lobes would induce an increase of the moment of inertia.Therefore, the observed increase of ω is most likely connected to the torque exerted by outgassing of the nucleus.The situation is different from epoch 3 to epochs 4 and 5 while Rosetta was in the vicinity of the nucleus.The angular velocity was seen to decrease and therefore could be the sign of either thermal dilatation or separation of the two lobes.
Separation of the lobes: considering the masses of the two lobes M 1 = αM and M 2 = (1−α)M, with M being the total mass and a = a 1 + a 2 the distance between the two lobes, the sum of the distances a 1 and a 2 of the center of gravity of each lobe to the center of gravity of the nucleus.The moment of inertia I can be described as Therefore, combining Eqs. ( 4) and ( 5), we have the relationship According to Jorda et al. (2015), a = 2600 m, yielding da = 0.92 m for epoch 4 and 2.8 m for epoch 5.This increase in separation could be confirmed (or refuted) by a careful analysis of landmarks on both lobes.Thermal dilatation: the evolution of the spin rate of comet 9P/Tempel over two perihelion passages was determined and discussed thoroughly by Belton et al. (2011).They found that the period decreased by 16.8 ± 0.3 min during the passage in 2000 and by 13.7 ± 0.2 min during the passage in 2005, suggesting a secular decrease in the net torque.These authors considered only the torque from outgassing.In addition to the main effect at perihelion, they also found a small deceleration before perihelion and a symmetrical acceleration after the perihelion.This somewhat mimics what can be expected from the effect of thermal dilatation or retraction (their Fig. 22), and it was tempting to also try this hypothesis with comet 67P.
For a homogeneous sphere, the moment of inertia is I = 2/5 M R 2 .An increase in length ∆ L due to thermal expansion with an increase in temperature ∆ T is ∆L = Lβ∆ T .The thermal dilatation coefficient β of a comet is totally unknown for the time being.As a guide, β = 12 × 10 −6 /K for steel, β = 50 × 10 −6 /K for water ice, and may reach 150 × 10 −6 /K for synthetic organic materials like Rilsan.If the whole comet were inflating homogeneously, we would have the relation Since we observe between epochs 3 and 4 dω ω = −7 × 10 −4 , with β = 100 × 10 −6 it would require a modest increase of the internal temperature ∆T = 2.1 K, or 18 K if β = 12 × 10 −6 /K.The radius R = 2000 m of the comet would have increased by 0.42 m in both cases between epochs 3 and 4 and by 1.3 m between epochs 3 and 5.But the actual situation is more complicated, since only an outside layer of the comet sees a change in temperature along its orbit.As discussed in McKay et al. (1986), one may define a depth below which the temperature is constant throughout the orbit.The thermal skin depth should be about 7 m for comet 67/P, and therefore only the dilatation of this external part could increase the moment of inertia.Indeed, if only the skin is affected by dilatation, it is found that the skin would have to move outward by about 30 m to change the inertia enough to induce the observed change in ω between epochs 3 and 4, which is absurd, and would certainly have been noticed in the navigation analysis of NAVCAM images.Distinguishing between the two options (separation between the lobes, and global thermal inflation) is probably possible by separately examining the inflation of both the body and the head, and of the total nucleus.

Torque effect of outgassing
Now we ignore the changes in the shape of the nucleus (dilatation or separation of the lobes) and only consider the effect of outgassing as exerting a torque on the nucleus, that can change the rotation state.The moment of inertia I is kept constant, and Eq. ( 2) applies.This equation can be integrated over time to relate a change in the rotation rate ∆ω during two epochs E i and E j to an integrated flux of ejected material during the same time.In this section we wish only to evaluate orders of magnitude, and some approximations are justified.The moment of inertia I is approximated as the one for a homogeneous sphere of radius R, I = 2/5 M R 2 .The force exerted by outgassing is f = Vdm/dt, where dm/dt is the mass ejection rate, and V is the velocity vector of ejection.The magnitude of the torque R X f is RV sin α dm/dt, where α is the angle between the radius vector R and the velocity V.It is generally assumed that the gas velocity is perpendicular to the local surface, but the expression is correct even if this assumption is not fulfilled.For simplicity we here only consider the component of the torque aligned with the angular momentum, which affects the rotation rate.Assuming a certain pattern of the outgassing (integrated over one rotation) to be stable over some time between epochs E i and E j , the geometrical part of the integral of Eq. ( 1) can be separated from the time integral of dm/dt, and Eq. ( 1) may be rewritten as where sin α is an average value of sin α over the nucleus.The angle α and its sine may be either positive or negative, increasing or decreasing the rotation rate, respectively.From Eq. ( 8) an estimate of sin α may be computed by placing numbers in the other terms that are known or estimated.For instance, if we were to find sin α larger than 1, this would demonstrate that the torque forces are insufficient to explain the observed change in rotation rate.Only the gas production rate must be taken into account in the integral of Q over time.This is because the solid material that is lifted off from the surface of the nucleus takes its energy and momentum from the gas drag.For an estimate of V, we rely on the detailed study of Davidsson et al. (2004) of the force exerted on a porous cometary material outgassing water vapor.In the frame of their layer energy absorption model (LEAM) that considers the ice sublimation inside a porous material, they computed that a force F of 0.18 newton per unit surface element would be exerted by a gas flow of Z = 2.5 × 10 −4 kg m −2 s −1 (their Figs.11 and 12 at maximum temperature).Therefore, their estimate of V is derived from F = Z V , which yields V = 720 m/s.Another approach (Zakharov, priv.comm.) would be to take the average of the thermal velocity components in the normal-direction averaged over only those molecules moving in the positive normal-direction: sqrt(2kT/πm H2O ) (which is one half of the averaged thermal speed), yielding 260 m/s for T = 230 K.In view of the actual Doppler velocity measurement of H 2 O line obtained by MIRO (Gulkis et al. 2015) of 680 m/s, we prefer to consider the 720 m/s value.
We may now derive sin α from ( 8), with M = 10 13 kg, R = 2000 m, observed ∆ω from Table 2, and an estimate of the integral of Q dt between two epochs.Between epochs 3 and 4 the length of period is 156 days.With an average production rate of Q(H 2 O) = 10 26 mol/s (or 3 kg/s) and V = 720 m/s, we find sin α = −1.7 × 10 −2 , a quite reasonable number.The corresponding angle is 0.95 • .The minus sign indicates that the torque is exerted against the existing angular momentum, yielding a decrease in the spin rate.If the nucleus were perfectly spherical and the outgassing perpendicular to the surface, α would be 0 and the torque would also be 0.But in view of the complicated shape of the nucleus (Sierks et al. 2015;Thomas et al. 2015), the direction perpendicular to the surface is rarely coincident with a direction toward the center of mass, therefore the angle α is often as large as several tens of degrees.
It is nonetheless quite striking that a modest production rate of 3 kg/s thrown almost vertically can still slow the rotation of a body of 10 billion tons in a measurable amount.Between epochs 1 and 2, which encompass the perihelion passage in 2009, our estimated gas production (line 9 × 2 for a full orbit+line 12 of Table 1) is approximately 3.5 × 10 9 kg; in this case, we find sin α = +1.2× 10 −2 , and a corresponding angle 0.7 • .This number is quite similar to the number derived from epochs 3 to 4, although it is of opposite sign.This change of sign is most likely linked to the geometry of the outgassing, dictated by a combination of the particular shape of the nucleus, the inclination of the rotation axis on the orbital plane, and the solar latitude around perihelion.While during the observations of epochs 3 to 5, the northern region was illuminated by the Sun and the southern region was not, the situation is different when integrated over a whole orbit.The Sun rapidly crossed the equator of the comet in May 2015 and will reach a maximum southern latitude of -52 • three weeks after perihelion.On 28 February 2015, we wrote: "Therefore, if the 2009 outgassing scenario repeats itself in 2015, we may expect that the rotation of the nucleus is going to continue to slow down for a while after February 2015, but then this slow down will stop when the Sun will explore more southern latitudes, the torque will reverse, and an acceleration will take place with a net result on the rotation rate change along one orbit similar to the one observed over the last orbit".With new measurements of the rotation rate provided by the ESOC Flight Dynamics team, Fig. 2 was updated with data up to 8 June 2015.The time derivative of the measured rotation rate is also very informative.It is expressed as being proportional to the outgassing rate Q(t) and sin α , as shown in Eq. ( 9) derived from (8): The rotation rate data of Fig. 2 was numerically time-derived to obtain the acceleration of Fig. 4.This acceleration shows a complex structure with strong periodic oscillations and six welldefined peaks: the first one in 2015 doy 41.6 (10 February), the last one in 2015 doy 121.8 (1 May).The period of the oscillation is determined with the time lapse of 80.27 days between the first and the last peak, yielding a period of 16.05 ± 0.15 days.We do not believe that these oscillations of dω/dt are caused by oscillations of Q(t) or oscillations of sin α ; for the time being, we interpret these oscillations as linked to the existence of a precession of the spin axis with this particular period, and we suspect that these oscillations are not real, but rather an artifact.The optical measurements of the NAVCAM camera could introduce a periodic artifact because the retrieval method of the rotation rate from images is by comparison with a model that probably does not include a precession.In fact, when the strong oscillations are smoothed out (by averaging in a box car of 32 days = two

Conclusions and perspectives
We have estimated the total mass of water outgassed at each of the three last perihelions of comet 67P.Taking into account some information collected during the early escort phase from Rosetta instruments (outgassing rate, contribution of CO and CO 2 , dustto-gas ratio, mean density, and total area of the nucleus), we were able to estimate the erosion rate per orbit in the present epoch, quantified by the thickness of the layer of lost material: about 1.0 m, within ±50%.This is an important parameter that must be taken into account for the interpretation of the numerous features that are observed at the surface of this fascinating comet.In particular, the erosion rate per orbit (0.1 to 0.3% of the total mass) is significant.The comet is peeling off at each orbit, supporting the idea that the composition of the material that is measured in the coma (gas and solid) is representative of the bulk material of the nucleus.
We have considered three potential mechanisms for a change in the rotation rate of the nucleus in two circumstances.The mechanisms are thermal dilatation, progressive separation of the two lobes, and the more classical outgassing-induced torque.The two circumstances are the changes associated with the perihelion orbit in 2009, and the approach in 2014 to early 2015, which is well documented by Rosetta ESOC measurements and analysis.For the perihelion orbit in 2009, gradual separation of the two lobes is excluded since it would have induced a slower rotation rate, while a higher rate was observed.Thermal dilatation is also excluded for this perihelion orbit, since this is a reversible phenomenon around perihelion, unless orbital changes would induce a complex history of the internal temperature profile, with a net cooling during the present orbit with perihelion at 1.3 AU.It would mean that in the recent past the orbit was nearer to the Sun than it is now.This contradicts the fact that before the encounter with Jupiter in 1959, the perihelion distance was larger at 2.7 AU.In contrast, we computed that an outgassing-induced torque is a plausible hypothesis in view of our previous SWAN/SOHO measurements of outgassing.The observed increase in the spin rate indicates that the net average angle α of outgassing is a modest 0.7 • .We calculated that the maximum spin rate ω L would be reached in 56 orbits around year 2375, meaning a rapid disruption of the nucleus at this time.
During the approach to perihelion in 2014 to early 2015, Rosetta accurately measured an increase in the period of 32 s up to 23 February, and 98 s up to 17 May, 2015.The implications were explored for the three mechanisms.A gradual separation of the two lobes by 0.92 m and 2.8 m for 23 February and 17 May, respectively, could explain the observed change in spin rate.A global inflation of 0.42 m and 1.3 m, respectively, of the nucleus could also explain the observation.However, if thermal dilatation were the cause, it would have to affect the whole nucleus, which contradicts what we understand of the thermal regime inside the nucleus.One could explore the stretching influence of the solar tide, but this is beyond the scope of this paper.A sublimation-induced torque is the most plausible explanation, requiring an average angle of α = 0.95 • , similar (but opposite in sign) to the one required for the integrated orbit in 2009.The change of sign can be explained by insolation geometry.
In addition, the refined analysis performed by ESOC Flight Dynamics team of Rosetta NAVCAM images (based on landmarks) allowed accurately measuring the decrease in the spin rate (Fig. 2).At end of February 2015, the spin rate decreased rather fast.We estimate that this fast decrease is difficult to explain with the mechanisms of thermal dilatation.If it is connected with the separation of the two lobes, it will soon be detectable in the refined analysis images, the distance between head and body landmarks should increase more than the distance between landmarks on the same lobe.In the case of the A38, page 9 of 10 most likely mechanism (outgassing), we wish to emphasize the importance of the spin rate change measurements.They provide a global picture of the activity of the comet, integrated over the whole nucleus and over a few days, which is not attainable by any other Rosetta measurements at present.This quantification of the activity is relative and not absolute for the time being.It is valid as long as the insolation geometry does not change significantly.To make it absolute, one would have to compare the various activity estimates from the various instruments to this spin rate deceleration.
We have indeed seen around 17 May, 2015 the reversal of torque that we predicted in our first version of this paper, which was written in February 2015.The acceleration of the rotation rate became positive and very likely will stay positive after early June, in conformity with the effect of the jet activity averaged over a full orbit, as was determined over the apparition of comet 67P in 2009 (Mottola et al. 2014).

Fig. 1 .
Fig. 1.Water production rates in comet 67P/Churyumov-Gerasimenko. They are derived from Lyman-α single-images collected by SWAN/SOHO at three last perihelion passages in 1996, 2002, and 2009.The error bars show the formal error resulting from instrument noise and model fitting.The error bars do not include a 30% systematic uncertainty.The change in the level of the production rate near perihelion from 1996 to 2009 is remarkable.However, the activity levels are more similar by 25 days after perihelion.These data are taken from Bertaux et al. (2014).The dashed line is an interpolation between a 2002 point at 18 days after perihelion and a 2009 point at 45 days after perihelion to complement the 2002 SWAN measurements after the last point at 28 days.The solid line is an interpolation between a point at 1.4 AU suggested by SWAN values in 2009 and the Rosina/COPS estimate at 2.6 AU.

Fig. 2 .
Fig. 2. Variations of the rotation rate measured since the rendezvous of Rosetta with comet 67P up to 8 June 2015.These data were obtained by the ESA navigation team at ESOC through a refined analysis of many landmarks on the nucleus observed by the ESA navigation camera NAVCAM.Bottom scale: day of year 2015.Top scale: days before 2015 perihelion.The perihelion is on August 13, 2015, or day 225.

Fig. 3 .
Fig. 3. Schematic representation of the torque exerted on the nucleus by the outgassing.The flow of gas V from one point makes an angle α with the line joining the center of mass of the nucleus and the point of outgassing.

Fig. 4 .
Fig. 4. Black triangles: deceleration of the spin rate since rendezvous with comet 67 P as a function of time from August 2014 to 6 May 2015.Bottom scale: day of year 2015.Top scale: days before 2015 perihelion.They were time-derived from ESOC data of Fig. 2, obtained by the ESOC navigation team from the fit of observed landmarks with a shape model and an assumed constant spin rate over a period of a few days.The strong oscillations with a period of 16 days are most likely not real and linked to the simplifying assumption of no precession (see text).The thick gray line is a smoothed average over 32 days.It shows that the spin rate change went through a minimum around April 10 and changed sign around mid-May 2015, a sign of outgassing pattern change linked to the solar illumination geometry change, passing to more southern latitudes.

Table 1 .
Gas-and dust-mass ejection per orbit.

Table 2 .
Rotation rates of the nucleus of 67P.
sid Rotation rate ω Change in rotation rate ∆ω