Free Access
Issue
A&A
Volume 510, February 2010
Article Number A24
Number of page(s) 18
Section Planets and planetary systems
DOI https://doi.org/10.1051/0004-6361/200912804
Published online 02 February 2010
A&A 510, A24 (2010)

First attempt at interpreting millimetric observations of CO in comet C/1995 O1 (Hale-Bopp) using 3D+t hydrodynamical coma simulations[*]

J. Boissier1 - D. Bockelée-Morvan2 - A. V. Rodionov3 - J.-F. Crifo4

1 - Institut de Radioastronomie Millimétrique, 300 rue de la piscine, Domaine universitaire, 38406 Saint-Martin-d'Hères, France
2 - LESIA, Observatoire de Paris, 5 place Janssen, 92195 Meudon Cedex, France
3 - TsNIIMASCh, Central Research Institute of Machine Building, Pionyerskaya St. 4, Korolev, Moscow region 141070, Russia
4 - CNRS and UVSQ, LATMOS, BP 3, 91371 Verrières-le-Buisson Cedex, France

Received 1 July 2009 / Accepted 3 October 2009

Abstract
Context. Millimetre line observations of comet C/1995 O1 (Hale-Bopp) close to perihelion, completed using the IRAM Plateau de Bure Interferometer, have detected temporal variations in the CO J(2-1) 230 GHz line shape, and in the position of its maximum emission brightness within the field-of-view, whose heuristic analysis has suggested the presence in the coma of a slowly rotating spiral-shaped enhancement of the CO density.
Aims. Here, we reanalyse these data using a physically consistent model of the coma.
Methods. We consider a large, rotating, icy nucleus with an arbitrarily aspherical shape and an adhoc rotation mode, and compute its tridimensional, time-dependent (``3D+t'') mixed CO + H2O coma, using a previously developed tridimensional hydrodynamical code (HDC). The line emission of CO is then computed using a molecular excitation and radiation transfer code (ERC). In the present, pioneering phase, the HDC and ERC both contain crude, and not thoroughly mutually consistent approximations. Several alternative CO surface flux distributions are considered, and the resulting CO 230 GHz line spectra and brightness maps are compared with observations.
Results. We find that when an uniform surface flux of CO is assumed, the spiral structures created by the nucleus asphericity in the CO coma are too faint to account for the observational data, whereas we confirm earlier conclusions based on a heuristic approach that the assumption of an area of suitable dimensions and localization with increased CO flux leads to results in agreement with a large subset of (but not all) the data. This suggests that the true CO coma production map may be more complex than the presently assumed rather simple-minded one. Refined and mutually consistent HDC and ERC are needed for a more satisfactory interpretation of the present and any similar future data.

Key words: comets: individual: C/1995 O1 (Hale-Bopp) - radio lines: planetary systems - techniques: interferometric

1 Introduction

Spectroscopic observations of comets at millimetre wavelengths allow the study of a number of molecules in the inner coma. Single-dish data provide limited spatial resolution, allowing the determination of only large-scale averages of the coma temperature and velocity fields, from which little information about the coma structure, and all the more about the nucleus outgassing pattern, can be obtained. With interferometric observations such as those discussed here, we are capable of imaging the emission of inner coma molecules with far higher spatial resolution than before, making possible to compare the data with ab-initio physical models of the coma. Such observations, however, require bright enough sources, and thus are not frequent.

The extremely bright comet C/1995 O1 (Hale-Bopp) offered in 1997 the opportunity to undertake observations with several interferometers including the Berkeley IllinoisMaryland Association (BIMA) array (Snyder et al. 2001; Veal et al. 2000; Wright et al. 1998), the Owens Valley Radio Observatory (OVRO) (Blake et al. 1999), and the Plateau de Bure Interferometer (PdBI) operated by the Institut de Radioastronomie millimétrique (IRAM) (Bockelée-Morvan et al. 2009; Altenhoff et al. 1999; Henry et al. 2002; Boissier et al. 2007; Wink et al. 1999). The outburst of the comet 17P/Holmes transformed this faint comet into a reasonable target for interferometric observations with the PdBI (Boissier et al. 2008) and SMA (Qi, personal communication). The Hale-Bopp observations remain the most suitable set of available interferometric data for the study of the coma structure experiencing normal cometary activity.

One of the most interesting characteristics of Hale-Bopp observations is the detection of rotating spiral structures in the visible light images of the coma obtained at the wavelengths of the OH, NH, CN, C2, and C3emission (Lederer et al. 2009; Lederer & Campins 2002; Lederer et al. 1999). Analysing the PdBI observations of the CO 230 GHz line, Henry et al. (2002) and Bockelée-Morvan et al. (2009) detected spiral structure of CO gas. A detailed analysis of the observations of sulfur compounds (CS, H2S, and SO) in this comet obtained at PdBI suggests that density structures were present for SO and CS (Boissier et al. 2007), and a similar conclusion was drawn from the images of the emission of HNC, DCN, and HDO acquired at OVRO. Structures are rarely resolved in a gas coma: one can mainly cite observations of this kind in comets 1P/Halley, C/1996 B2 (Hyakutake), and 8P/Tuttle (e.g., Woodney et al. 2008; Crifo et al. 2005). The completion of their physical (i.e., gasdynamical) modelling is an important challenge. To date, only comets Halley (Szegö et al. 2002) and Hyakutake (Rodionov et al. 1998) have been studied in this way.

Section 2 summarizes the main characteristics of the CO J(2-1) (230 GHz) line observations acquired at PdBI in comet Hale-Bopp in March 1997, and presented in detail in Bockelée-Morvan et al. (2009) and Henry (2003). Because comet Hale-Bopp displayed a high CO production rate close to perihelion ($\sim$ $2 \times 10^{30}$ s-1, Bockelée-Morvan et al. 2009; Biver et al. 1999a), the observations of the J(2-1) 230 GHz line are of excellent signal-to-noise ratio. It was thus possible to study the time evolution of the recorded signals.

Bockelée-Morvan et al. (2009) analyzed these data with the help of a phenomenological model assuming a spherical rotating nucleus with uniform CO flux except within a region with a Gaussian flux enhancement producing a spiralling jet. Their best-fit solution was obtained with a jet half aperture of $\sim$18$^{\circ}$ located close to the nucleus equator (latitude $\sim$20$^{\circ}$) and accounting for $\sim$35% of the total CO production.

However, it is extremely unlikely that the true nucleus is spherical. In addition, the inferred very high (a factor $\simeq$14) local CO surface-flux increase seems surprising, as well as the simultaneously rather strong assumption that the CO flux is uniform over the rest (98%) of the surface.

In a phenomenological model, the molecules are assumed to be distributed in the coma in a mathematically simple manner, regardless of whether this distribution is compatible or not with the laws of fluid mechanics. For instance, for the sake of simplifying the radiative transfer computations, the above-mentioned model assumes a constant gas coma temperature and velocity, despite the presence of a density increase inside the spiral structure. It is easy to show that these assumptions imply non-conservation of the flow momentum. This raises the question of which physical process(es) could produce such an effect hitherto never advocated in the outer coma (see e.g., Rodionov et al. 2002).

In contrast, a physical model computes all gas flow parameters from explicitly defined physical processes (or, at least, aims to achieve this). The very high activity of the comet makes it possible to use a computationally very efficient gasdynamical approach - the 3D+t Euler equations. However, these equations must include terms allowing for all non-adiabatic physical processes taking place in the coma (such as dust loading, photoionization, photodissociation, molecular radiative emission, and absorption), and the computation of these terms requires 3D radiation transfer codes: (1) to evaluate the net budget of the molecular absorption and emission, the radiative transfer of the far-IR H2O spectrum must be computed; (2) since the dust IR emission enters into play there, the radiative transfer of the solar spectrum at visible wavelengths in the dust coma must be computed; (3) to evaluate the rate of dissociation of H2O, the radiative transfer of the far-UV solar spectrum must be computed; (4) to evaluate the heating rate due to the dissociation products H and OH, the transfer of the superthermal H atoms must be computed. This takes enormous effort, far beyond the present state of the art. In the present work, the required 3D transfer codes are replaced by approximate treatments. Our present modelling goal thus need not be a perfect fit to the data, but results exhibiting the main observed features, with correct orders of magnitude. However, one should keep in mind that, in this highly active comet, very little direct information about the nucleus has been retrieved from the observations at any wavelength, leaving most of the properties of Hale-Bopp's nucleus entirely open to speculation.

Rodionov & Crifo (2006) developed a time-dependent 3-D fluid model of the coma around a rotating aspherical nucleus producing CO and H2O, using parameters applicable to comet Hale-Bopp close to 1 AU. They showed that, even when the nucleus is assumed to be homogeneous in surface icy content, and have a constant surface CO flux, rotationally induced coma gas structures appear, due to the nucleus asphericity, demonstrating that the existence of gas spirals does not necessarily require surface heterogeneity. In the present paper, by computing the resulting CO coma emission, we investigate whether this effect is able to explain the structures identified in the CO interferometric data obtained at the Plateau de Bure in March 1997. The code of Rodionov & Crifo (2006) is slightly modified to allow a localized increase in the CO production inspired from the above-mentioned phenomenological model, and we similarly investigate whether this can explain the above-mentioned structure.

Section 3 contains a short description of the 3D+t HDC of Rodionov & Crifo (2006), and of the various simulations carried out with it here. In Sect. 4, we present the ERC used to simulate the observations on the basis of the results of the fluid model results. The results of the comparison between the observed and simulated data are presented in Sect. 5.

2 Observations

\begin{figure}
\par\includegraphics[angle=-90,width=8.5cm,clip]{12804f1.ps}
\end{figure} Figure 1:

ON-OFF CO J(2-1) (230.538 GHz) observations in comet Hale-Bopp on March, 11, 1997: a) seven representative individual ON-OFF spectra (1 min of integration ON source) recorded from 4.42 to 13.76 UT, (out of the 13 observed, Bockelée-Morvan et al. 2009; Henry et al. 2002). The intensity unit is the main beam brightness temperature. The UT time is given in the upper right corner; b) average of the 13 spectra recorded from 4.42 to 13.76 UT; c) line area as a function of UT time. The dashed line indicates the mean value; d) line velocity shift as a function of time. The dashed line indicates the value measured on the average spectrum.

Open with DEXTER

The CO J(2-1) (230.538 GHz) line observations discussed here were carried out with the IRAM Plateau de Bure interferometer on March 11, 1997 (Bockelée-Morvan et al. 2009; Henry 2003). The heliocentric distance of the comet was 0.99 AU, and its geocentric distance 1.37 AU. The observations included both single-dish observations (4.4-13.8 h UT) and interferometric observations (4.5 to 12.5 h UT). At that time, the Plateau de Bure comprised 5 antennas and was in a compact configuration, with baselines (the spacing between two antennas) ranging from $\sim$20 to $\sim$150 m. Data were acquired with spectral channels of 0.78 kHz separation, corresponding to $\sim$0.10 km s-1 in velocity scale.

The detailed description of the data and of their phenomenological analysis is given in Henry (2003) and Bockelée-Morvan et al. (2009). A short report was presented in Henry et al. (2002). We discuss here only the features that have been interpreted in terms of surface distribution of the CO flux at the nucleus surface - i.e., uniform flux production over the surface, apart from a small area where the flux is much stronger.

\begin{figure}
\par\includegraphics[width=8.5cm,clip]{12804f2.ps}
\end{figure} Figure 2:

Time evolution of the velocity shift measured on ON-OFF and visibility spectra. The antenna pairs (i.e., baselines) are indicated and sorted by increasing lengths (e.g., 3-4: 32 m; 2-4: 88 m, 3-5: 101 m; 2-5 : 147 m). Observations are shown by filled circles with error bars, and the fitted sinusoids are shown by dashed lines. Open squares and solid curves are the result of the phenomenological model fit. Figure adapted from Bockelée-Morvan et al. (2009).

Open with DEXTER

2.1 Single-dish data

In the so-called ON-OFF mode, the antennas are used as independent single-dish telescopes, and position-switching is used to remove the sky background. The full width at half maximum (FWHM) of the primary beam of the 15-m Plateau de Bure antennas is 20.9 $^{\prime \prime }$ at 230 GHz, which corresponds to $\sim$20 000 km (diameter) on the comet.

ON-OFF 1-min CO spectra were taken at different times during the observing period. Within each of them, the CO J(2-1) line is strong enough to be detected with a high ($\sim$50) signal-to-noise ratio. It is thus possible to study the time evolution of the line. Its profile changes drastically (Fig. 1a), the dominant emission swinging from positive to negative relative velocities. This can be conveniently characterized by the evolution of the mean velocity shift of the line, $\Delta v$(i.e., first moment with respect to the nucleus rest velocity), shown in Fig. 1d. Figure 2 (top left panel) shows that it can be reproduced well (over nearly one period) by a sine-like curve whose period is that of the nucleus rotation (11.35 h, Jorda et al. 1999). In contrast, the average profile (over the measuring period) is almost symmetrical (Fig. 1b), and the line area does not show significant variations (Fig. 1c). Given the orientation of the comet spin axis, which is close to the plane of the sky (the angle between the spin axis and the line of sight, i.e. aspect angle $\theta_{\omega}$ being $68^{\circ}$ based on the spin orientation of Jorda et al. 1999), these characteristics are expected if there is enhanced and sustained (night and day) CO production from a source close to the equator.

2.2 Interferometric data

In the interferometric mode, the signals received by the different antennas are correlated with each another. In this mode, the direct output of the observations are complex visibilities for each couple of antennas. We define for each combination of antennas (i.e., baseline) the baseline vector ${\vec \sigma}$ with components (u, v) in the Fourier plane. For each spectral channel i, of bandwith $\delta v$ (km s-1), the complex visibility $\mathcal{V}_i({\vec \sigma})$ (W m-2 Hz-1 or jansky) is defined by (see e.g., Thompson et al. 1991, Chap. 4, Sect. 4.1)

\begin{displaymath}
\mathcal{V}_i ({\vec \sigma}) = \frac{c}{\nu \delta v} \int_...
...\pi\nu}{c} {\vec \sigma} \cdot {\vec s}\right) {\rm d} \Omega,
\end{displaymath} (1)

where ${\vec s}$ is a vector in the sky plane with components (x,y)in radian units, $F_i({\vec s})$ (W m-2 sr-1) is the brightness distribution on the sky at frequency $\nu$ contributing to channel i, $A_{\rm pp}$ is the normalized (dimensionless) power pattern of the antennas (defined by the ${\it FWHM}$ of the primary beam), and $d\Omega$ is an element of solid angle on the sky.

Interferometric maps are obtained by computing the inverse Fourier transform of the available visibilities, and accounting for only a partial sampling of the Fourier plane being obtained (see e.g., Thompson et al. 1991). When all the CO 230 GHz Hale-Bopp Plateau de Bure data are considered, the ${\it FWHM}$ of the synthesized beam is 2.00 $^{\prime \prime }$ $\times$ 1.38  $^{\prime \prime }$ ( $1990 \times 1370$ km).

Visibility spectra $\mathcal{V}_i$ were recorded for each pair of antennas (i.e., for 10 baselines). As for the ON-OFF spectra, the velocity shift of the CO 230 GHz line in these visibility spectra was measured. The time evolution of the velocity shift measured for the 10 baselines is shown in Fig. 2. The velocity shift is most accurately measured for the shortest baselines (3-4 (32 m), 1-5 (48 m), 1-4 (55 m), 1-3 (60.8 m), and 2-4 (88 m), where the numbers label the different antennas) which sample typically roughly rectangular coma regions of widths 3500 to 10 000 km, and lengths 20 000 km (Bockelée-Morvan et al. 2009). The time interval over which data exist is only 2/3 of a rotation period. For convenience, the data were fitted with sine waves of the rotation period, with a degree of confidence that decreases with the signal-to-noise ratio. One sees that both the amplitude and the relative phase shift of the sines vary from spectrum to spectrum. As explained in Bockelée-Morvan et al. (2009), small differences can be expected because the different baselines probe different regions of the coma, of different sizes. Indeed, the elapsed time between the detection of the rotating structure by one baseline and its detection by another baseline causes some delay. The best-fit phenomenological model of Bockelée-Morvan et al. (2009) can reproduce fairly well two of the curves (baselines 3-4 and 2-4), but is in clear disagreement with the other ones, either in the position of the extrema (4-5), or in the magnitude of the variation (1-3, 3-5) or both (Fig. 2). The discrepancies are observed for baselines of field of view perpendicular to the spin vector (Bockelée-Morvan et al. 2009).

2.3 Emission brightness maps

Considering five successive subsets of the interferometric data, one can produce five interferometric maps of the CO 230 GHz emission brightness that have a signal-to-noise ratio allowing their individual study. We note that this brightness is not proportional to the CO column density because the CO line is optically thick. For comparison with the computed brightness maps, one can study the time evolution of astrometric position of the photometric center (i.e., position of the brightness maximum) with respect to the mean position measured on the overall map (bottom right frame of Fig. 3). The astrometric precision of the peak position, estimated by dividing the interferometric beam size by the signal-to-noise ratio, is 0.07 $\hbox {$^{\prime \prime }$ }$. The successive observed positions differ significantly. Their displacement is not random but follows approximately the direction perpendicular to the nucleus spin axis (the position angle of the projection of the spin axis onto the plane of the sky is $pa_{\omega} = 211^{\circ}$, based on the spin orientation of Jorda et al. 1999). This behavior is expected for a CO jet close to the nucleus equator in an optically thin coma: when the jet changes its direction, the whole brightness distribution is modified and the brightness peak moves close to the equatorial plane. As the present computations - which incorporate optical thickness effects - demonstrate, this remains true in the present optically thick case.

\begin{figure}
\par\includegraphics[angle=-90,width=9cm,clip]{12804f3.ps}
\end{figure} Figure 3:

Successive CO 230 GHz line integrated maps obtained from the data subsets covering: 1/ 4 h30-5 h30 UT, 2/ 6 h00-7 h00 UT, 3/ 7 h00-8h00 UT, 4/ 8 h50-9 h50 UT, and 5/ 10 h45-11 h45 UT. Twenty-nine velocity channels were averaged. The first positive contour and the contour spacing is 0.2 mJy (2-$\sigma $). The peak intensity at the center of the maps in units of line area is 6.76, 5.96, 5.40, 5.46, and 3.81 Jy km s-1 for maps 1 to 5, respectively. The beam shape and size is shown in the bottom left corner. The (RA, Dec) = (0 $^{\prime \prime }$, 0 $^{\prime \prime }$) coordinates correspond to the mean photometric center determined from the whole data set and coincides approximately to the nucleus position (Bockelée-Morvan et al. 2009). The projected Sun direction is indicated by an arrow. The offset between the brightness peak position measured on these maps and the position measured on the map that includes the whole data set is shown in the bottom right frame. The positions were obtained by fitting the Fourier transform of a Gaussian to the visibilities (i.e., the fit is performed in the Fourier uv plane). The positional accuracy is about 0.1 $\hbox {$^{\prime \prime }$ }$. The arrow on the lower-right panel indicates the direction of the nucleus north pole, projected on the sky. In the same panel, the grey line represents $\Delta O$, the distance between the most distant individual points. Figure adapted from Henry et al. (2002) and Bockelée-Morvan et al. (2009).

Open with DEXTER

2.4 Observational summary

To make the comparison between observations and model computations easier, we introduce the following seven quantities:

  • $v_{\min}$ and $v_{\max}$ (in km s-1) define the limiting projected velocities, in the nucleus rest frame, where CO emission is present in the averaged ON-OFF spectrum shown in Fig. 1b.

  • The velocity shift of an individual spectrum is $\Delta v =
\sum_{i} T_{i} v_{i}$/ $\sum_{i} T_{i}$, where Ti and vi are, respectively, the intensity and the velocity of channel i.

  • The amplitude of variation in the velocity shift during the observing period is $A = \Delta v_{\max}-\Delta v_{\min}$ (in km s-1), where $\Delta v_{\max}$ and $\Delta v_{\min}$ are, respectively, the highest and the lowest $\Delta v$ measured in the individual ON-OFF spectra.

  • $\langle \Delta v \rangle$ (in km s-1) is the mean velocity shift of the line measured in the averaged ON-OFF spectrum.

  • $\Delta O$ (in arcsec) is the extent of the displacement of the photometric center of CO emission during the observing period (Fig. 3).

  • $F_{\rm SD}$ and $F_{\rm Int}$ (in Jy km s-1) are, respectively, the line integrated ON-OFF and interferometric fluxes, the latter corresponding to the peak brightness in the interferometric maps. $F_{\rm SD}$ is measured in the averaged ON-OFF spectrum. $F_{\rm Int}$ is measured in the interferometric map obtained using the whole data set.
Table 1 provides the values of these quantities deduced from the observational data, together with the values deduced from the five alternative models which we describe below.

Table 1:   Main characteristics of ON-OFF spectra and interferometric maps for observed and simulated data (Sect. 2.4).

3 Gas production and outflow model

Rodionov & Crifo (2006) developed a 3D+t HDC to compute the CO-H2O coma created by a non-spherical, rotating nucleus. The code is a time-dependent, Eulerian version[*] of the 3D code described extensively in Rodionov et al. (2002). This code was developed to treat near-nucleus regions; for dealing with Hale-Bopp observations, its computational domain was extended up to 300 000 km (to allow for correct line-of-sight integrations), which implied taking into account new processes specific of the outer coma. The computed distribution of H2O and CO number densities, as well as of the mixture velocity and temperature, are stored at 72 equispaced times during the rotation. A thorough description of the code and the solutions discussed hereafter will be given elsewhere (Rodionov et al. in preparation). Here, we provide only an overview of the code and draw attention to the free parameters and approximations that it uses. There are two groups of parameters and approximations: one concerning the nucleus and its gas production mechanism, the other concerning the coma.

3.1 Nucleus and gas production

\begin{figure}
\par\includegraphics[width=9cm,clip]{12804f4.ps}
\vspace*{3mm}
\end{figure} Figure 4:

Nucleus used in the 3-D+t hydrodynamical simulations (Rodionov & Crifo 2006). The shape is that of comet Halley, smoothed and scaled-up to match the estimated size of comet Hale-Bopp. The illumination of the nucleus (grey scale) by direct Solar radiation is shown at two different rotation phases (separated by one third of the period). The rotation axis is vertical, the Sun is to the right.

Open with DEXTER

The arbitrarily assumed nucleus shape, shown in Fig. 4, is that of comet Halley, smoothed to eliminate small details, and scaled-up to the Hale-Bopp mean radius of 25 km, in agreement with the compromise ($30 \pm 10$ km) suggested by Fernández (2002) on the basis of several independent studies. The assumed rotation mode is also that of P/Halley, but with a period decreased to 11.4 h (Jorda et al. 1999).

Table 2:   Model parameters for the different solutions.

Following all authors (e.g., Prialnik et al. 2004), we assume CO to be produced by diffusion across the nucleus outer layers. In the present study, we further assume that the CO flux is distributed over the surface according to a pattern independent of the phase of the rotation, so that the total CO production rate is constant during the rotation. This is consistent with observational constraints: the line area of the CO ON-OFF spectra obtained with the PdBI does not show any significant temporal variation that can be attributed to a rotational modulation (Fig. 1c). The model postulates the surface flux of CO to be either constant over the entire nucleus (at a value $F_0 = Q_{\rm CO}/{\mathcal A}$ where ${\mathcal A}$ is the external area of the nucleus), or to be constant apart from a ``circular spot'' within which it is higher, producing a coma structure referred to as a ``jet'' in the following. The characteristics of this spot were chosen to match approximately those determined by Bockelée-Morvan et al. (2009): its axis is radial (i.e., passes across the nucleus center of mass C) and its CO flux excess at a point M is $\Delta
F_0 {\rm e}^{-[(\delta/20)^2 ln2]}$ where $\delta$ (in degrees) is the angle between ${\vec {CM}}$and the spot axis, and $\Delta F_0$ is such that the flux inside the spot, $\Delta Q_{\rm CO}$, represents about 30% of the total CO production $Q_{\rm CO}$. The physical size of the surface from which the excess flux (higher than half the jet maximum flux) emerges is 7 km (in radius). This corresponds to 1.7% of the nucleus surface. It is here important to notice that, in an aspherical nucleus, the longitude $\phi$ of the spot is also an important parameter, since the outflow direction is not the cone axis, as in the phenomenological model: the gas mixture leaves the surface in the direction of its ``average normal'', and gradually becomes radial, but in a final direction that depends on the latitude, on $\phi$, and on the local time (owing to the strong local time dependence of H2O). Here, $\phi$ was chosen arbitrarily.

In all cases, the adopted value of $\Delta F_0$ is such that the surface-integrated total CO production rate $Q_{\rm CO}$ (given in Table 2) is in order of magnitude agreement with the value $2.1 \times 10^{30}$ molecules s-1 derived by Bockelée-Morvan et al. (2009) assuming an isothermal coma at the temperature of 120 K, and a constant expansion velocity of 1.05 km s-1derived from the width of the radio lines (Biver et al. 1999a).

To compute the water production, the nucleus surface is assumed to be an intimate mixture of inert dust and ice, the latter producing H2O under solar illumination. The surface icy fraction f is assumed to be constant over the nucleus - including or excluding inside the CO spot (see below). The ice temperature $T_{\rm I}$ and net H2O flux $Z = f Z_{\rm I}(T_I,M_0)_{\rm I}$are computed at each point from the correct ice sublimation energy-budget equation

\begin{displaymath}(1-A_V) \; {\max}[\; \cos{z_\odot}, F_{s+c}\;] c_\odot
= \epsilon \sigma T_I^4 + Z_{\rm I}(T_{\rm I},M_0) L_V
\end{displaymath} (2)

where AV is the ice albedo, $\epsilon$ its emissivity, LVits sublimation heat, $c_\odot$ is the solar constant, $z_\odot$the solar zenith angle, $\sigma $ is the Stefan constant, and M0 is the initial flow Mach number (to be found by integration of the EE), and the dimensionless parameter Fs+c < 1, representing diffuse illumination and heat exchange by conduction, is discussed below.

We were unsure whether we should consider the CO production to be correlated or anticorrelated with the ice abundance, or to be unrelated in any way. This affects: (1) the outflow direction as mentioned above; and (2) the distant photodissociative heating discussed below. In the following, when a CO spot is assumed, only two cases are considered. In the first case, the icy fraction f inside the spot is increased by an amount $\Delta f = (1 - f) {\rm e}^{-[(\delta/20)^2
ln2]}$, corresponding to pure ice at the spot center; this produces from the spot (under peak solar illumination) a production rate $\Delta Q_{\rm H_2O}$ equal to about 5% of the total water production rate $Q_{\rm H_2O}$. In the second case, f has the same value over the whole surface, including inside the CO spot. The reader should observe that many alternative assumptions could be made with respect to the H2O-CO mixture production, for instance that only part of the surface is icy, and that CO diffuses negligibly, partly or totally from the icy area. At this time, the effect of these assumptions on the results has not yet been studied.

The total H2O production rate $Q_{\rm H_2O}$ depends on the value of f. To obtain values in order of magnitude agreement with the observations (e.g., Biver et al. 1999a; Dello Russo et al. 2000), we adopted f = 0.15. The resulting values of $Q_{\rm H_2O}$ are given in Table 2. One sees that, when Fs+c is small, $Q_{\rm H_2O}$ varies by about a factor of 2 during the rotation - a consequence of the adopted shape. As Fs+c increases, the rotational modulation of $Q_{\rm H_2O}$ decreases. The analysis of the Biver et al. (1999a) observations of OH indicated - surprisingly - a negligible rotational modulation, thus favouring a high value of Fs+c; we return to this point below.

The present observations are extremely sensitive to the coma CO temperature and velocity fields, which themselves depend on the nucleus surface temperature field and on the H2O photodissociative heating. The surface dust temperature $T_{\rm d}$ of the sunlit areas is much greater than the ice temperature $T_{\rm I}$, so that most of the CO molecules leave the surface at $T_{\rm d}$, but it can be shown that, owing to the dominance of the H2O flux at these points, the H2O-CO mixture common initial velocity and temperature are determined by $T_{\rm I}$ alone. In the shadowed areas, the latent heat exchanges are small, and $T_{\rm d} \simeq T_{\rm I}$, hence the previous conclusion remains true, even though the H2O flux is not dominant. Thus, the ratio $T_{\rm I}$(sunlit)/$T_{\rm I}$(shadow) is expected to strongly influence the observed day-to-night CO velocity asymmetry. In the shadow, $T_{\rm I}$ depends on conductive heat exchanges with the nucleus interior, and radiative exchanges with the coma, the latter being dominated by dust-scattered solar light. We represent these effects by an ad hoc term in the surface energy budget Eq. (2): the net flux of energy input from these effects, expressed as a fraction of the subsolar solar flux. Several values were used in alteration, in the range Fs+c = 0.05-0.5. This value drastically affects the ice sublimation rate in the parts of the nucleus that are not directly exposed to sunlight, hence controls the day-to-night asymmetry of both the temperature $T_{\rm I}$ and the H2O flux distribution. For all solutions discussed here, the parameters f and Fs+care adjusted to obtain a mean total production rate $Q_{\rm
H_2O} \simeq 1.3 \times 10^{31}$ molecules s-1, in order of magnitude agreement with observations (e.g., Biver et al. 1999a; Dello Russo et al. 2000).

3.2 Hydrodynamic model of the coma


The H2O-CO mixture flow in the unusually dense Hale-Bopp coma was computed by solving 3D+t Euler equations (EE). These EE should in principle be solved with additional equations representing all significant physical processes known to occur inside the gas coma (i.e., momentum and energy exchanges with the dust, photodissociation, photoionization, and radiative heating or cooling of the molecules) and/or should include terms allowing for these effects (Rodionov et al. 2002). In 3D, this would already be a formidable task, far beyond the present state-of-the-art development of available codes, since it would require, for instance, (1) computing locally the radiative cooling by H2O rotational emission in an optically thick coma; (2) computing locally the photolytic heating caused by the non-local thermalization of the fast hydrogen atoms produced by photodissociation; and (3) computing in a self-consistent manner the coma attenuation of the far UV solar radiation driving the photodissociation. Using a 3D+t code, it would be even worse, because the computations should be performed at each time step. Thus we were forced to limit ourselves to an heuristic computation of the photodissociation of H2O and its resulting coma heating, known to severely affect the coma temperature and velocity fields (at 1 AU, the CO photodissociation scale length is $\sim$1 500 000 km, much larger than the simulated coma, and the excess energy of the dissociated atoms is much less than in the case of H2O, hence the effect can be neglected).

3.2.1 Heuristic photodissociation model

Classical photodissociation conservation equations for H2O, H, O, and OH were coupled with the EE governing the molecular mixture, and the local photodissociation rate $\beta(r,\phi_\odot)$ was assumed to depend on both the cometocentric distance r and the local solar zenith angle $\phi_\odot$. Using precise photodissociation sections, we compute the subsolar distance $R_{\rm UV} \simeq 900$ km, where the dissociation rate has decreased by 1/e, and then assume that $\beta(r,\phi_\odot) = 0$ inside a ``finger-shaped'' volume consisting, on the day side, of a half sphere of radius $R_{\rm UV}$ and, on the night side, of its cylindrical shadow volume. Outside this volume, we assume that $\beta(r,\phi_\odot) = 1.3 \times 10^{-5}$ s-1. For the resulting heating rate, we write
                     $\displaystyle \dot{E}(r,\phi_\odot)$ = $\displaystyle E_{\rm H} \left(1-{\rm e}^{- \frac{R_{\rm H}}{r}(\alpha \cos{\phi_\odot} +1 -\alpha)}\right)$  
  $\textstyle \quad +$ $\displaystyle E_{\rm OH} \left(1-{\rm e}^{-\frac{R_{\rm OH}}{r}(\alpha \cos{\phi_\odot} + 1 -\alpha) }\right),$  

an equation inspired by the 1-D Monte Carlo simulation of the H and OH cooling in spherical geometry (
Bockelée-Morvan & Crovisier 1987), in which $E_{\rm H}$ and $E_{\rm OH}$ are the excess energy carried by the dissociation fragments, the distances $R_{\rm H}$ and $R_{\rm OH}$ control the variation of the photolytic heating with r, and the parameter $\alpha$ controls its variation with  $\phi_\odot$. Only 3-D modelling of fast H decelerating inside the coma could yield precise values of $\alpha$. Hence, we have run solutions with different assumed values of $\alpha$. Small values yield a similar heating of the day and night coma, whereas higher values yield a dominant dayside heating.

4 Molecular excitation and radiative transfer

The model is detailed in Boissier et al. (2007) and was first used to interpret the PdBI observations of sulfur-bearing species in comet Hale-Bopp. We describe here its application to the observations of the CO J(2-1) 230 GHz line. In summary, it requires: 1) the computation of the population of the J = 1 and J = 2 CO rotational levels using an excitation model; 2) the computation of the brightness distribution of the CO 230 GHz line using a radiative transfer model and the output of the hydrodynamical simulations; and 3) the simulation of the PdBI observations (ON-OFF spectra and interferometric data).

4.1 Excitation model

We used the excitation model described in Biver et al. (1999b), which accounts for neutral and electron collisions, and for radiative pumping of the CO vibrational bands by the solar radiation. The total cross-section for de-excitation by neutral collisions $\sigma_{\rm c}$ is assumed to equal $2 \times 10^{-14}$ cm2. The electron density is estimated from that derived from the ``Giotto'' measurements in comet 1P/Halley, applying physically justified scaling laws to allow for its dependence on the water production rate and heliocentric distance (Biver et al. 1999b). In other words, the $x_{\rm ne}$ parameter introduced by Biver et al. (1999b) to account for uncertainties in the absolute value of the electron density is taken to equal 1.

This 1-D excitation model integrates the level population budget equations along radial directions in the coma, which is justified because the present hydrodynamical simulations do not exhibit significant deviations from a radial flow in the regions sampled by the interferometric and ON-OFF beams (Figs. 5-7). Most of the CO molecules in the field of view are found to have their J = 1 and J = 2 rotational levels populated according to a local thermal equilibrium (Bockelée-Morvan et al. 2009) as a result of high collision rates combined with long radiative lifetimes of the CO rotational levels. Both $\sigma_{\rm c}$ and $x_{\rm ne}$ are not critical parameters in these calculations.

We did not use the 3-D distribution of the gas kinetic temperature given by the hydrodynamical model, for the following reasons: (1) the computed values (150-250 K at 1000-10 000 km) are significantly higher than implied by the observational constraints; (2) it would require the development of a 3-D excitation model, which is beyond the scope of this paper; (3) last but not least, it is very unlikely that the temporal variations that we attempt to reproduce are caused only by excitation effects.

Instead of the 3-D temperature, we assumed a constant kinetic temperature of 120 K, in agreement with the values derived for radial distances $r \sim 10~000$ km from the relative intensities of millimetre molecular lines, including CO lines, assuming an isothermal coma (Bockelée-Morvan et al. 2009; Biver et al. 2002). This temperature is also consistent with rotational temperatures measured from infrared ro-vibrational spectra of several molecules (Brooke et al. 2003; Disanti et al. 2001). Nevertheless, infrared observations imply lower temperatures ($\sim$70-80 K) in the inner ( $r \sim 1000$ km) regions (Disanti et al. 2001), as do the present gasdynamical results. These results also predict day-to-night asymmetries in the gas temperature, as well as lower temperatures in the spiral CO structure (by 30% with respect to the ambient gas in the outer regions; Figs. 67). The gas kinetic temperature affects the population of the CO rotational levels, and the opacity of the CO lines. Both the populations of the J = 1, 2 levels and the opacity of the CO J(2-1) line increase with decreasing kinetic temperature. In optically thin conditions, a lower temperature in gas structures (e.g., by 20% as found for solutions #2 to #5, Sect. 5) would result in relatively higher amounts of CO emission (by 20%) being produced by these structures than expected for an isothermal coma at 120 K. However, for optically thick conditions applying in the inner (1000 km) coma, this is counterbalanced by more significant opacity effects.

\begin{figure}
\par\includegraphics[%
height=11cm,width=6.4cm,clip]{12804f5a-NEW...
...cludegraphics[%
height=11cm,width=6.4cm,clip]{12804f5d-NEW.eps}\\\end{figure} Figure 5:

Solution #1: H2O and CO number density (molecules m-3), velocity and temperature in the plane perpendicular to the rotation axis at 6 different rotation phases (the corresponding time is given in the top left corner). The Sun direction is the +X axis. The white lines with arrows are flow lines.

Open with DEXTER

\begin{figure}
\par\includegraphics[%
height=11cm,width=6.3cm,clip]{12804f6a-NEW...
...%
height=11cm,width=6.3cm,clip]{12804f6d-NEW.eps}\\\vspace*{4mm}
\end{figure} Figure 6:

Solution #2: H2O and CO number density (molecules m-3), velocity and temperature in the plane perpendicular to the rotation axis at 6 different rotation phases (the corresponding time is given in the top left corner). The Sun direction is the +X axis. The arrow is the projected axis ${\vec {CM}}$ of the spot. The white lines with arrows are flow lines.

Open with DEXTER

\begin{figure}
\par\includegraphics[height=11cm,width=7cm,clip]{12804f7a-NEW.eps...
...ics[height=11cm,width=7cm,clip]{12804f7d-NEW.eps}\\\vspace*{5mm}
\end{figure} Figure 7:

Solution #5: H2O and CO number density (molecules m-3), velocity and temperature in the plane perpendicular to the rotation axis at 6 different rotation phases (the corresponding time is given in the top left corner). The Sun direction is the +X axis. The arrow is the projected axis ${\vec {CM}}$ of the spot. The white lines with arrows are flow lines.

Open with DEXTER

4.2 Gridding and integration

In computing the brightness distribution, we divided the plane of the sky into $256~\times~ 256$ square pixels of 0.5 $\hbox {$^{\prime \prime }$ }$ width (362 km in the coma). To adapt the grid to the variations in the local parameters, the pixels were divided into subpixels whose number depends on the distance to the map center: $50 \times 50$ for the central pixel, $20 \times 20$ for the pixels at less than 1000 km, $5 \times 5$ for the pixels at less than 5000 km, and $2 \times 2$ for the outermost pixels. These numbers were empirically determined to obtain enough precision in a reasonable CPU time. Each subpixel defines a line of sight along which we solved the radiative transfer equation (self-absorption and stimulated emission included) for 40 velocity channels i separated by 0.1 km s-1, corresponding to the channel spacing in acquired spectra. A synthetic brightness spectrum Fi is obtained for each pixel after summation of the fluxes in the subpixels.

The ON-OFF synthetic spectrum is obtained by convolving the brightness distribution with the primary beam pattern of the 15-m PdBI antennas, described by a 2D Gaussian. Visibilities were calculated from the Fourier transform of the brightness distribution (Eq. (1)).

4.3 Time variations

The Plateau de Bure interferometer recorded ten values of visibility (one per pair of antennas) every 45 s-scan. A full simulation of the observations would require us to compute for each scan the expected visibilities corresponding to the true state of the coma. This would require the computation of the brightness distribution on the plane of the sky using the 3D model output corresponding to the rotational phase of the nucleus at the time the scan was recorded. Following the method outlined by Bockelée-Morvan et al. (2009), we used only 12 snapshots to represent an entire nucleus rotation. We thus assumed that the brightness distribution does not evolve much within 1 h, which seems reasonable given the observed time variations.

4.4 Geometrical considerations

The results of the physical model depend on the relative positions of four important directions: (1) the enhanced flux mean direction; (2) the comet-Sun axis; (3) the rotation axis of the nucleus; and (4) the line of sight.

In the phenomenological model of Bockelée-Morvan et al. (2009), the Gaussian enhanced flux emanates from a spherical cap with best-fit axis 20 degrees above the equator. In the physical model, (a) the surface is neither a spherical cap, nor orthogonal to its axis; (b) the flux vector is at each point orthogonal to the local surface element; (c) hence, the flow does not originate in a virtual point source, as in the phenomenological model; and (d) its mean direction needs not coincide with the axis of the Gaussian. Thus, there are (unavoidable) significant differences between the postulated phenomenological model source and the postulated physical model source.

\begin{figure}
\par\includegraphics[width=6.7cm,clip]{12804f8a-NEW.eps}\hspace*{...
...ace*{2mm}
\includegraphics[width=6.7cm,clip]{12804f8d-NEW.eps}\\
\end{figure} Figure 8:

Solution #3: H2O and CO number density (molecules m-3), velocity and temperature in the plane perpendicular to the rotation axis at 6 different rotation phases (the corresponding time is given in the top left corner). The Sun direction is the +X axis. The arrow is the projected axis ${\vec {CM}}$ of the spot. The white lines with arrows are flow lines.

Open with DEXTER
\begin{figure}
\par\includegraphics[width=6cm,clip]{12804f9a-NEW.eps}\includegra...
...f9c-NEW.eps}\includegraphics[width=6cm,clip]{12804f9d-NEW.eps}\\
\end{figure} Figure 9:

Solution #4: H2O and CO number density (molecules m-3), velocity and temperature in the plane perpendicular to the rotation axis at 6 different rotation phases (the corresponding time is given in the top left corner). The Sun direction is the +X axis. The arrow is the projected axis  ${\vec {CM}}$ of the spot. The white lines with arrows are flow lines.

Open with DEXTER

In the physical model, the Sun direction and the spin axis are orthogonal whereas at the time of observations and in the phenomenological model their angle is 140$^{\circ}$, based on the rotation axis parameters derived by Jorda et al. (1999) from optical observations of the dust coma. When computing simulated data from the outputs of the physical model, we fixed the line of sight to be along a direction close to the true phase angles (46$^{\circ}$ for the Sun, and 80$^{\circ}$ for the angle between the spin axis and the line of sight), so that the velocity projections along the line of sight are consistent with observations. Therefore, simulated and observed spectra can be compared. As a consequence of the above-mentioned geometrical differences, the projected directions of all the axes of the physical model on the plane of the sky do not match the true directions. Hence, the brightness maps derived from the physical model are expected to exhibit built-in differences from the observed ones.

5 Results

In the phenomenological model, different solutions were attempted in fitting the enhanced source geometry and intensity to obtain a best-fit solution; here, the geometry remained unchanged, but several different parameter combinations were used: we investigated the effect of the presence or not of a localized CO flux increase  $\Delta F_0$, the presence or not of a bright diffuse illumination Fs+c, the presence or not of a day-night asymmetry in coma heating (playing with $\alpha$), and the presence or not of an increase in the water flux ($\Delta f$) coincident with  $\Delta F_0$.

5.1 The various computed physical solutions

Five different 3D+t gasdynamical solutions are discussed here; they differ in terms of their values of the parameters $F_{s+c}, \alpha$ and the CO flux and icy area fraction assumed inside the spot for increased production rates, resulting in different values of the spot CO production  $\Delta Q_{\rm CO}$ and spot H2O excess production for normal illumination  $\Delta Q_{\rm H_2O}$. The corresponding values, as well as the total production rates, are given in Table 2.

\begin{figure}
\par\includegraphics[angle=-90,width=17cm,clip]{12804f10.ps}
\end{figure} Figure 10:

Time evolution of the synthesized spectrum obtained with solutions #1 to #5. The presented spectra were computed at rotational phases of the nucleus comparable to those of the observed spectra, presented in Fig. 1 and reproduced here in the first column. The intensity of the spectra from solution #1 have been divided by a factor of 3 to make the comparison with other solutions easier.

Open with DEXTER

Solution #1 - The assumed dust light-scattering coefficient Fs+c is small, leading to strong diurnal variations in the surface temperature and H2O production. The production of CO is homogeneous (no spot). The coeffiicient $\alpha$ is 0.45, producing a dominantly sunward coma photodissociative heating. This is the solution described by Rodionov & Crifo (2006).

Solution #2 - Same as #1, except that 30% of $Q_{\rm CO}$ is released in a small area representing 2% of the total surface; outside of it, the production is uniform.

Solution #3 - Same as #2, except that Fs+c is large, resulting in a small day-night asymmetry in the surface temperature and H2O production.

Solution #4 - Same as #3, except that there is some H2O production increase inside the spot with increased CO production.

Solution #5 - Same as solution #4, except that $\alpha$is assumed to be small, producing a more efficient gas heating of the night side of the coma.

\begin{figure}
\par\includegraphics[angle=-90,width=16.5cm,clip]{12804f11.ps}
\end{figure} Figure 11:

Synthetic average ON-OFF profiles obtained with the 5 solutions compared to the observed spectrum, shown in the upper-left frame. The simulated intensities in the case of solution #1 have been divided by a factor of 3 to make the comparison with other solutions easier.

Open with DEXTER

5.2 General properties of the physical solutions

A full physical description of the solutions #1 to #5 will be presented - together with additional solutions - in a future publication. Here, we only sketch their properties. The solutions #1, #2, and #5 are presented, respectively, in Figs. 5-7, which show, at six times separated by 1/6 of the rotation period, the distribution of the H2O and CO number density, the coma gas velocity, and the coma gas temperature in the equatorial plane. Similar figures for solutions #3 and #4 are available as online material (Figs. 8 and 9). However, it is possible to characterize them by refering to #2 and #5.

Focusing first on the H2O density, we see a clear dominance of the local time over the effect of rotation in all solutions. A two orders of magnitude day-to-night difference is evident for #1 and #2, which reduces to a factor $\simeq$2-3 in #3, #4, and #5, except in the vicinity of the antisolar axis, where the night side density is greater. The first effect is caused by the higher Fs+c, and the second effect is explained below.

The flow velocity distribution is much more uniform (at any distance and in any solution) than the density: there is at most a factor of two between the lowest and the highest values. This being said, we see that it is smaller on the night side, and that its variation appears to reflect the photodissociative heating in the coma, increasing with increasing distance, with decreasing heating anisotropy (compare #2 and #5), and having small values in the antisolar ``UV extinction shadow cylinder''. Because flux is preserved during the nearly radial outflow, the antisolar velocity minimum produces the above-mentioned antisolar H2O density maximum for solutions #3 to #5.

The CO density is clearly predominantly shaped by the rotation (i.e., by the nucleus topography and/or inhomogeneity of production), as is obvious for all solutions. When a homogeneous CO production is assumed, the density patterns (produced by the surface topography and temperature distribution) evolve in a non-intuitive way (solution #1). Spiral components exist, which are too faint to be visible on the density isocontours, but appearing clearly on column density maps if a radial structure enhancement algorithm is applied (Rodionov & Crifo 2006).

As soon as the spot with increased CO flux is postulated (solutions #2 to #5), a spectacular spiral density structure appears, which rotates with the nucleus. The density enhancement in this structure is comparable to that postulated for the spot surface flux.

Finally, we consider the temperature field. As expected, it is sensitive to all parameters. All solutions contain a cold ``UV extinction shadow cylinder'' and display a strong radial increase produced by photodissociative heating, whose azimuthal extent is controlled by $\alpha$. Simulations with solution #1 do not show significant temporal modulations of the coma temperature (Fig. 5). An unexpected feature of solutions #2 to #4 is a pronounced temperature minimum (an approcimately 20% decrease) inside the spiral CO structure (we note that no similar signature exists in the velocity field). This effect is probably caused by the photolytic heating rate (in Watt per molecule) being practically constant over a distance comparable to the spiral structure width, while the gas density strongly increases when the spiral ``passes by''. This will be analyzed further in a future gasdynamical publication.

Regarding the order of magnitude of the computed temperatures, in the inner region where photolytic heating does not dominate it appears to be in rough agreement with the previously mentioned estimates based on IR spectra (about 70 K near 1000 km). At greater distances where the heating is more dominant, all presently computed dayside temperatures are about twice as great as than estimates (about 120 K at 10 000 km) inferred from single-dish observations of millimeter molecular lines, including CO lines. Since these estimates assume an isothermal coma, it would worthwhile to use physical 3D models to obtain a clearer interpretation of the temperature data. This being said, the agreement at small distances is encouraging, since the present model does not allow for gas-dust energy exchanges, even though these exchanges surely occur close to the nucleus in this very dusty comet, either directly (convective exchanges) or indirectly (by absorption of the thermal emission of the dust in the IR followed by collisional de-excitation). As to the factor of 2 discrepancy at larger distances, it is not surprising, given that the present gas energy budget model is, as stated above, primitive. The photolytic heating is modelled by a simple heuristic formula, in the absence of a 3D dissociation fragments cooling model. Molecular radiative cooling occurs, caused by the rotational emission of water and possibly other species. Because the water rotational lines are optically thick, the efficiency of this process must be evaluated using a 3D radiative transfer code and by taking into account the coupling between gas thermodynamics and excitation (Bockelée-Morvan & Crovisier 1987).

\begin{figure}
\par\includegraphics[angle=-90,width=15.5cm,clip]{12804f12.ps}
\end{figure} Figure 12:

Time evolution of the velocity shift $\Delta $v deduced from the observed ON-OFF spectrum ( top left panel, also shown in Figs. 1 and 2), and computed from outputs of solutions #1 to #5. For solution #1, the evolution of $\Delta $v is drawn with two scales. The scale in grey color corresponds to a zoom by a factor of 20 with respect to the main scale.

Open with DEXTER

5.2.1 CO line fluxes

The ON-OFF ( $F_{\rm SD}$) and interferometric $F_{\rm Int}$ fluxes measured in the synthetic spectra and maps differ from the observed values (Table 1). This is because the total CO production rate in the different solutions overestimates or underestimates the true value. For solutions #2 to #5, the averaged (over velocity) opacity of the CO line at the center of the interferometric map is 0.7, whereas the ON-OFF flux is not significantly affected by opacity effects. Opacity effects are less significant for solution #1. Based on the observed ON-OFF flux, the total CO production rate deduced from the 3D+t gasdynamical solutions #3 to #5 is in the range $1.91{-}2.05 \times 10^{30}$ molecules s-1 and consistent within 10% with the value of $2.1 \times 10^{30}$ molecules s-1deduced by Bockelée-Morvan et al. (2009). The small discrepancy can be attributed to differences in the CO velocity field between these solutions and the adopted constant value of 1.05 km s-1 in Bockelée-Morvan et al. (2009). For solution #2, the retrieved $Q_{\rm CO}$ is lower ( $1.65 \times 10^{30}$ molecules s-1), which is due to the lower CO velocities in the night side for this solution.

5.2.2 Line profiles

Figure 10 shows the individual ON-OFF spectra computed at different times for solutions #1 to #5. These synthetic spectra are analogous to the observed spectra shown in Fig. 1a. The corresponding rotation-averaged spectra (analogous to the observed average spectrum in Fig. 1b) are shown in Fig. 11. Visual comparison of the figures clearly indicates that solutions #1 and #2 are inadequate. Interestingly, one has a uniform CO production, whereas the other has a localized strong production enhancement. On the other hand, they have similarly strong day-to-night contrasts in nucleus surface temperature. As a result, though the spectra obtained with solutions #1 and #2 differ, they both peak permanently at a positive Doppler velocity. In contrast, profiles obtained with solutions #3 to #5 show temporal variations qualitatively consistent with the measurements. The time evolution of the velocity shift obtained for the different solutions is presented in Fig. 12. We now compare the quantities defined in Sect. 2.4 and presented in Table 1.

\begin{figure}
\par\includegraphics[angle=-90,width=16cm,clip]{12804f13.ps}\vspace*{3mm}
\end{figure} Figure 13:

Successive positions labelled 1 to 5 (relative to the mean position) of the brightness center of the observed interferometric maps ( top left panel, also shown in Fig. 3), and computed from solutions #1 to #5. The arrow indicates the direction of the north pole of the comet, projected onto the sky. The grey line represents $\Delta O$, one of the quantities used to compare the model results to the observations.

Open with DEXTER

Amplitude A. Rotation-induced velocity-shift variations are present for solution #1 (homogeneous production of CO), but have a very weak amplitude of variation: $A_{\rm 1} = 0.009$ km s-1 compared to $A_{\rm obs} \sim 0.4$ km s-1 (hereafter, the subscripts refer to either the observations (obs) or the solutions). For the solutions with a localized enhancement in CO production (solutions #2 to #5), the computed amplitudes of the velocity shift variation are in the range $A_{\rm 2{-}5} = 0.25{-}0.31$ km s-1 and are $\sim$30% smaller than the observed value.

Velocity shift on the averaged spectrum $\langle \Delta v \rangle$. Different values are obtained for the different solutions (Table 1). For solutions #1 and #2, where strong diurnal variations in nucleus surface temperature are present, $\langle \Delta v \rangle$ has an incorrect sign. The other solutions, which provide a small day-to-night temperature contrast, yield $\langle \Delta v \rangle$ values of the correct sign and order of magnitude. The closest agreement is however obtained for solutions #3 and #4, where a larger day-to-night difference in photolytic heating is present than for solution #5.

Velocity range ($v_{\min}$, $v_{\max}$). We can see from Table 1 that $v_{\min}$ and $v_{\max}$ are strictly identical in all solutions, except for solution #5, where the gas heating on the night side of the coma is higher than for other solutions. The synthetic ON-OFF spectra obtained with solution #5 cover a velocity range of 4 km s-1, which is consistent with that of the observed spectra.

5.2.3 Motion of the image brightness center

The time evolution of the peak brightness position provided by solutions #1 to #5 is presented in Fig. 13. Rotation-induced motions are present in solution #1 with an ``homogeneous CO production'', but have a small amplitude ( $\Delta
O_{\rm 1} = 0.2\hbox{$^{\prime\prime}$ }$ versus $\Delta O_{\rm obs} \sim 0.8 \hbox{$^{\prime\prime}$ }$). The positions obtained for solutions #2 to #5 are similar to first approximation with an amplitude $\Delta O_{\rm 2{-}5} =
1.1{-}1.5\hbox{$^{\prime\prime}$ }$. These positions cannot be blindly compared to the observed values since the geometry of the model differs from the true geometry (as mentioned in Sect. 4.4). The important points are that: (1) the overall amplitude of the motion is comparable to the measurements; and (2) the main axis of the motion is perpendicular to the rotation axis, as requested.

5.2.4 Interferometric spectra

Figure 14 presents the computed time evolution of the velocity shift measured in the visibility spectra obtained with the ten baselines; this is to be compared with the corresponding ten curves of Fig. 2. Solutions #2 to #5 yield curves that are roughly in phase, as does the phenomenological model (Bockelée-Morvan et al. 2009, Fig. 2). Hence, they do not reproduce the observed phase shifts. Their amplitudes of variation are, for some baselines, consistent with the observed values, but, for the other baselines, are about twice as great as the observed values. As opposed to this, solution #1 provides curves that are not in phase (curves observed in phase correspond to parallel or almost parallel baselines), but whose amplitudes are about one half to one fourth of the observed ones. Because of their different orientations, the different baselines probably do not probe the fine shock structures spiralling into the inner coma equivalently. On the other hand, the strong CO structure included in solutions #2 to #5 governs the time variations and does not produce significant phase shifts.

\begin{figure}
\par\includegraphics[angle=0,width=12.5cm,clip]{12804f14.ps}
\end{figure} Figure 14:

Time evolution of the velocity shift measured in the synthesized visibility spectra for all baselines. The antenna pairs (i.e., baselines) are indicated in the upper-right corner, and sorted by increasing lengths.

Open with DEXTER

6 Discussion

We have simulated the Plateau de Bure observations of the CO 230 GHz line in comet Hale-Bopp using five different time-dependent hydrodynamical calculations. It is clearly not easy to identify one solution that could be considered the closest description of the observations. As we have seen, all solutions have their merits and deficiencies. We however observe that, no attempt has been made in this work to optimize the parameters of the CO jet that affect the magnitude of the A, $\langle \Delta v \rangle$, and $\Delta O$ quantities (Bockelée-Morvan et al. 2009). Most importantly also, even though the present coma model is already fairly advanced, improvements are clearly needed in the computation of several physical coma processes. In principle, it is even conceivable to adopt other assumptions for the surface production of CO and even of H2O.

The present results show that, although a rotating realistic (i.e., aspherical) nucleus with a uniform surface gas production of CO does create spiralling structures in the coma, they are too faint to explain the observed rotational modulation of the ON-OFF spectra and the observed motion of the brightness peak in the interferometric maps. We confirm that the assumptions of a constant CO flux over most (98%) of the surface of the nucleus, the rest of the surface (a narrow cone of 18 degree half-angle low latitude) having a much higher flux representing $\sim$30% of the total CO production, can explain a large part of the observations. The data are more accurately reproduced when the day/night asymmetry of the models in terms of H2O production and gas heating is weak. This implies little rotational modulation of $Q_{\rm H_2O}$ (Table 1). Interestingly, short-term temporal variations in the water production rate have not been reported in the literature. The lines of several molecules (CS, H2S, HCN) with volatilities intermediate between those of CO and H2O were thoroughly monitored over almost one nucleus period in March 1997, using the PdBI (Bockelée-Morvan et al. 2005; Boissier et al. 2007). No significant intensity variations were observed.

The failure to explain the relative phase shifts between the interferometric velocity offset curves is puzzling. From a thorough analysis of the CO PdBI data, Bockelée-Morvan et al. (2009) identified a second moving structure, possibly associated with an outburst, and suggested that this structure might have severely affected the signals recorded by some baselines of specific orientation, causing the observed delays.

When dealing with coma-integrated observations, it is relatively easy to infer global nucleus properties (e.g., molecular production rates) using conservative laws. But here, one deals, for the first time, with spatially and rotationally resolved observations of a parent molecule. The objective of extracting from them spatial information about the nucleus surface flux can be respected. However, it appears that, at the same time, an advanced understanding of the coma must be achieved, an objective of substantial interest since, in contrast to the case of the nucleus, a large set of observational constraints result from the large variety of coma data. This requires us to develop physically robust methods for modelling the coma, because of the non-intuitive and non-trivial nature of fluid flows. This is illustrated, for instance, by the extreme sensitivity of the computed radiative properties to the photodissociation - the only non-adiabatic effect allowed for inside the model. This sensitivity is itself a consequence of the high sensitivity of the rotational spectra to the temperature. Hence, not only would a thorough 3-D modelling of the excitation and radiation transfer be required itself, but also advanced allowance for all non-adiabatic effects operating in Hale-Bopp's coma.

Acknowledgements
We are grateful to Vladimir Zakharov for his help during the writing of this paper. This work has been supported by the Programme national de planétologie of Institut national des sciences de l'univers.

References

  1. Altenhoff, W. J., Bieging, J. H., Butler, B., et al. 1999, A&A, 348, 1020 [NASA ADS] [Google Scholar]
  2. Biver, N., Bockelée-Morvan, D., Colom, P., et al. 1999a, Earth Moon and Planets, 78, 5 [Google Scholar]
  3. Biver, N., Bockelée-Morvan, D., Colom, P., et al. 2002, Earth Moon and Planets, 90, 5 [NASA ADS] [CrossRef] [Google Scholar]
  4. Biver, N., Bockelée-Morvan, D., Crovisier, J., et al. 1999b, ApJ, 118, 1850 [Google Scholar]
  5. Blake, G. A., Qi, C., Hogerheijde, M. R., Gurwell, M. A., & Muhleman, D. O. 1999, Nature, 398, 213 [NASA ADS] [CrossRef] [PubMed] [Google Scholar]
  6. Bockelée-Morvan, D., & Crovisier, J. 1987, In Diversity and Similarity of Comets, ESA-SP 278, 235 [Google Scholar]
  7. Bockelée-Morvan, D., Gunnarsson, M., Biver, N., et al. 2005, in Asteroids, Comets, and Meteors, Book of Abstracts (Ithaca: Cornell Univ.), 111 [Google Scholar]
  8. Bockelée-Morvan, D., Henry, F., Biver, N., et al. 2009, A&A, 505, 82 [Google Scholar]
  9. Boissier, J., Bockelée-Morvan, D., Biver, N., et al. 2007, A&A, 475, 1131 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  10. Boissier, J., Bockelée-Morvan, D., Biver, N., et al. 2008, LPI Contributions, 1405, 8081 [NASA ADS] [Google Scholar]
  11. Brooke, T. Y., Weaver, H. A., Chin, G., et al. 2003, Icarus, 166, 167 [NASA ADS] [CrossRef] [Google Scholar]
  12. Crifo, J. F., Fulle , M., Kömle , N. I., et al. 2005, In Comets II, ed. M. Festou, H. U. Keller, & H. Weaver (Tucson: Univ. of Arizona Press), 471 [Google Scholar]
  13. Dello Russo, N., Mumma, M. J., DiSanti, M. A., et al. 2000, Icarus, 143, 324 [NASA ADS] [CrossRef] [Google Scholar]
  14. Disanti, M. A., Mumma, M. J., Dello Russo, N., et al. 2001, Icarus, 153, 361 [NASA ADS] [CrossRef] [Google Scholar]
  15. Fernández, Y. R. 2002, Earth Moon and Planets, 89, 3 [Google Scholar]
  16. Henry, F. 2003, Ph.D. Thesis, Thèse de doctorat de l'Université Pierre et Marie Curie, Paris [Google Scholar]
  17. Henry, F., Bockelée-Morvan, D., Crovisier, J., et al. 2002, Earth Moon and Planets, 90, 57 [Google Scholar]
  18. Jorda, L., Rembor, K., Lecacheux, J., et al. 1999, Earth Moon and Planets, 77, 167 [Google Scholar]
  19. Lederer, S. M., & Campins, H. 2002, Earth Moon and Planets, 90, 381 [Google Scholar]
  20. Lederer, S. M., Campins, H., Osip, D. J., et al. 1999, Earth Moon and Planets, 78, 131 [NASA ADS] [CrossRef] [Google Scholar]
  21. Lederer, S. M., Campins, H., & Osip, D. J. 2009, Icarus, 199, 484 [NASA ADS] [CrossRef] [Google Scholar]
  22. Prialnik, D., Benkhoff, J., & Podolak, M. 2004, Comets II, 359 [Google Scholar]
  23. Rodionov, A. V., Crifo, J.-F., Szego, K., Lagerros, J., & Fulle, M. 2002, Planet. Space Sci., 50, 983 [NASA ADS] [CrossRef] [Google Scholar]
  24. Rodionov, A. V., Jorda, L., Jones, G. H., et al. 1998, Icarus, 136, 232. [NASA ADS] [CrossRef] [Google Scholar]
  25. Rodionov, A. V., & Crifo, J.-F. 2006, Adv. Space Res., 38, 1923 [NASA ADS] [CrossRef] [Google Scholar]
  26. Snyder, L. E., Veal, J. M., Woodney, L. M., et al. 2001, ApJ, 121, 1147 [Google Scholar]
  27. Szegö, K., Crifo, J.-F., Rodionov, A. V., et al. 2002, Earth Moon and Planets, 90, 435 [CrossRef] [Google Scholar]
  28. Thompson, A., Moran, J., & Swenson, G. 1991, Interferometry and Synthesis in Radio Astronomy, ed. Krieger Publishing Company [Google Scholar]
  29. Veal, J. M., Snyder, L. E., Wright, M., et al. 2000, ApJ, 119, 1498 [Google Scholar]
  30. Wink, J., Bockelée-Morvan, D., Despois, D., et al. 1999, Earth Moon and Planets, 78, 63 [Google Scholar]
  31. Woodney, L. M., Schleicher, D. G., & Bair, A. N. 2008, LPI Contributions, 1405, 8316 [NASA ADS] [Google Scholar]
  32. Wright, M. C. H., de Pater, I., Forster, J. R., et al. 1998, ApJ, 116, 3018 [Google Scholar]

Footnotes

... simulations[*]
Based on observations carried out with the IRAM Plateau de Bure Interferometer. IRAM is supported by INSU/CNRS (France), MPG (Germany) and IGN (Spain).
... version[*]
Both 3D and 3D+t Godunov-based hydrodynamic codes use 3D+t governing equations; only the methods of solution (and the size of the solution!) differ.

All Tables

Table 1:   Main characteristics of ON-OFF spectra and interferometric maps for observed and simulated data (Sect. 2.4).

Table 2:   Model parameters for the different solutions.

All Figures

  \begin{figure}
\par\includegraphics[angle=-90,width=8.5cm,clip]{12804f1.ps}
\end{figure} Figure 1:

ON-OFF CO J(2-1) (230.538 GHz) observations in comet Hale-Bopp on March, 11, 1997: a) seven representative individual ON-OFF spectra (1 min of integration ON source) recorded from 4.42 to 13.76 UT, (out of the 13 observed, Bockelée-Morvan et al. 2009; Henry et al. 2002). The intensity unit is the main beam brightness temperature. The UT time is given in the upper right corner; b) average of the 13 spectra recorded from 4.42 to 13.76 UT; c) line area as a function of UT time. The dashed line indicates the mean value; d) line velocity shift as a function of time. The dashed line indicates the value measured on the average spectrum.

Open with DEXTER
In the text

  \begin{figure}
\par\includegraphics[width=8.5cm,clip]{12804f2.ps}
\end{figure} Figure 2:

Time evolution of the velocity shift measured on ON-OFF and visibility spectra. The antenna pairs (i.e., baselines) are indicated and sorted by increasing lengths (e.g., 3-4: 32 m; 2-4: 88 m, 3-5: 101 m; 2-5 : 147 m). Observations are shown by filled circles with error bars, and the fitted sinusoids are shown by dashed lines. Open squares and solid curves are the result of the phenomenological model fit. Figure adapted from Bockelée-Morvan et al. (2009).

Open with DEXTER
In the text

  \begin{figure}
\par\includegraphics[angle=-90,width=9cm,clip]{12804f3.ps}
\end{figure} Figure 3:

Successive CO 230 GHz line integrated maps obtained from the data subsets covering: 1/ 4 h30-5 h30 UT, 2/ 6 h00-7 h00 UT, 3/ 7 h00-8h00 UT, 4/ 8 h50-9 h50 UT, and 5/ 10 h45-11 h45 UT. Twenty-nine velocity channels were averaged. The first positive contour and the contour spacing is 0.2 mJy (2-$\sigma $). The peak intensity at the center of the maps in units of line area is 6.76, 5.96, 5.40, 5.46, and 3.81 Jy km s-1 for maps 1 to 5, respectively. The beam shape and size is shown in the bottom left corner. The (RA, Dec) = (0 $^{\prime \prime }$, 0 $^{\prime \prime }$) coordinates correspond to the mean photometric center determined from the whole data set and coincides approximately to the nucleus position (Bockelée-Morvan et al. 2009). The projected Sun direction is indicated by an arrow. The offset between the brightness peak position measured on these maps and the position measured on the map that includes the whole data set is shown in the bottom right frame. The positions were obtained by fitting the Fourier transform of a Gaussian to the visibilities (i.e., the fit is performed in the Fourier uv plane). The positional accuracy is about 0.1 $\hbox {$^{\prime \prime }$ }$. The arrow on the lower-right panel indicates the direction of the nucleus north pole, projected on the sky. In the same panel, the grey line represents $\Delta O$, the distance between the most distant individual points. Figure adapted from Henry et al. (2002) and Bockelée-Morvan et al. (2009).

Open with DEXTER
In the text

  \begin{figure}
\par\includegraphics[width=9cm,clip]{12804f4.ps}
\vspace*{3mm}
\end{figure} Figure 4:

Nucleus used in the 3-D+t hydrodynamical simulations (Rodionov & Crifo 2006). The shape is that of comet Halley, smoothed and scaled-up to match the estimated size of comet Hale-Bopp. The illumination of the nucleus (grey scale) by direct Solar radiation is shown at two different rotation phases (separated by one third of the period). The rotation axis is vertical, the Sun is to the right.

Open with DEXTER
In the text

  \begin{figure}
\par\includegraphics[%
height=11cm,width=6.4cm,clip]{12804f5a-NEW...
...cludegraphics[%
height=11cm,width=6.4cm,clip]{12804f5d-NEW.eps}\\\end{figure} Figure 5:

Solution #1: H2O and CO number density (molecules m-3), velocity and temperature in the plane perpendicular to the rotation axis at 6 different rotation phases (the corresponding time is given in the top left corner). The Sun direction is the +X axis. The white lines with arrows are flow lines.

Open with DEXTER
In the text

  \begin{figure}
\par\includegraphics[%
height=11cm,width=6.3cm,clip]{12804f6a-NEW...
...%
height=11cm,width=6.3cm,clip]{12804f6d-NEW.eps}\\\vspace*{4mm}
\end{figure} Figure 6:

Solution #2: H2O and CO number density (molecules m-3), velocity and temperature in the plane perpendicular to the rotation axis at 6 different rotation phases (the corresponding time is given in the top left corner). The Sun direction is the +X axis. The arrow is the projected axis ${\vec {CM}}$ of the spot. The white lines with arrows are flow lines.

Open with DEXTER
In the text

  \begin{figure}
\par\includegraphics[height=11cm,width=7cm,clip]{12804f7a-NEW.eps...
...ics[height=11cm,width=7cm,clip]{12804f7d-NEW.eps}\\\vspace*{5mm}
\end{figure} Figure 7:

Solution #5: H2O and CO number density (molecules m-3), velocity and temperature in the plane perpendicular to the rotation axis at 6 different rotation phases (the corresponding time is given in the top left corner). The Sun direction is the +X axis. The arrow is the projected axis ${\vec {CM}}$ of the spot. The white lines with arrows are flow lines.

Open with DEXTER
In the text

  \begin{figure}
\par\includegraphics[width=6.7cm,clip]{12804f8a-NEW.eps}\hspace*{...
...ace*{2mm}
\includegraphics[width=6.7cm,clip]{12804f8d-NEW.eps}\\
\end{figure} Figure 8:

Solution #3: H2O and CO number density (molecules m-3), velocity and temperature in the plane perpendicular to the rotation axis at 6 different rotation phases (the corresponding time is given in the top left corner). The Sun direction is the +X axis. The arrow is the projected axis ${\vec {CM}}$ of the spot. The white lines with arrows are flow lines.

Open with DEXTER
In the text

  \begin{figure}
\par\includegraphics[width=6cm,clip]{12804f9a-NEW.eps}\includegra...
...f9c-NEW.eps}\includegraphics[width=6cm,clip]{12804f9d-NEW.eps}\\
\end{figure} Figure 9:

Solution #4: H2O and CO number density (molecules m-3), velocity and temperature in the plane perpendicular to the rotation axis at 6 different rotation phases (the corresponding time is given in the top left corner). The Sun direction is the +X axis. The arrow is the projected axis  ${\vec {CM}}$ of the spot. The white lines with arrows are flow lines.

Open with DEXTER
In the text

  \begin{figure}
\par\includegraphics[angle=-90,width=17cm,clip]{12804f10.ps}
\end{figure} Figure 10:

Time evolution of the synthesized spectrum obtained with solutions #1 to #5. The presented spectra were computed at rotational phases of the nucleus comparable to those of the observed spectra, presented in Fig. 1 and reproduced here in the first column. The intensity of the spectra from solution #1 have been divided by a factor of 3 to make the comparison with other solutions easier.

Open with DEXTER
In the text

  \begin{figure}
\par\includegraphics[angle=-90,width=16.5cm,clip]{12804f11.ps}
\end{figure} Figure 11:

Synthetic average ON-OFF profiles obtained with the 5 solutions compared to the observed spectrum, shown in the upper-left frame. The simulated intensities in the case of solution #1 have been divided by a factor of 3 to make the comparison with other solutions easier.

Open with DEXTER
In the text

  \begin{figure}
\par\includegraphics[angle=-90,width=15.5cm,clip]{12804f12.ps}
\end{figure} Figure 12:

Time evolution of the velocity shift $\Delta $v deduced from the observed ON-OFF spectrum ( top left panel, also shown in Figs. 1 and 2), and computed from outputs of solutions #1 to #5. For solution #1, the evolution of $\Delta $v is drawn with two scales. The scale in grey color corresponds to a zoom by a factor of 20 with respect to the main scale.

Open with DEXTER
In the text

  \begin{figure}
\par\includegraphics[angle=-90,width=16cm,clip]{12804f13.ps}\vspace*{3mm}
\end{figure} Figure 13:

Successive positions labelled 1 to 5 (relative to the mean position) of the brightness center of the observed interferometric maps ( top left panel, also shown in Fig. 3), and computed from solutions #1 to #5. The arrow indicates the direction of the north pole of the comet, projected onto the sky. The grey line represents $\Delta O$, one of the quantities used to compare the model results to the observations.

Open with DEXTER
In the text

  \begin{figure}
\par\includegraphics[angle=0,width=12.5cm,clip]{12804f14.ps}
\end{figure} Figure 14:

Time evolution of the velocity shift measured in the synthesized visibility spectra for all baselines. The antenna pairs (i.e., baselines) are indicated in the upper-right corner, and sorted by increasing lengths.

Open with DEXTER
In the text


Copyright ESO 2010

Current usage metrics show cumulative count of Article Views (full-text article views including HTML views, PDF and ePub downloads, according to the available data) and Abstracts Views on Vision4Press platform.

Data correspond to usage on the plateform after 2015. The current usage metrics is available 48-96 hours after online publication and is updated daily on week days.

Initial download of the metrics may take a while.