The extraordinary composition of the"blue comet"C/2016 R2 (PanSTARRS)

We present a multi-wavelength study of comet C/2016 R2 (PanSTARRS). This comet was observed on 23-24 January 2018 with the IRAM 30m telescope, and in January to March 2018 with the Nan\c{c}ay radio telescope. Visible spectroscopy was performed in December 2017 and February 2018 with small amateur telescopes. We report on measurements of CO, CH3OH, H2CO and HCN production rates, and on the determination of the N2/CO abundance ratio. Several other species, especially OH, were searched for but not detected. The inferred relative abundances, including upper limits for sulfur species, are compared to those measured in other comets at about the same heliocentric distance of about 2.8 AU. The coma composition of comet C/2016 R2 is very different from all other comets observed so far, being rich in N2 and CO and dust poor. This suggests that this comet might belong to a very rare group of comets formed beyond the N2 ice line. Alternatively, comet C/2016 R2 (PanSTARRS) could be the fragment of a large and differentiated transneptunian object, with properties characteristic of volatile-enriched layers.


Introduction
Comets are the most pristine remnants of the formation of the solar system 4.6 billion years ago. Investigating the composition of cometary ices provides clues to the physical conditions and chemical processes at play in the primitive solar nebula. Comets may also have played a role in the delivery of water and organic material to the early Earth (see Hartogh et al., 2011, and references therein). The latest simulations of the evolution of the early solar system (Brasser & Morbidelli , 2013;O'Brien et al., 2014) suggest a more complex scenario. On the one hand, ice-rich bodies formed beyond Jupiter may have been implanted in the outer asteroid belt, participating in the supply of water to the Earth, or, on the other hand, current comets coming from either the Oort Cloud or the scattered disk of the Kuiper belt may have formed in the same trans-Neptunian region, sampling the same diversity of formation conditions. Understanding the diversity in composition and isotopic ra-⋆ Based on observations carried out with the IRAM 30m telescope. IRAM is supported by INSU/CNRS (France), MPG (Germany) and IGN (Spain). ⋆⋆ The radio spectra are available at the CDS via anonymous ftp to cdsarc.u-strasbg.fr (130.79.128.5) or via http://cdsweb.u-strasbg.fr/cgi-bin/qcat?J/A+A/ tios of the comet material is therefore essential for the assessment of such scenarios (Altwegg & Bockelée-Morvan, 2003;Bockelée-Morvan et al., 2015).
Recent years have seen significant improvement in the sensitivity and spectral coverage of millimetre receivers enabling sensitive spectral surveys of cometary atmospheres and simultaneous observations of several molecules. We report here observations of a very peculiar comet, C/2016 R2 (PanSTARRS), with the 30-m telescope of the Institut de radioastronomie millimétrique (IRAM). This object is a long-period, dynamically old, Oort cloud comet that passed close to the Sun 21 600 years ago and has an orbit inclination of 58 • . It is very peculiar in the sense that in autumn 2017, while approaching the Sun (perihelion was on 9.6 May 2018 UT at 2.60 AU), it exhibited a deep blue coma and tail, due to the presence of strong CO + lines in the optical spectrum and very little dust or other emission lines. In contrast, usual comets show a dust tail and a coma of neutral or yellowish colour, resulting from the scattering of solar radiation by dust and the emission of C 2 Swan bands, if any at this heliocentric distance, which are not seen here. It seems to belong to a category of comets of which we know only very few examples: C/1908 R1 (Morehouse) (de la Baume Pluvinel & Baldet , 1911) or C/1961 R1 (Humason) (Greenstein , 1962), none having been observed with modern astronomical facilities. As in those two comets that were however observed closer to the Sun, optical spectroscopy of comet C/2016 R2 revealed unusually strong N + 2 lines (Cochran & McKay , 2018). Here we present millimetre observations of CO, CH 3 OH, H 2 CO, HCN, CO + and upper limits on several other molecules obtained with the IRAM 30-m, complemented by observations of the OH radical with the Nançay radio telescope, and amateur observations of CO + , N + 2 and dust continuum. From the millimetre data we derive constraints on the gas temperature and outgassing pattern. Production rates and their time evolution on 23-24 January, as well as abundances relative to CH 3 OH or CO are derived. We then compare the abundances in comet C/2016 R2 with those measured in eight other comets observed at a similar heliocentric distance.

Observations conducted at IRAM 30-m
Comet C/2016 R2 (PanSTARRS) was observed with the IRAM 30-m radio telescope on two consecutive evenings, 23.8 and 24.8 January 2018 UT, under good weather (mean precipitable water vapour of 1 mm). We targeted the CO(J=2-1) line at 230538.000 MHz first, before the HCN(J=3-2) line at 265886.434 MHz and the H 2 S(1 10 −1 01 ) line at 168762.762 MHz (see CDMS catalogue, Müller et al., 2005, the line frequencies uncertainties (¡ 1 kHz) being much lower than the finest sampling of the spectrometer used (20 kHz)), with other molecular lines, especially of methanol, in the band or the other side band. We used the EMIR 1 mm and 2 mm band receivers (Carter et al., 2012) with the fast Fourier transform spectrometer (FTS) set to a 200 kHz sampling covering ∼ 2 × 8 GHz in two polarizations simultaneously.
A log of the observations is given in Table 1. We used the secondary mirror wobbler, with a 180 ′′ throw at a frequency of 0.5 Hz to cancel the sky background. The comet was tracked with JPL#14 orbital solution. Coarse mapping of the strong CO line shows that the peak intensity was shifted by 2-3 ′′ to the west northwest (WNW) (Figs. 1 and 2). We adjusted the pointing offset to integrate closer to the peak of intensity. All pointing offsets mentioned henceforth in Tables 1 and 6 are relative to the position of this peak of intensity.
The comet immediately revealed itself as very different from most other comets observed so far. The CO(2-1) line is one order of magnitude stronger than expected (based on a correlation between visual magnitudes and CO production rates from Biver (2001)). Meanwhile the narrow and blueshifted line shape (Fig. 3) strongly resembles the one of distant comet 29P/Schwassmann-Wachmann 1 (Gunnarsson et al. , 2008), whose activity is dominated by a continuous, asymmetric and large outgassing of CO (≈ 4 × 10 28 molec. s −1 ). In total, 15 methanol lines are clearly detected; two lines of formaldehyde are also detected and the HCN(3-2) line is marginal. The sum of the three lines of HNCO covered by the observations that are expected to be the strongest is close to 3σ. Lines of neutral species (Figs.3-7) extend from -0.7 to +0.5km s −1 on the axis of Doppler velocities relative to the comet rest frame, which is due to the projection of the expansion velocity vector of the gas. In contrast, the two CO + lines at 236062.553 and 235789.641 MHz show a very marginal feature of ∼ 2 km s −1 in width superimposed on a more pronounced broad red-shifted line of 30 km s −1 . No other species are detected.
Sample spectra are shown in Figs. 1-7 and 8, and the line intensities and derived production rates are given in Table 5.  Colour-coded map of the CO(2-1) line integrated intensity of comet C/2016 R2 (PanSTARRS) obtained with the IRAM 30-m telescope on 24.85 January 2018. The (0,0) offset position corresponds to the JPL ephemeris position, with a ∼ 2 ′′ pointing uncertainty. The arrow points towards the projected direction of the Sun (phase angle of 17 • ).

Observation of OH radical with the Nançay radio telescope
In order to assess the water production rate of the comet, the OH lines at 1665 and 1667 MHz were searched for with the Nançay radio telescope between 2 January and 31 March 2018 (observation circumstances for selected periods are presented in Table 2). The telescope tracked the comet around transit for about 1 h on average every two days.

Complementary optical observations
In support of these observations conducted with large facilities, amateur astronomers followed the activity of the comet (which exhibited a total visual magnitude between 10 and 11) with their own telescopes. We report here on results obtained Notes. (a) Residual pointing offset (relative to the position of peak intensity measured on coarse maps) for the "On" nucleus pointings. (b) "On" nucleus pointings. But the first tuning of each day includes 15 and 22 min (23.8 and 24.8 January, respectively) spent at 5 to 14 ′′ pointing offsets when coarse mapping was carried out.  Fig. 3. Average on-nucleus (offset 1.7 ′′ ) spectrum of the CO(2-1) line at 230.538 GHz obtained with the IRAM 30m telescope on 24.80 January 2018. The vertical scale is the main beam brightness temperature and the horizontal scale is the Doppler velocity in the comet rest frame. The simulated profile with Q CO = 5 × 10 28 molec. s −1 in a 0-60 • cone at v exp = 0.56 km s −1 and Q CO = 5.5 × 10 28 molec. s −1 in a 60-120 • cone at v exp = 0.50 km s −1 is superimposed in green dotted line.
by Etienne Bertrand (St Sordelin), Hubert Boussier (IAU station K21) and François Kugel (IAU station A77). The first two astronomers obtained moderate-resolution (λ/δλ ∼ 500) visible spectra (λ=390-600 nm) of the comet, while the latter two also provided A f ρ values of the comet (measured outside the range of CO + lines). The data and logs of the observations are also accessible from the comet observation database hosted by the Paris Observatory . Hubert Boussier took spectra on 29.9 November, and 18.9 and 22.9 December 2017. In this paper we present an analysis of the spectrum from the latter of the three (Figs. 9,10) obtained between 22:01 and 22:42 UT with a 28-cm Schmidt-cassegrain telescope working at F/D=7.0 equipped with a LISA spec- Fig. 4. Average off nucleus spectra of the CO(2-1) line obtained with the IRAM 30m telescope on 23.78-24.85 January 2018: left: 12.1 ′′ tail-ward (∆RA> 0), right: 7.9 ′′ sun-ward. Scales are as in Fig. 3. The simulated profiles with Q CO = 5 × 10 28 molec. s −1 in a 0-60 • cone at v exp = 0.56 km s −1 and Q CO = 5.5×10 28 molec. s −1 in a 60-120 • cone at v exp = 0.50 km s −1 are superimposed in green dotted line.
trometer (Cochard and Thizy, , Shelyak instruments) and Atik 314L CCD detector. The slit width and integration window was 2.4 × 14 ′′ . Sky-glow and light pollution lines were removed using the signal in 2.4 × 33 ′′ windows at -38 and +72 ′′ along the slit. The effective spectrum resolution is 1 nm (5.1 Å per pixel). Etienne Bertrand obtained his spectrum (Fig. 10) on 22 February 2018 between 19:09 and 20:09 UT during a one-hour integration with a 20-cm Schmidt-cassegrain telescope working at F/D=6.3 equipped with an Alpy600 spectrometer (Cochard and Thizy, , Shelyak instruments) and Atik 414EX CCD detector. The slit width and integration window was 3.8 × 53 ′′ . Sky-glow and light pollution lines are removed using the signal in 3.8 × 37 ′′ windows at ±177 ′′ along the slit. The spectral resolution is also λ/δλ ∼ 500, that is, on the order of 1 nm. The average aper-  tures used for these observations (equivalent diameter of 8 ′′ and 14 ′′ for the 22 Dececember and 22 February observations, respectively) are similar to the IRAM 30-m beam, which is useful for comparing both types of observation.
The optical spectra are presented in Figs. 10 and 11. The correction for the system and sky transmission was done using a nearby A-type reference star and an out-of-atmosphere synthetic spectrum. Only a relative calibration was obtained (no absolute calibration of the fluxes in e.g. Wm −2 Å −1 ). Due to poor transmission of the atmosphere and optical system; and limited sensitivity of the detector towards the blue end of the spectrum, especially for the LISA spectrometer, the noise and uncertainty both increase below 400-430 nm. The wavelength calibration was performed using a reference internal source.  . Each line has been divided by its statistical weight (9/15 and 5/15 for the F = 5/2 − 3/2 and F = 3/2 − 1/2 transitions respectively) before averaging. Scales are as in Fig. 3. The line shows a narrow peak close to the zero velocity and a broader component mostly red-shifted indicative of antisunward acceleration of CO + ions by the solar wind.
The spectra are clearly dominated by the CO + doublets contrarily to most optical spectra of comets, dominated by CN and C 2 Swan lines (e.g. Fig.1 of Feldman et al., 2004). The N + 2 line at 391 nm is also detected in both spectra, but in the noisiest part of the spectrum; we assumed that there is no contribution from the nearby CN line at 388 nm, both because it was not detected by Cochran & McKay (2018) and because the abundance of HCN is very low (see following section).

IRAM data
Thanks to the simultaneous detection of two to six methanol lines and to the very high S/N obtained on the CO line, we were  The spectrum was extracted from the central 10 (binned) pixels rows. Green filled lines represent the simulated CO + spectrum with respective average g-factors. Blue filled lines show the N + 2 spectrum for a N + 2 /CO + column density ratio of 0.06. The vertical intensity scale has been normalized. The dotted line is the sum of N + 2 and CO + synthetic spectra. The signal towards the red end (¿500nm) of the spectrum is likely under-corrected because the LISA spectrometer is not optimised to cover the full wavelength range and the comet was slightly out of focus in that part of the spectrum. As a consequence, the signal of the lines beyond ∼510 nm is lower than expected and some lines are not detected. This part of the spectrum is not used in the analysis. able to constrain the excitation conditions of the molecules in the coma and the outgassing pattern of the comet. This is important to derive accurate outgassing rates and relative abundances. Table 3 provides the derived rotational temperatures of methanol lines and implied gas temperature. A rotational diagram of the 165 GHz lines is presented in Fig. 12. All groups of lines (at 165 GHz,(213)(214)(215)(216)(217)(218)(219)(220)242 GHz,(252)(253)(254)(255)(256)(257)(258)(259)(260)(261)(262)(263)(264)(265)(266) provide rotational temperatures T rot in very good agreement in the 18-24 K range. Due to the radiative decay of the rotational levels within the ground vibrational state, we expect colder rotational temperatures for the 242 GHz lines. A higher collisional rate would limit this decay and observations suggest that indeed the collision rate might be higher than modelled, possibly due to the outgassing in narrow jets resulting in higher local densities. Nevertheless we adopted a gas temperature T gas = 23 K to derive all the molecular production rates. We note also that if we use the "jet" part (velocity interval -0.8 to -0.2 km s −1 ) of the 165 GHz lines observed at high spectral (40 kHz) res- The spectrum was extracted from the central 50 pixels rows. Green filled lines show the simulated CO + spectrum with respective average g-factors. Blue filled lines indicate the N + 2 spectrum for a N + 2 /CO + column density ratio of 0.08. The vertical intensity scale has been normalized. Dotted lines is the sum of N + 2 and CO + synthetic spectra. We can also constrain the kinetic temperature T gas from the CO(2-1) line profile. Assuming that the width of the narrow component of the CO line is due to thermal broadening, its full width at half maximum (FWHM) of 0.237 ± 0.003 km s −1 implies T gas ≤ 34 ± 1 K. This is compatible with T gas = 23 K derived from the rotational temperatures. Modelling of the CO line shape (following section) yields information on the actual dispersion of velocities.

Outgassing pattern and expansion velocity
After determining the gas temperature and implied thermal broadening of the lines, we modelled the CO line profile to infer the gas expansion velocity and outgassing pattern. For T gas = 20 − 23 K the CO line shape asymmetry is well explained assuming an expansion velocity v exp = 0.56 ± 0.01 km s −1 on the Earth side and considering that most of the outgassing is restricted to a cone with a half-opening angle of 60 • . Figure 13 shows the CO production rate per unit solid angle (assuming symmetry around the Comet-Earth axis) as a function of the Earth-Comet-Vgas or colatitude angle (φ) inferred from the line profile. For this profile inversion (i.e. conversion of each channel intensity into a production rate per solid angle, with the channel velocity v i converted into a colatitude angle φ i from the relationship v i = −v exp cos(φ i )), we assumed that v exp (0.56 km s −1 ) does not depend on φ. In reality, we expect the expansion velocity to be lower at higher phase angle due to lower surface temperature. This is consistent with the smaller blue-shift (by 0.06 km s −1 ) observed in spectra acquired on the anti-sunward side (Fig. 4).
In any case, line profile inversion with a lower velocity down to 0.3 km s −1 still produces a cut-off of production around the colatitude angle φ= 120-130 • . To approximately simulate the profile of gas production that was to be used to derive total production rates, we assumed constant production rates per solid angle q 1 and q 2 in the colatitude ranges φ 1 = 0 − 60 • and φ 2 = 60 − 120 • , respectively. The expansion velocity was set to 0.56 km s −1 (φ 1 = 0 − 60 • ) and 0.50 km s −1 (φ 2 = 60 − 120 • ), respectively. We adjusted q 1 and q 2 to account for the observed mean Doppler shift of the line. q 2 is about half q 1 and the total production rates in these two regions (Q i = 2π φ q i dφ) are in the ratio Q 1 :Q 2 =10:11. The profile of the total production rate per solid angle q is plotted in Figs. 13 and 14. The resulting line profiles are superimposed to the observed line profiles in Figs. 3 and 4.
Although the agreement between simulation and observation is relatively good for the central position (Fig. 3), small differences appear at offset positions (Fig. 4). This shows that the actual outgassing pattern deviates somewhat from the axisymmetric description done in this study. Indeed the line intensity spatial distribution (Figs. 1 and 2) points at some asymmetry in the outgassing.
A similar profile inversion was undertaken for the average of the four strongest methanol lines (J=1 to 4) observed at 165 GHz. The result (Fig. 14) is similar to the one from CO, although with a lower S/N, suggesting that CH 3 OH and CO have similar production patterns from the nucleus. However the peak of the methanol lines appears at slightly more negative velocities, possibly linked to the uncertainty on the absolute line frequencies. Indeed the shift is ∼ 50kHz, which is on the order of the differences between available frequencies in 2016 and previously published values (Müller et al., 2005, and references therein). If we expect the same line shape for CH 3 OH and CO, the inferred CH 3 OH line frequencies are 165050.258 ± 0.016, 165061.153±0.023, 165099.220±0.017 and 165190.544±0.010 MHz, for the 1 1 − 1 0 E, 2 1 − 2 0 E, 3 1 − 3 0 E, and 4 1 − 4 0 E lines, respectively.

Molecular production rates
In summary, to analyse all spectra, we used a gas temperature of 23 K and the outgassing pattern described above. A simpler description of the outgassing pattern assuming hemispheric outgassing at v exp = 0.56 km s −1 yields production rates that are only 5% higher. Production rates are provided in Tables 5 and  6 and plotted in Fig. 16. They are calculated as described in Biver et al. (2016) and previous papers. Results from the maps are provided in Table 6 which gives average intensities at several radial offsets from the peak of intensity determined from the maps (Figs. 1 and 2) for the strongest lines and the derived production rates. There is only a marginal trend of increasing production rate with offset for CO (∼ +20%), not significant for the other molecules, which could also be related to instrumental effects (error beam and distortion of the beam shape at higher elevations). For the three molecules, observations are compatible with most of the production coming from the nucleus. Fig. 13. Inversion of the CO(2-1) line profile from Fig. 3, assuming a constant expansion velocity of 0.56 km s −1 . The production rate per steradian is derived as a function of the colatitude angle (φ=0 • on the Earth side, 180 • on the opposite), assuming symmetry along the Comet-Earth line. The dots correspond to the values derived for each spectral channel and the connected line with error bars is the running average taking into account thermal broadening. The grey shaded region corresponds to the production rate profile used to determine the total production rates, with v exp = 0.56 km s −1 and 0.50 km s −1 in the dark and light grey regions, respectively (see text).

Search for OH
In order to constrain the water production rate of the comet, the 18-cm OH lines were observed with the Nançay radio telescope. From 2 January to 26 March, 2018, 48 individual daily spectra of about one-hour integration were secured 1 . The observing and reduction procedures were as explained in Crovisier et al. (2002). During the observing period, the OH maser inversion varied from −0.31 to −0.15. Neither the individual spectra nor their average show any significant signal. For the 2 January -9 February average (corresponding to a mean maser inversion of −0.31), the 3-sigma upper limit on the 1667 MHz line area, integrated over a 3 km s −1 width, is 12 mJy km s −1 . According to the model of Crovisier et al. (2002), this corresponds to an upper limit of 1.1 × 10 28 molec. s −1 for the OH production rate, that is, 1.2×10 28 molec. s −1 for the water production rate. Observations and resulting upper limits are summarised in Table 2.

N +
2 /CO + from visible spectra Following the work presented by Cochran & McKay (2018), we estimate here the N + 2 /CO + column density ratio from the visible spectra presented in Sect 2.3. We used line frequencies from Lofthus and Krupenie (1977), Kepa et al. (2004), and Magnani and A'Hearn (1986). Fluorescence efficiencies (L/N or g-factors) of N + 2 are from Lutz et al. (1993) and those for CO + are from Magnani and A'Hearn (1986). We use those pertaining to the heliocentric velocity range -10 to +25 km s −1 (see sect. 6.3).
Synthetic fluorescence spectra were computed assuming a width of ∼ 0.4 nm for each v ′ − v ′′ band and an instrumental resolution of 1 nm. A N + 2 /CO + abundance ratio of 6-8% was used to obtain the simulation plotted in Figs. 10 and 11. Band areas were computed to estimate the relative intensities (Table 4). Due to non-perfect relative calibration across the visible band and also variable uncertainties due to removal of atmospheric lines (mostly light pollution lines of Hg-I and high pressure sodium (Na-I) broad emission above 540 nm, especially for the spectrum of 22.9 December, 2017, in Fig. 9), we took the weighted average of the column densities derived from each detected line. The band observed around 427 nm is the sum of the N + 2 (0-1) band at 427.4 nm and the CO + (2−0)Π 3/2 band at 427.2 nm (band head, Magnani and A'Hearn , 1986). In order to evaluate the possible contribution of the N + 2 (0-1) band, we computed the signal of the CO + (2 − 0)Π 3/2 using the average of the column density derived from the nearby CO + (2 − 0)Π 1/2 band at 424.9 nm and the column density derived from all CO + bands. The 424.9 nm band should be less affected by relative calibration issues, but the average of all lines should average out local fluctuations of the baseline and uncertainties on the g-factors.
The derived N + 2 /CO + column density ratio is 0.09 ± 0.03 for 22.9 December 2017 and 0.06 ± 0.02 for 22.8 February 2018, giving an average of 7 ± 2%. If we only take into account lines below 460 nm to limit biases with wavelength, the values be-  (Lutz et al. , 1993), and for CO + (Magnani and A'Hearn , 1986) for (-10, +30 km s −1 ) heliocentric velocity interval. (b) in relative units, normalized to the first CO + band: integrated band intensity (in ADU × Å) on the spectrum divided by L/N. (c) B : CO + (2 − 0)Π 1/2 blue wing and R : CO + (2 − 0)Π 3/2 red wing. (d) The intensity of the CO + (2 − 0)Π 3/2 band is derived from this range of column density estimated from the nearby Π 1/2 band and the average of all bands. (e) The intensity of the N + 2 (1 − 2) band is inferred from the subtraction of the expected contribution of the CO + (2− 0)Π 3/2 . ( f ) Average column density ratio, including a 20% additional uncertainty to each column density due to the use of approximate values of the L/N parameter which varies with the heliocentric velocity of the ions (Magnani and A'Hearn , 1986). come 0.08 ± 0.02 and 0.07 ± 0.03, respectively, yielding an average of 8 ± 2%. Since ionization efficiencies are similar for the two molecules (β i,N 2 = 3.52 × 10 −7 s −1 and β i,CO = 3.80 × 10 −7 s −1 at 1 AU for quiet Sun (Huebner et al., 1992)), this translates to a N 2 /CO production rate ratio of 8 ± 2%, consistent with the value of 6% determined by Cochran & McKay (2018).

Dust production
Our analysis on the dust production of the comet is solely based on optical data, as the detection of continuum emission in the infrared or in the radio has not been reported so far. The optical spectra present only a weak continuum (Fig. 9). Images show a "CO + " blue coma, and do not reveal any dust tail (Cochran & McKay , 2018;Noël et al. , 2018), suggesting a low dust production. From coma photometry, F. Kugel and H. Boussier ) report A f ρ values from 530 to 830 cm, with a mean value of 670 cm for the period from 8 February to 8 March, 2018 (at a nearly constant phase angle of 21 • ). The A f ρ parameter introduced by A' Hearn et al. (1984) is proportional to the dust loss rate multiplied by the cross-section of the dust particles within the aperture; it can be used as a proxy for the dust production. According to A'Hearn et al. (1995), a value of 1000 cm corresponds to a dust production rate of about 1 ton per second. For comet C/2016 R2 this would suggest that the dust production is about one order of magnitude less than the total gaseous production, in mass (A f ρ/Q gas ≈ 0.1 cm.kg −1 vs. 2-5 for other comets at 3 AU (Table 7) and A f ρ/Q H 2 O = 0.5 cm.kg −1 on average from A' Hearn et al. (1995) for comets closer to the sun). The conversion of the A f ρ parameter into dust production rate is highly dependent on assumptions made on the dust size distribution and maximum size. If released dust con-sists predominantly of large grains, then the dust production rate could be much higher than the estimation given above. These large grains could form a dust trail which is not seen in the images, including those acquired in mid-December 2017 when the Earth crossed the plane of the orbit of the comet. As discussed in the following section, the A f ρ value for C/2016 R2 is comparatively much lower than values measured in comets showing comparable gaseous activity.

Discussion
The derived CO production rate (10.6 × 10 28 molec. s −1 , equivalent to ∼ 5 tons per second) corresponds to a large outgassing rate, only surpassed by comet Hale-Bopp and by 17P/Holmes during its massive outburst (Table 7). This was unexpected given the small brightness of the comet (total visual magnitude m 1 around 10-11). Indeed, from the correlation between visual magnitudes and CO production (Biver , 2001) we obtain a CO production rate 20 times lower. This indicates a dust-poor coma, consistent with the very low A f ρ values (670 ± 110 cm) measured by F. Kugel and H. Boussier during this period. At the same heliocentric distance for comet Hale-Bopp the A f ρ was 100 times higher  (at a similar phase angle) for a CO production only twice higher. We note that the visual magnitude is difficult to compare to that of other comets as the optical brightness is dominated by CO + "tail" emission lines.

Molecular abundances
Pending the assessment of the abundance of other possible major molecules (CO 2 , O 2 ), for which no ionic emission has been reported so far, the major species in the coma of comet C/2016 R2 (PanSTARRS) at 2.8 AU from the Sun are CO and N 2 with an abundance ratio on the order of 100 : 8. Observations and searches for NH 3 or NH 2 and hydrocarbons are needed to complete the inventory, but these species are not expected to be as abundant as the previous species. Water is not detected with an upper limit H 2 O/CO < 0.1. Comparing with comets observed at similar heliocentric distances (Table 7), methanol is relatively abundant relative to water, but deficient compared to CO. Hydrogen cyanide and sulfur species are strongly depleted. Figure 15 compares the abundances of CO and HCN relative to methanol in comets observed at 2.3-3.3 AU from the Sun. Molecular production rates measured in these comets are listed in Table 7.

Upper limits on minor species
The wide frequency coverage (Table 1) has enabled us to look for many other species previously detected in comets (Biver et al. , 2015). The most significant upper limits are provided in Table 5. For species where several lines of similar expected intensities are sampled, we provide the combined upper limit (weighted rms) derived from all lines. In general the upper limits in comparison to the reference molecule (HCN for CN-species, CH 3 OH for CHO-species and H 2 S for S-species) are not very constraining (Biver et al. , 2015). The marginal detection of HNCO in C/2016 R2 requires confirmation as it corresponds to a HNCO/CH 3 OH ratio that is between three and ten times higher than in other comets and a HNCO/HCN two orders of magnitude higher. H 2 S is clearly under abundant (< 6% relative to CH 3 OH) in comparison to other comets observed at ∼3 AU from the Sun (Table 7: H 2 S/CH 3 OH = 30-70%). We also looked for the CO isotopologues 13 CO and C 17 O, but they were not detected. The best S/N is obtained considering the blue-shifted part of the lines, and gives 12 C/ 13 C> 54 (82 ± 42 at 2 σ), which is compatible with the terrestrial value (90).

Ions observed with IRAM-30m
The rotational lines of CO + were observed for the first time in comet Hale-Bopp at the CSO (Lis et al. , 1997) and the IRAM 30-m (Henry et al., 2001). The CO + lines at 236 GHz are marginally detected, red-shifted and broader than the CO line, as would be expected for ions (cf Hale-Bopp spectra of CO + and HCO + (Lis et al. , 1997)). The two strongest CO + (N = 2 − 1) transitions at 236062.553 and 235789.641 MHz show a marginal 3 − 2 σ peak of FWHM=1 km s −1 close to the zero velocity in the comet frame and a broader component of 25 km s −1 in width between -5 and +30 km s −1 more clearly detected (5 and 4 σ, respectively). This asymmetry towards larger red-shifted velocities is due to the acceleration of ions by the solar wind in the anti-solar direction. Figure 8 shows the combined spectrum (average of the two lines). Column densities cannot be derived in a straightforward manner as the properties of this profile show evidence of the acceleration of CO + ions in the beam due to interaction with the solar wind. A rough estimate for an expansion velocity of 10 km s −1 (mean Doppler shift of the lines), and a rotational population at 23 to 200 K is 0.4-1.7×10 12 cm −2 , but we have no precise idea of the rotational population of CO + ions which are in a collisionless regime. We assumed that CO + comes from the photo-ionization of CO (scale-length on the order of 12 × 10 6 km) which yields optically thin lines in any case. The marginal narrow component would yield a column density about five times lower. However, this profile can be used to obtain more appropriate g-factors to analyse the optical spectra taken in a similar aperture: given the heliocentric velocity of the comet and the small phase angle of the observations, this means that the heliocentric velocity of the CO + ions contributing the most of the signal is in the -10 to +25 km s −1 range. Hence we used the corresponding L/N values from Magnani and A'Hearn (1986) to interpret the optical spectra.
HCO + is not detected. It has been detected in comets with weaker productions of CO, but at closer distances to the Sun where water production and protonation of CO in the coma are much more efficient (Milam et al. , 2004).

Variation with time
Since the comet exhibited a unexpectedly large outgassing of CO, with some optical images showing rapid changes in the CO + structures, we looked for possible evidence of short-term variations due to either the rotation of the nucleus or a transient outburst phase. Figure 16 shows the production rates from Table 5 and their evolution over the two days of observations. On this short time scale, generally the variations were less than 20% and possibly more related to pointing and calibration uncertainties. Indeed, at high elevations (above 60-70 • ) the beam efficiency of IRAM-30m degrades. We made some modelled corrections of the beam efficiency but could not precisely track its variation. This may account for ∼10% variations. No significant variations of the CO(2-1) line Doppler shift (−0.28 ± 0.02 km s −1 ) that could be correlated with a variation in the production pattern -and independently of any calibration issue -are observed either. So, on a timescale of several hours, we do not see any significant variation of the activity.

Summary and conclusion
We performed a multi-wavelength (optical, millimetre, centimetre) compositional study of comet C/2016 R2 (PanSTARRS). This comet has the following characteristics.
-A very large CO production, on the order of 10 29 molec. s −1 (i.e. 5 tons/s) at r h =2.8 AU from the Sun, only a factor of two below that of C/1995 O1 (Hale-Bopp) at same r h . -Very low dust production, based on the Afρ which is lower than values measured in comets of similar gaseous activity at ∼ 3 AU from Sun by more than a factor of 15. -Unusual optical spectrum showing essentially CO + and N + 2 lines.
-Abundant N 2 in the coma (N 2 /CO∼0.08), with N 2 dominating the nitrogen budget. -A strong depletion of H 2 O, CH 3 OH, H 2 CO, HCN, H 2 S relative to CO (by more than one order of magnitude) compared with other comets observed at a similar heliocentric distance. -A depletion of HCN and sulfur species relative to methanol.
The origin of the huge production of CO of C/2016 R2 could be related to the size of its nucleus. The nucleus of comet Hale-Bopp has a radius estimated to R n ∼ 37 km (Altenhoff et al., 1999), making this comet one of the largest ever observed. Distant comet 29P/Schwassmann-Wachmann 1, which also shows a large CO production on its circular 6-AU orbit (∼ 4 × 10 28 molec. s −1 ), has a radius of ∼23 km (Stansberry et al., 2004). Unfortunately, measurements of the nucleus size of C/2016 R2 have not yet been reported. Given the pure CO ice sublimation rate at 2.8 AU Z CO (2.8) ∼ 1.6 × 10 22 molec.m −2 .s −1 (Sekanina , 1991), the outgassing rate of C/2016 R2 could be supplied by a pure CO ice object of 3 km in diameter. However, detailed thermal and structural modelling is necessary to provide valuable constraints on the nucleus size from the CO production rate.
The depletion of H 2 O, CH 3 OH, H 2 CO, HCN, H 2 S relative to CO, together with the low dust production, suggest that the large CO activity reflects a CO-rich ice composition and large ice/dust ratio within the nucleus material of C/2016 R2. Indeed, due to low latent heat, sublimation of CO ice can release smaller particles and larger aggregates than the outgassing of H 2 O ice (Gundlach et al., 2015). Therefore, if C/2016 R2 had a dust/ice ratio comparable to that of other comets, significant dust production should have been observed. Here, we assume that the thermal properties of the nucleus surface and subsurface layers of C/2016 R2, controlled by the heat conductivity and porosity of the material, are overall comparable to other comets. This seems a reasonable assumption. Indeed, the illumination-driven CO outgassing indicated by the CO line profile (Sect. 3.1.2) shows that the CO sublimation front is close to the surface, consistent with a low-porosity material. Other comets observed at r h > 3 AU from the Sun (Hale-Bopp and 29P) display a similar CO outgassing pattern (e.g. Gunnarsson et al., 2003). In comet Hale-Bopp, which had a higher A f ρ/gas ratio and higher abundance of H 2 O and HCN at 2.8 AU (Table 7), icy grains were also found in the coma (Lellouch et al. , 1998). We can speculate that following the low dust abundance, icy grains are also deficient in the coma of this comet, which could explain the low abundance of water and HCN, which could be significantly produced from the sublimation of icy grains in other comets like Hale-Bopp at such heliocentric distances.
The other particularity of comet C/2016 R2 is the large abundance of N 2 , dominating other N-bearing species.
Measured N 2 /HCN is < 0.006 for Hale-Bopp (derived from Cochran & McKay , 2018, using HCN/CO=0.01), and ∼0.2 in 67P at r h =3 AU (from Rubin et al. , 2015;Le Roy et al. , 2015). For C/2016 R2, N 2 /HCN = 2000, that is, four orders of magnitude higher. Since there is no hint of NH 2 lines in optical spectra, and other N-bearing species (Table 5) searched for at IRAM are less abundant than N 2 by at least two orders of magnitude, we can conclude that contrary to other comets, most of the nitrogen escaping the nucleus is in the form of N 2 .
The strong nitrogen deficiency in some comets, both in the ice and refractory phases, was revealed during the space investigation of 1P/Halley, and confirmed with the Rosetta mission (Rubin et al. , 2015;Fray et al., 2017). Its interpretation remains elusive. According to Owen and Bar-Nun (1995), the trapping of N 2 by amorphous water ice in the cooling solar nebula was inefficient within Neptune's orbit, resulting in the formation of planetesimals deficient in N 2 . Laboratory experiments show that the trapped N 2 /CO ratio is depleted by a factor of ∼ 20 at 24 K, with respect to the gas phase value, the depletion factor being strongly dependent on temperature (Bar-Nun et al., 2007). With a protosolar ratio N/C = 0.29 and assuming that all C and N is in the form of CO and N 2 , the N 2 /CO ratio in the solar nebula gas phase is 0.15, resulting in a trapped N 2 /CO of a few times 10 −3 at ∼ 25 K, relatively consistent with the value measured for 67P (Rubin et al. , 2015). If this interpretation is correct, the very low N 2 /CO values measured in, for example, comet Hale-Bopp would indicate a formation in warmer regions of the solar nebula.
Considering the gas-phase species detected in comet C/2016 R2, the N/C ratio for this comet is ∼ 0.15, that is, close to the solar value of 0.29±0.10 (Lodders et al., 2009). To explain this property, a possibility is that C/2016 R2 agglomerated from grains formed at an extremely low temperature, favouring the trapping of high quantities of N 2 both as trapped gas and in condensed form. However, the high CO content of C/2016 R2 ices together with the low dust/ice ratio (assuming our extrapolation from coma to nucleus composition properties is correct) may suggest another scenario. Models examining the thermal evolution of the relatively large Kuiper Belt object show that due to radiogenic heating, the most internal layers reach high temperatures (Prialnik et al., 2008;Sarid & Prialnik, 2009); released gases migrate towards colder regions where they refreeze. The more volatile ices refreeze closer to the cold surface than the less volatile, so that the pristine dust/ice mixture becomes enriched in volatile ices such as CO and N 2 (Prialnik et al., 2008;Sarid & Prialnik, 2009). Dynamical studies of the transneptunian population argue for a rich collisional history in the Kuiper Belt (Morbidelli & Rickman, 2015, and references therein). In the second scenario, comet C/2016 R2 would be a fragment of the disruptive collision of a large Kuiper Belt object, with properties representative of volatile-enriched layers.
Comet C/2016 R2 is representative of a family of comets that we observe only rarely each century. Besides C/1908 R1 (Morehouse) and C/1961 R1 (Humason) (de la Baume Pluvinel & Baldet , 1911;Greenstein , 1962), other candidates are comets 29P/Schwassmann-Wachmann 1 and C/2002 VQ 94 (LINEAR) (Korsun et al., 2008) which both showed optical spectra dominated by strong emissions of CO + and N + 2 , characteristics of abundant CO and N 2 production and high N 2 /CO ratio. The diversity of the dust/gas ratios seen in these comets (as judged from the Afρ values, which indicate that 29P and C/2002 VQ 94 are dust-rich in comparison to the other comets in the sample) may favour the second scenario in which these comets are collisional fragments of differentiated transneptunian objects.