Low-frequency waves at comet 67P/Churyumov-Gerasimenko (cid:63) Observations compared to numerical simulations

Context. A new type of low-frequency wave was detected by the magnetometer of the Rosetta Plasma Consortium at the comet during the initial months after the arrival of the Rosetta spacecraft at comet 67P / Churyumov-Gerasimenko. This large-amplitude, nearly continuous wave activity is observed in the frequency range from 30mHz to 80mHz where 40mHz to 50mHz is the dominant frequency. This type of low frequency is not closely related to the gyrofrequency of newborn cometary ions, which di ﬀ ers from previous wave activity observed in the interaction region of comets with the solar wind. Aims. This work aims to reveal a global view on the wave activity region using simulations of the comet-solar wind interaction region. Parameters, such as wavelength, propagation direction, and propagation patterns, are within the focus of this study. While the Rosetta observations only provide local information, numerical simulations provide further information on the global wave properties. Methods. Standard hybrid simulations were applied to the comet-solar wind interaction scenario. In the model, the ions were described as particles, which allows us to describe kinetic processes of the ions. The electrons were described as a ﬂuid. Results. The simulations exhibit a threefold wave structure of the interaction region. A Mach cone and a Whistler wing are observed downstream of the comet. The third kind of wave activity found are low-frequency waves at 97 mHz, which corresponds to the waves observed by Richter et al. (2015, Ann. Geophys., 33, 1031). These waves are caused by the initial pick-up of the cometary ions that are perpendicular to the solar wind ﬂow and in the interplanetary magnetic ﬁeld direction. The associated electric current becomes unstable. The simulations show that wave activity is only detectable in the + E hemisphere and that the Mach cone and whistler wings need to be distinguished from the newly found instability driven wave activity.


Introduction
When a comet approaches the Sun, insolation leads to sublimation of ices and the emanating gas forms an extended cometary atmosphere, called the coma. The neutrals are ionised by means of photoionisation and/or charge exchange and then picked up by the solar wind and its magnetic field. This triggers an interaction that finally results in the formation of boundaries and structures, such as cometary bow shocks, diamagnetic cavities, and plasma tails (cf. Cravens & Gombosi 2004;Coates & Jones 2009;Koenders et al. 2013Koenders et al. , 2015Goetz et al. 2016).
In recent decades several cometary spacecraft missions, which were equipped with plasma instruments, visited comets and confirmed the presence of these structures and boundaries (cf. Coates & Jones 2009). Observations at previous comets, however, have been carried out when these comets were close to their maximum activity, that is close to the Sun, and had a large gas production rate Q ≈ 5 × 10 27 s −1 to 10 30 s −1 . This large gas production rate causes a vast comet-solar wind plasma interaction region that extends up to millions of kilometres upstream of the comet.
The movie associated to Fig. 10 is available at http://www.aanda.org During the arrival months (August to November 2014) of the Rosetta spacecraft at comet 67P/Churyumov-Gerasimenko (Glassmeier et al. 2007a), its activity is much lower than those of the previously visited comets (Hässig et al. 2015;Bieler et al. 2015): Q ≈ 1 × 10 26 s −1 . By this, the size of the main interaction region is expected to be much smaller than those at the more active comets. Moreover, the expected sizes are even smaller than the scales of the ion motion. Its size is determined by the gyroradius of the cometary ions r gyr,CI = m CI u SW,⊥ /q CI B, where m CI , q CI , and u SW,⊥ denote the mass and charge of the cometary ions and the solar wind velocity perpendicular to the magnetic field B. For typical solar wind parameters around 3 AU, the gyroradius of the cometary ions is of the order of 4 × 10 4 km, which exceeds the size of the obstacle by far. Thus, a different type of comet-solar wind plasma interaction is expected, and consequently, the measurements by the Rosetta Plasma Consortium (RPC; Carr et al. 2007) explored an uncharted type of plasma region during these first months.
The first indications of cometary plasma activity at comet 67P/Churyumov-Gerasimenko were obtained by the RPC Ion Composition Analyser (RPC-ICA). Nilsson et al. (2015a,b) found cometary ions that were accelerated perpendicular to the solar wind flow and parallel to the solar wind electric field, which did not yet populate a full ring in the velocity space, as predicted by numerical simulations (cf. Bagdonat & Motschmann 2002b;Rubin et al. 2014). In addition, it was found that the solar wind flow direction differs from the Sun direction close to the comet. Similar findings were made by Broiles et al. (2015), who, based on the RPC Ion Electron Sensor (RPC-IES) data, found that the deflection can reach up to 60 • and that changes in the magnetic field, as measured by the RPC Magnetometer (RPC-MAG; Glassmeier et al. 2007b), change the angle of the solar wind deflection. These observations clearly indicate that during its low-activity phase the interaction region can be divided into two regions depending on the solar wind electric field. In the following, we denote these two regions as +E hemisphere and −E hemisphere. The cometary ions are accelerated into the +E hemisphere and perform their cycloidal motion in it (Rubin et al. 2014). In contrast the solar wind is deflected into the −E hemisphere.
As the particle density in the cometary interaction region is extremely low, the interaction between the solar wind protons and cometary ions is not mediated by collisions. It is rather mediated by plasma waves, which have large amplitudes as observed at previous cometary spacecraft missions, e.g. comets 21P/Giacobini-Zinner ), 1P/Halley (Glassmeier et al. 1989), and 26P/Grigg-Skjellerup (Glassmeier & Neubauer 1993). At these comets the frequency in the spacecraft frame of reference was close to the local gyrofrequency of the cometary water group ions. Ringbeam and non-gyrotropic pick-up ion instabilities are responsible for the observed wave activity (e.g. Winske & Gary 1986;Brinca & Tsurutani 1987;Motschmann & Glassmeier 1998). Similar waves were expected at comet 67P/Churyumov-Gerasimenko, but with amplitudes that are much below 1 nT during its low-activity phase (Volwerk et al. 2013b,a). In contrast, RPC-MAG observed quasi-harmonic, large-amplitude (δB/B ∼ 2), and nearly continuous low-frequency wave activity in the vicinity of the comet 67P/Churyumov-Gerasimenko . These waves are neither purely compressional nor purely transverse, and their frequency is in the range from 30 mHz to 80 mHz, but groups around 40 mHz to 50 mHz. A clear dependence of the wave frequency on the local magnetic field does not exist, which strongly suggests that these observed low-frequency waves are different from the previously observed type of waves. Richter et al. (2015) suggest that a cross-current instability (Chang et al. 1990;Sauer et al. 1998) is generating this new type of cometary waves. Richter et al. (2015) report the first observation of the new type of waves at times around Rosetta's encounter with 67P/Churyumov-Gerasimenko in August 2014, and observe the waves everywhere in the inner coma within about 100 km distance to the nucleus. An example of this newly detected type of wave is presented in Fig. 1, which shows the three components of the magnetic field during a five minute interval on October 16, 2014. The spacecraft was at (0.1; 6.3; −7.7) km in the Cometcentred solar equatorial (CSEQ) frame of reference, about 10 km away from the nucleus. That day the magnetometer sampled with 20 Hz. A linear trend is removed and a low pass butterworth filter with a cut-off frequency of 2 Hz of order four is applied to the data to remove spacecraft noise and reaction wheel signatures. The bottom panel of Fig. 1 shows the trace of the spectra of the magnetic field between 02:00 UT and 03:00 UT on that day. A spectral enhancement at about 77 mHz is clearly visible. This frequency is in the range reported by Richter et al. (2015) as typical for the new type of wave observed. For comparison, we note that the local water group ion frequency was about 1 mHz. To further elucidate the properties of this new type of wave and its spatial structure within the overall solar wind-comet interaction region during the weak activity phase of 67P/Churyumov-Gerasimenko, we set up numerical simulations to gain a more global view of the oscillations. As a result of the small size of the interaction region, about 200 km, single fluid models that neglect the particle motion are not applicable to this type of interaction. This implies that more complex models, such as hybrid models that describe ions as particles and electrons as a fluid, have to be applied to this scenario. In this paper, we aim to conduct a qualitative comparison of the RPC measurements with the numerical plasma simulations and a study of the global properties of the newly detected cometary wave activity.
for example to weakly active comets (Bagdonat & Motschmann 2002b;Rubin et al. 2014) as well as to strongly active comets ).
Since it is not possible to describe all ions individually, the A.I.K.E.F. code gathers ions as macro-particles. These have the same charge-to-mass ratio as the real ions, whereby the amount of ions that are represented by a macro-particle is not fixed. During each timestep in the simulation the A.I.K.E.F. code solves the Newtonian equation of motion for each macro-particle i, where q i , m i , x i , and u i denote the charge, mass, position and velocity of the macro-particle. The electric field E is calculated with the momentum equation of the electrons. Under the assumption of massless electrons, quasi-neutrality, and a slowly changing electric field, it follows that where c,i , u i , and p e denote the charge density of the ions, bulk velocity of all ions, and electron pressure. The electron pressure is computed by solving the electron pressure equation where f e and u e denote the number of degrees of freedom for the electron fluid and bulk velocity of the electrons. Following Bagdonat & Motschmann (2002b), the number of degrees of freedom for the electrons is set to two. The right-hand side of the equation, ∂ t , denotes sources and sinks because of the different processes in the cometary atmosphere described later. The observed features in this simulation also occur if an adiabatic state equation for the electron fluid is used (Koenders et al. 2013). The magnetic field is obtained using the Faraday law The A.I.K.E.F. code is a particle-in-cell code. Thus, the macroparticles, at least 125 per cell in the presented simulations, can basically move to every position within the simulation box, whereas the electromagnetic fields and the moments, c,i and u i , are stored and computed on the nodes of a numerical mesh. The A.I.K.E.F. code uses a rectangular mesh without mesh refinement in this study. The mesh size in the presented simulation is about ∆x = 5 km = 0.025 x 0 , which is much below the scales of the features of interest, and the time step is dt = 2.6 × 10 −3 s, which fulfil the Courant-Friedrichs-Lewy criteria. The simulation box is a cube with an edge length of L = 1167 km = 6 x 0 ; x 0 denotes the ion inertia length in the undisturbed solar wind (Table 1). The moments and extrapolation of the fields to the macro-particles are performed by means of the cloud-in-cell technique. We solve the field equations by means of the current advance method and a cycling leapfrog (cf. Matthews 1994;Bagdonat & Motschmann 2002a). We apply a numerical smoothing to stabilise the simulation. Despite the fact that first measurements at 67P/Churyumov-Gerasimenko reveal that the coma is highly inhomogeneous  Hansen et al. (2007). The ion inertia length and the inverse gyrofrequency are given for the undisturbed solar wind.
and that the nucleus is not spherical (cf. Sierks et al. 2015;Hässig et al. 2015), we use the classic Haser model (Haser 1957) to reduce the complexity of the resulting structures and to simplify the analysis. In the Haser model, the neutral gas density is given by The neutral gas velocity is set to u ng = 800 m s −1 and r is the radial distance to the centre of the simulation box. Owing to the coarse grid size and small nucleus, the nucleus, which would be located in the centre of the simulation, is not modelled. The gas production rate in this scenario is set to Q = 8 × 10 25 s −1 and, for the sake of simplicity, only water is used as neutral gas species. The photoionisation ν of this neutral gas is modelled within the A.I.K.E.F. code by inserting new macro-particles into the simulation box. Each time step the new macro-particles are inserted equally distributed into the simulation box. However, the amount of ions which are represented by a new macro-particle is given by where ν = 1.37 × 10 −7 s −1 is the photoionisation rate (Hansen et al. 2007) and n ng the neutral gas density. The mass of the cometary water group ions is m ci = 19 amu = 19 · m sw ; m sw is the solar wind proton mass. In addition to ionisation, the simulation also considers collisions between ions and the neutral gas. This can lead to chargeexchange and a drag force on the plasma. Collision rates of k ci,ng = 1.7 × 10 −15 m 3 s −1 , for the cometary ions, and k sw,ng = 10 −14 m 3 s −1 , for the solar wind protons, are used according to Gombosi et al. (1996), and Kriegel et al. (2014). Furthermore, the simulation considers dissociative recombination; for details we refer to Koenders et al. (2015) because these collisions only have a minor impact on the plasma interaction of weakly active comets, as shown by Rubin et al. (2014).
The photoionisation is described by the first three terms. Here, only a mean energy of E ph ng = 12 eV is used to model the ionisation (Huebner et al. 1992;Rubin et al. 2014). The fourth term considers the recombination. Elastic ion-neutral collisions, according to Schunk & Nagy (2009), are modelled in the last term to cool the electron fluid.
In order to avoid problems of numerical stability, the neutral gas density is limited to the density value, which is reached at 2∆x. In addition, within a radius of 5∆x, the electron pressure is set constant to the mean value on a shell with a radius of 6∆x, which proved to be a suitable value.
The parameters of the simulation are summarised in Table 1. The set-up used here is identical with that chosen by Rubin et al. (2014). The simulation is representative for the situation of the plasma environment at 67P/Churyumov-Gerasimenko in the first months after the arrival of the Rosetta spacecraft. Actual parameters exhibit significant variations, which requires a scenario to be selected as computational resources are limited. Hence, a systematic parameter study cannot be performed.

Global structure of the interaction region
For this study we simulated the processes in the solar windcomet interaction region for a time interval of 1116 s. This corresponds to about twice the cometary ion gyroperiod and about 380 passage times, that is the time a solar wind parcel needs to pass through the simulation box. It took about 53 s until the simulation reached its quasi-stationary state. The initial state of the simulation is a box only filled with solar wind particles. With ongoing time, cometary ions are inserted and the interaction evolves until the quasi-stationary state is reached. For our analysis we use the remaining 1063 s. Figures 2 and 3 present the state of the simulation after 264 s. This time is selected arbitrarily.
In its quasi-stationary state, the solar wind density remains constant upstream of the centre of the neutral gas cloud (Figs. 2a and b). Downstream both a small thin enhancement with values up to twice the undisturbed solar wind density and a depletion is visible. This small enhancement indicates the presence of a Mach cone as described in other simulations (Bagdonat & Motschmann 2002b;Rubin et al. 2014). The orientation of the solar wind velocity remains the same. Only in the region of the depletion and close to the centre of the comet the solar wind velocity is slightly more oriented towards +z. cometary ion density drops significantly and the speed of the ions increases. In contrast to strongly active comets, which reveal a large anti-sunwards plasma tail, the simulation for the present scenario shows only a small anti-sunwards plasma tail.
In addition, at a few tens of kilometres the cometary ions are mainly moving in -z direction, i.e. in the direction of the solar wind convective electric field, with only a minor anti-sunwards motion occurring (Figs. 2e and f). This supports our separation of the interaction region into a +E and a −E hemisphere. Besides these patterns, which are familiar from previous simulation studies (cf. Bagdonat & Motschmann 2002b;Rubin et al. 2014), a new pattern is visible upstream and beneath the cometary nucleus in the region with z < 0. In the +E hemisphere, several thin and nearly parallel lines of equal ion density appear. Similar lines are present in the x-y plane (Figs. 2d). In the first instance we refer to these structures as fan-like structures as discussed later in further detail.
The configuration of the magnetic field on the y = 0-and z = 0 cross section is presented in Figs. 3a and b. The magnetic field, which is oriented towards the y direction in the undisturbed solar wind, is slightly bent around the comet. In addition, a small magnetic pile-up to field strengths up to 6 nT ≈ 3B IMF occurs upstream of the comet. Plasma structures, such as bow shocks or diamagnetic cavities are not present, as anticipated for the low-activity case studied here (cf. Koenders et al. 2013). Instead magnetic signatures of the Mach cone and the pick-up ion tail are visible (cf. Rubin et al. 2014).
The electric field (Figs. 3c and d), is orientated towards the -z direction in most parts of the simulation box. Thus, the field in these regions is controlled by the convective electric field of the solar wind: E con = −u i × B, where the mean ion velocity is about the velocity of the solar wind u sw . The electric field is strongly reduced only within about a 20 km radius to the nucleus, since the mean velocity is strongly slowed down by the high number of slow newborn ions and by ion-neutral collisions. In this small region, the simulation predicts a radial flow of the cometary ions away from the nucleus, which is deflected into the small plasma tail. Because of the limited resolution and the small spatial scale of this region, however, this study indicates the presence of two separated plasma regions at the weakly active comet.
The total electric current density (Figs. 3e and f) is determined via ∇ × B. On the z = 0 cross section, the current distribution exhibits fan-like structures, which are also visible in the cometary ion density distribution.
As our model comet is only weakly active the velocity-space configuration at (0, −200, −200) km, that is in the +E hemisphere ( Fig. 4) differs to a large extent from that observed at more active comets (cf. Neugebauer et al. 1989;Coates et al. 1993). The distribution is elongated in the z direction and no ring or shell distribution is apparent. The maximum velocity is about 80 km s −1 in the -z direction. These fast ions were ionised above the comet and, therefore, travelled a longer distance into the observation point and experienced a longer acceleration by the convective electric field. Ions with a very low speed are ionised closer to the observation point. The distribution reveals a kink at v z ≈ −40 km s −1 , where the ions have a much larger negative v y component than the other ions. Particle tracking shows that the ions in the kink originate from regions close to the nucleus. The velocity-space distribution of the solar wind protons is not shown; we merely note that the macro-particles representing the solar wind are located in a cloud that is centred around their bulk velocity of 394 km s −1 . The slowest solar wind particle still has a velocity of about 320 km s −1 . Thus, there are no protons moving in a sunwards direction.
As already mentioned, in general, the global structure of the interaction region in our simulations is similar to the structures presented in Bagdonat & Motschmann (2002b) and in Rubin et al. (2014). Our simulation reveals, however, the additional presence of a fan-like structure and a kink in the velocity space configuration of the cometary ions. In the low-activity stage, the convective electric field in the solar wind-comet interaction region is mostly driven by the solar wind, and it dominates the Hall-and the ∇p e -terms in Eq. (3). In addition and in contrast to the strongly active stage, the draping of the magnetic field is much weaker. This also holds for the pile-up of the magnetic field ahead of the comet. Consequently, most of the cometary ions are picked up into a cycloidal motion perpendicular to the magnetic field. The height of these cycloids is about h = 2r gyr,ci = 8 × 10 4 km. Hence, the cometary ions move mainly towards the -z direction in the simulation. Richter et al. (2015) reported the presence of a low-frequency wave activity in the environment of the comet 67P/Churyumov-Gerasimenko during its weakly active phase (Fig. 1). To compare observed and simulated wave properties, we select four different points within the simulation box, in which we consider time series and determine wave properties such as frequency and polarisation. These points are A = (0, 0, −50) km, B = (0, −200, −200) km, C = (0, +200, −200) km, and D = (0, 0, +200) km. Point D is located in the −E hemisphere, while the other points are located in the +E hemisphere. Point B is located in the negative y hemisphere, while C is located in the positive hemisphere. The temporal variations of relevant parameters are shown for points A and B in Figs. 5 and 6, respectively. At point A the magnetic field and cometary ion density data exhibit irregular low-frequency oscillations, while the solar wind density is almost not varying (Fig. 5). Prominent field oscillations in the low-frequency regime, however, are clearly visible at points B (Fig. 6) and C (simulated data not shown). The magnetic field is characterised by variations between 1 nT pp and 2 nT pp , amplitudes, which is comparable to the ambient magnetic field magnitude at the two selected points. In case of the cometary ion density, the variations are even larger than the quasi-stationary values. The solar wind density remains nearly constant as already seen at point A. At point D the variations are negligible or very small. The magnetic field oscillates with an amplitude of about 0.04 nT pp and the cometary ions with about 0.1 cm −3 . The solar wind density remains constant, which  Notes. The symbols B, E, j, n sw , n ci , u sw , and u ci denote the mean values of the magnetic field, electric field, current density, solar wind and cometary density, and the bulk velocities of the solar wind and cometary ions. By means of a minimum variance analysis of the magnetic field, the orientation of the normalised wave vectork was determined. The ratio of the eigenvalues is given by λ B,1,2,3 . λ E,1,2,3 denotes the ratio of the eigenvalues of a minimum variance analysis of the electric field andk E,max denotes the orientation of the maximum variance of the electric field. f gyr,e , f gyr,sw , and f gyr,ci denote the local gyrofrequency of solar wind protons and cometary ions.

Local wave analysis
is similar to observations at the other points. Mean values of the various quantities are listed in Table 2.
Results of a spectral analysis of the time series derived from the simulation at points A and B are presented in Figs. 7 and 8. At point A, which is very close to the nucleus, neither the magnetic nor the solar wind and cometary ion density exhibit any clear spectral peak indicating the absence of major plasma oscillations. The absence of coherent oscillations is also indicated by the low coherence between the magnetic field and plasma oscillations ( Fig. 7 bottom panel). Further away from the nucleus, at point B = (0, −200, −200) km, however, the magnetic field reveals a clear spectral peak at about 95 mHz. A similar spectral peak is also observed in the cometary ion density; no spectral enhancement is present in the solar wind proton density. The spectral peak of the simulated magnetic field oscillations corresponds to the spectral peak in the actual observations made at 67P/Churyumov-Gerasimenko .
For a coherence analysis the magnetic field and cometary ion velocity components are rotated into the magnetic field minimum variance system (k min ,k int ,k max ). The coherencies between the magnetic field components in the maximum and medium variance directions and the cometary ion density reach values of about 0.7 at the peak frequency. Because of ∇ · B = 0, a minimum variance analysis is a suitable tool to estimate wave propagation directions. At point A, which is rather close to the nucleus, the wave propagation direction is not well defined. The ratios of the eigenvalues of the magnetic field oscillations in the minimum, medium, and maximum variance direction are 1:1.9:2.2. At point B, however, the ratios are 1:10:18. Therefore the minimum variance analysis of the magnetic field data provides the well-defined wave propagation direction k/|k| = (0.8, 0.6, 0.0). This direction corresponds to the direction of the normal of the fan-like structures at the observation point B. Moreover, the eigenvalue ratios reveal that the magnetic field is elliptically polarised in accordance to Richter et al. (2016) Contrary to this, the electric field at point B is polarised linearly, as can be seen in Fig. 6, and the ratio of the eigenvalues of 1:2:27 result from a minimum variance analysis of the electric field variations. The angle between wave propagation and the electric field oscillation directions is out 8 • at point B. At point C similar results are obtained. All values are listed for completeness in Table 2.
Transformation of the magnetic field into a local mean magnetic field-aligned system furthermore allows us to determine polarisation properties of the waves. The hodograms are shown in Fig. 9 and indicate that at points B and C, that is at the two  corresponding points of the positive and negative y hemisphere, the waves are neither purely compressional nor purely transversal. The polarity of the magnetic field is right handed.
Our local wave analysis reveals the existence of wave activity at around 95 mHz, at frequencies far above the local cometary ion gyrofrequency. The waves are large-amplitude waves. Their magnetic polarisation is neither purely transverse nor purely compressional. Cross-coherence is observed for magnetic and cometary ion density oscillations. The solar wind protons contribute only margninally to the oscillations. Wave properties in the x-y-plane are almost symmetric with respect to the x-axis. Propagation direction is along the normal of the fan-like structures and only in the x-y-plane. With respect to wave frequency, amplitude, and magnetic polarisation, the simulated wave properties agree with those observed at 6P/Churyumov-Gerasimenko .

Wake structure
The global structure of the interaction region and local wave analysis suggests that the interaction is a superposition of two  types of structures: a stationary wake and a time-varying wake (the fan-like structures). In order to prepare for the analysis of the fan-like structure in the next section, this sections aims to separate the fan-like structure and the signatures in the wake. For this purpose Fig. 10 shows the z component of the magnetic field. Panels a, b, and c give details of the x-y-plane in the −E hemisphere, the z = 0 km-plane, and the +E hemisphere, respectively. Panel d provides information on the structure in the x = 0 cross section.
In the z component of the magnetic field the fan-like structure, i.e. the thin, nearly parallel lines, can easily be identified. These lines, which are present on the dayside, are merely observed in the +E hemisphere, but not in the −E hemisphere. In addition to the fan-like structure, a wake structure on the nightside is visible. This structure consists of wedge-shaped regions in the z = const. cross sections of the B z component. In these wedges, the sign of B z is opposed to the wedge on the opposite y hemisphere. This pattern is only present on the x-y-plane, thus, the wake is aligned to the magnetic field.
Since we are able to control various processes in the simulations, we are also able to isolate the wake structure from the fan-like structure. To do so, a simulation, hereafter called Simulation B, in which the Hall-and the ∇p e terms in the electric field calculation are switched off, was conducted in a first instance. This means that the electric field in this simulation is only determined by the convective electric field. All other parameters are equal to the default simulation, called Simulation A. The final state of the cometary ion density and the cometary ion bulk velocity is presented in Figs. 11a and b. In contrast to Simulation A, the cometary ions are much more focussed in the +E hemisphere. Moreover, hardly any anti-sunwards motion of  In the next step, we conducted a simulation in which the Hall term and the ∇p e terms in the electric field calculation are reactivated. In exchange the treatment of the cometary ions in this simulation is manipulated. This simulation uses the final configuration of the cometary ions in Simulation B as the initial configuration, i.e. their density and velocity distribution, but now the further movement and acceleration of the cometary ions are switched off, i.e. the ions remain at their position and their velocity is not updated. Cometary ions are still considered in the calculation of the moments. Thus, the cometary ions in the simulation box act as a stationary current source that interacts with the impinging solar wind. In contrast to this very special current configuration, a simple line current in -z direction could also be used that would lead to similar results (cf. Thompson et al. 1996). However, we used this special current configuration to enhance the comparability of the resulting structures. As one can see from the resulting magnetic field configuration, a wing-type wake forms behind the comet (Fig. 11d). This wake structure is similar to the wake of Simulation A (Fig. 10c). Consequently, we can deduce that the wake is caused by the current caused by the cometary ions, moving in -z direction. Since this is about constant in time, the wing-typed wake is a quasi-stationary structure in the environment. In this simulation no field fluctuations in the magnetic field can be observed except for the initial field changes at the start of the simulation. We state that under conditions inhibiting the excitation of waves, no fan-like structures are observed, and merely a classical interaction picture exists.

Global wave analysis
The discussion of the global structures in the plasma interaction region and the observed low-frequency wave activity supports the assumption that the fan-like structure is linked to the wave activity. This section aims to analyse the properties of these fanlike structures.
The fan-like structures are not present in the entire simulation box (Fig. 10). In the z = 0-cross section, the fan-like structures start to be visible (Fig. 10b), but these structures are most pronounced in panel c. Since they are not visible on the cross section 200 km above the comet, the fan-like structures are only present in the direct vicinity of the nucleus and the +E hemisphere, which agrees with the observations carried out in the local wave analysis of the low-frequency wave.
The opening angle of the fan-like structures, the angle with respect to the x-axis, decreases from 70 • to about 45 • closer to the nucleus. This exceeds the opening angle of the Mach cone, of about 8 • in the far field, and that of the wing-typed wake. Similar to the Mach cone, which is a standing wave caused by the fast magnetosonic wave, we anticipate that the fan-like structure is an interference pattern of a wave. Thus, the opening angle of the fan-like structure is given by the local flow speed and the phase velocity of the wave. Since these parameters change rapidly in the direct environment of the comet, the opening angle changes as well. In addition, the distances between the maxima or minima of the fan structure change from small values upstream to larger distances downstream. At the point (0, −200, −200) km the distance is 55(21) km along k direction. The figure shows that the phase between the two y hemispheres seems to be shifted. On the y = 0 cross section, the enhancements are nearly parallel to the z-axis, but starting at about z = 100 km.
The temporal evolution of the fan-like structure is presented in a movie. In the movie the z component of the magnetic field on the z = −200 km cross section is shown for a time interval of about 240 s. It is clearly visible that the fan-like structure is moving downstream and that the size, i.e. the distance between two nearly parallel lines, increases towards +x direction. In addition, it is visible that the fan-like structures split at several places during the sequence. Similar effects can be observed in Fig. 12 A similar increase can be observed in Fig. 12b, which is indicted by a shorter interval. We estimate the local velocity of the fan-like structures at that point, using a graphical approach and find a mean velocity of 6.8(16) km s −1 . Using the distance between the maxima of about 55 km, we get a frequency of f = v/λ = 123 mHz, which is close to the frequency observed at point B in the local wave analysis. Thus, the fan-like structures are linked to the low-frequency wave activity discussed in the local wave analysis. As one can see from that figure as well, however, some of the fan-like structures split during the simulation and also small variations in the amplitude are visible. It can be deduced from Fig. 12c, which shows the B z component along the y-axis, that the signatures in the magnetic field are more chaotic within a distance of 200 km to the y = 0-line than at larger distances.

Parameter studies
The simulations presented here are rather time consuming. Therefore, a full parameter study is not possible at the moment. This mainly effects the determination of the wave frequency on parameters such as interplanetary magnetic field strength or cometary activity. Such a study would require us to run the simulations for an extended time, much beyond the stage where quasi-stationary situations are reached. Parameter studies, however, are possible with respect to the spatial properties. Table 3 lists all simulations performed together with the changed parameters. The simulation described so far is used as the default simulation and labelled as Simulation A. Simulations B and C have already been discussed in Sect. 5. Simulations D and E are similar to Simulation A, but use different directions of the interplanetary magnetic field. The Parker angle in simulation D is θ = 45 • (Fig. 13), while it is θ = 135 • in Simulation E. The latter result is not shown since it is similar to the result of Simulation D. The structures are mirrored at the y = 0 line and the sign of the z component of the magnetic field is changed. All other parameters of these simulations are equal to those of Simulation A. As a result of the changed IMF orientation the structure of the wing-typed wake also changed. In comparison to simulation A with θ = 90 • the fan-like structure is no longer symmetric to the Sun-comet line.
In analogy to a bow shock structure, a quasi-parallel and quasi-perpendicular regime can be identified, both with respect to the orientation of the background magnetic field and the Table 3. Overview of different simulation set-ups used in this study.

Simulation Properties
A default simulation with Parker angle 90 • B as A but without Hall-and ∇p e -terms in field calculation C as A but the movement of the cometary ions and their acceleration is suppressed. The cometary ions are considered in the calculation of the moments. As initial state of this simulation the final state of Simulation B is used. D as A but Parker angle is 45 • E as A but Parker angle is 135 • F 67P/C-G at 2.3 AU normal of the fan-like structures. In the quasi-perpendicular regime, the phase fronts are compressed, i.e. the wavelength of the corresponding waves is reduced. On the quasi-parallel side, a larger distance between the fan-like structures can be observed, which corresponds to larger wavelengths. In addition, in comparison to the default Simulation A, phase fronts appear early, i.e. at smaller x-positions on the quasi-parallel side. Moreover, in comparison to the fan-like structure in Simulation A, their wave amplitudes are reduced in Simulations B and C. This is probably due to the reduced convective electric field as E = −u i × B. Finally, Simulation F covers the plasma environment at comet 67P/Churyumov-Gerasimenko in February 2015, with solar wind conditions corresponding to the standard conditions at 2.3 AU (Parker 1958;Hansen et al. 2007). In contrast to Simulation A Simulation F uses a Parker angle of 66 • . The density in the neutral gas atmosphere is given by Haser (1957), as in Simulation A. The total neutral gas production rate is set to Q = 5 × 10 26 s −1 . Figure 14a shows the solar wind density and solar wind velocity in the y = 0 cross section. In this scenario the plasma around the comet is much more disturbed than in the scenario of Simulation A. Here a strong density enhancement on the +z hemisphere is present. In addition, a stronger deflection of the solar wind protons is visible. In contrast to Simulation A, the cometary ions (Fig. 14b) reveal a significant velocity component in the x direction, besides the velocity in the -z direction. The fan-like structure is also present in this scenario (Figs. 14c and d). The opening angle is smaller than in Simulation A, however, and the structures are not longer parallel to the z-axis. Instead, it seems that the fan structure is aligned parallel to the velocity of the cometary ion velocity, which is also true for Simulation A.

Discussion
The global structure of our simulated solar wind-comet interaction exhibits a threefold wave structure that is caused by the presence of cometary ions at first instance. In the weakly active stage, most of the newborn cometary ions are accelerated by the solar wind convective electric field E = −u i × B. The acceleration by the magnetic force, qu × B, at the beginning is negligible. Therefore the cometary ions first move into −z direction, i.e. perpendicular to the solar wind flow and to the magnetic field. For comparison, a cometary ion test particle injected into undisturbed solar wind moves about 500 km in −z direction and gains a speed of about 62 km s −1 within 15.9 s, which is close to the actual speed of the cometary ions in the simulation at z = −500 km. The motion in the x direction is still small at this distance, about 5 km s −1 . As shown by Rubin et al. (2014) the height of the cycloid is about 72 000 km, which exceeds the size of our simulation box. These results match results known from previous simulations of weakly active comets (cf. Bagdonat & Motschmann 2002b;Rubin et al. 2014). This has also been reported by Nilsson et al. (2015a,b) and Behar et al. (2016), who used actual RPC-ICA measurements from the coma of 67P/Churyumov-Gerasimenko. Furthermore, our simulation predicts an elongated distribution of the cometary ions in the velocity space, which is clearly separated from the solar wind protons. This is in agreement to the observations made by the RPC-ICA and the RPC-IES instrument (cf. Nilsson et al. 2015a;Broiles et al. 2015). However, the velocity space configuration found differs to a large extend from the observations made at more active comets, where ring-beam or non-gyrotropic ion distributions dominate (Neugebauer et al. 1989;Coates et al. 1993;Motschmann & Glassmeier 1993).
Following the simulation and observations, a separation of the plasma environment into a +E hemisphere and an −E hemisphere is favourable. While the newborn cometary ions move into the +E hemisphere parallel to the solar wind electric field, the solar wind protons are deflected into the −E hemisphere. This is why the slow cometary ions reduce the mean ion velocity and the convective electric field in the interaction region. Consequently, the magnetic force is no longer balanced by the electric force and the solar wind protons are deflected in the −E hemisphere. Broiles et al. (2015) and Behar et al. (2016) also came to the same conclusions using actual observations. In the simulation the deflected protons trigger the first interaction signature, the Mach cone. The associated wave handling this proton velocity perturbation is a standing fast magnetosonic wave. We denote this wave as the first wave structure. Its opening angle in the far field is given by the magnetosonic speed and the flow speed, which is about 8 • in the current scenario. More details on the Mach cone are given by Bagdonat & Motschmann (2002b) and Rubin et al. (2014).
A second interaction signature, the wing-typed wave, which is described in more detail in Sect. 5, is caused by the cross-field current density of the initial pick up of cometary ions. In simulations of the cometary interaction this type of wake has not been discussed yet. In Simulation A (Figs. 10b and c) this wake is superimposed by the fan-like structure, which is discussed later. This fan-like structure was removed in Simulation C, where the motion and acceleration of the cometary ions are suppressed but the cometary ions are considered in the calculation of the moments. Thus, the cometary ions contribute to a current, which is considered in the field calculation. The current is large in the vicinity of the nucleus and in the +E hemisphere. In this region the wake signature is present, while it is much weaker in the −E hemisphere, where the cometary ion density is low and the current density is very low. The wake signature in the magnetic field is not accompanied by a density variation and, as demonstrated by Simulation D and E, the wake structure changes its orientation with changes in the magnetic field orientation. Thus, the wing-typed wave propagates close to the magnetic field and the solar wind transports the wave downstream.
A similar structure has been predicted for magnetised asteroids (Baumgärtel et al. 1994(Baumgärtel et al. , 1997Gurnett 1995b;Omidi 2002) or a constant current source (Thompson et al. 1996) in the solar wind. These authors found that a Whistler wave is triggered by the small obstacle instead of a bow shock. Since the Whistler waves mainly propagate close to the magnetic field (about 19 • ) they are confined to a cone around the magnetic field line passing the obstacle (Gurnett 1995a). Because the wave is superimposed by the solar wind, it is transported downstream by the solar wind and creates a wedge-shaped region that can be accessed by the waves. This is the case at least for a perpendicular magnetic field configuration. Similar characteristics can be observed at the wake structure in our simulation. Since the current system is extended into the +E hemisphere, however, the wake is also extended towards this hemisphere. Nevertheless, we can conclude that a weakly active comet triggers a Whistler wing structure.
As a third type of interaction signature we consider the low-frequency wave activity at a frequency of about 95 mHz in the frame of a resting observer. The waves are observed in the +E hemisphere, while the waves are barely present in the −E hemisphere. As reported above, the simulation suggests that the frequency is nearly constant in the interaction region although the local plasma properties change drastically in the simulation box. The polarisation of the simulated magnetic field fluctuations is neither transverse nor purely compressional, but of mixed nature. This is also reported by Richter et al. (2015Richter et al. ( , 2016 for the waves observed at 67P/Churyumov-Gerasimenko. This third wave type creates a complex interference pattern, the fan-like structure. This is triggered by an active process in the plasma, rather than a simple obstacle, as the steady deflection of the solar wind or the cometary ion pick-up current. As the observed wave frequency is about 50 times higher than the local cometary ion gyrofrequency, this wave activity differs to a large extend from waves observed at various other comets (cf. Smith et al. 1986;Neubauer et al. 1986;Tsurutani et al. 1989;Glassmeier et al. 1989;Motschmann & Glassmeier 1993;Neubauer et al. 1993;Richter et al. 2011). At these comets ion cyclotron waves with frequencies (in the spacecraft or cometary frame of reference) close to the cometary ion gyro-frequency were identified. These waves are excited by an ion ring-beam instability (cf. Wu & Davidson 1972;Volwerk et al. 2013b) and grow on spatial scales of several cometary ion gyro-radii, which are scales much larger than the scale of the simulation box. Richter et al. (2015) suggest a possible wave source mechanism: an excitation by a cross-field electric current. This type of excitation mechanism has been studied by Chang et al. (1990) and Sauer et al. (1998) and leads to large growth rates when the current is perpendicular to the ambient magnetic field and becomes unstable. The current is initially caused by the acceleration of the cometary ions perpendicular to the interplanetary magnetic field. Our simulations confirm that the cometary ions move approximately perpendicular to the magnetic field (Fig. 15) because the dynamics of the cometary ions are dominated by the solar wind convective electric field. Thus, these cometary ions contribute to a current that is perpendicular to the magnetic field. We can also report that in case the current system is set to be constant, the fan-like structure and subsequent low-frequency waves vanish, as demonstrated in Simulation C. In case the motion of the cometary ions is not suppressed, the cometary ion current can vary and may become unstable. Thus, the low-frequency wave is linked to the variations in the cometary ion motion and the subsequent cometary ion current. In their first approach Richter et al. (2015) assess the variations in the unstable cometary ion current to be of the order of a few nA m −2 . Using magnetic field variations of the order of 1 nT and Ampere's law, a transverse wavelength of 500 km was estimated. In the simulation we see that the current driven by the cometary ion j ci = en ci u ci is not constant and varies (Fig. 15). For example, at point (0, −200, −200) km the cometary ion current varies by approximately 2.5 nA m −2 and the magnetic field fluctuates by about 0.5 nT. This leads to a transverse wavelength of about 1250 km. By using the variations at point (0, 0, −50) km, we gain δ j ci ∼ 10 nA m −2 , δB ∼ 1 nT, and a wavelength of 250 km. Both estimated wavelengths are similar to the wavelength estimated by Richter et al. (2015) but differ from the distance between the maxima in the fan-like structure of about 55 km. Since this fan-like structure is caused by an interference with the impinging solar wind, however, the wavelength cannot be compared to the distances of the fan-like structures.
A spatially distributed source, such as a region of an unstable cometary ion current, likely explains the complex shape of the fan-like structure. Contrarily, a localised wave source that continuously emits a wave in a supersonic flow leads to a single standing wave. A sketch of this interference is shown in Fig. 16a. In contrast to the classical picture we use a wave of a more complex character, which leads to an enhancement on one side and a reduction on the other. An example for this type is the Mach cone, which is triggered continuously by the deflection of the solar wind. As discussed, it leads to an enhancement of the magnetic field and density on one side and a reduction of these quantities on the other. In case waves are emitted from several locations, the resulting interference pattern should consist of several standing waves (Fig. 16b). If the wave excitation is in opposite directions, our simple superposition scheme creates a pattern that is similar to the fan-like structure we observe in the simulation. In a) a wave is triggered continuously and transported with the solar wind flow. Owing to superposition a standing wave is created. Instead of a single wave source several sources are used in b), which leads to the formation of several standing waves. In case c) the sources triggers the waves in opposite directions, which creates a interference pattern that is similar to the fan-like structures.
This simple scheme is only applicable if the resulting perturbations are small and the flow is not strongly modified, which is true in the case of the fan-like structure. Such an extended wave source, which emits waves in opposite directions, is the unstable current, which is present in the larger vicinity of the nucleus. Consequently, it should be able to create an interference pattern similar to the fan-like structure. Since this current is embedded in the flow of the cometary ions, it moves into the +E hemisphere. Consequently, on each z = const. cross section a fan-like structure is triggered, which leads to the enhancements nearly parallel to the z-axis (Fig. 10d).
The first approach to the wave excitation and the simple explanation of the fan-like structure need to be investigated in more theoretical detail. A more detailed linear instability analysis of the cross-field current situation is presented by Meier et al. (2016).

Conclusions
This work studies the global state of the plasma environment of comet 67P/Churyumov-Gerasimenko during its weakly active phase by the use of hybrid plasma simulations. In this activity stage the classical structures known from comets, such as bow shocks or diamagnetic cavities, are not present. Instead, our simulations reveal that the newborn cometary ions create three different types of wave structures, which dominate the plasma environment.
The presence of slow newborn cometary ions leads to a deflection of the solar wind protons. This perturbation triggers a standing fast magnetosonic wave, a Mach cone, which was also observed by investigators in previous studies (cf. Bagdonat & Motschmann 2002a;Rubin et al. 2014). We denote this wave as the first wave structure.
In the simulation the newborn cometary ions are accelerated by the solar wind convective electric field perpendicular to the solar wind flow and the interplanetary magnetic field. This was also found by Nilsson et al. (2015a,b), Broiles et al. (2015), Behar et al. (2016) using in situ observations. The moving cometary ions in the environment of the comet also represents a cross-field current density, which causes the second type of wave, a Whistler wing.
As the third type of wave, we consider low-frequency wave activity. In the simulation this wave activity peaks at a frequency of 95 mHz and its polarisation is neither purely transverse nor purely compressional, but of mixed nature. Similar observations were made by Richter et al. (2015Richter et al. ( , 2016 who found a nearly continuous wave activity between 20 and 100 mHz. In the simulations this low-frequency wave activity causes a complex interference pattern, i.e. a fan-like structure, which can be found in the +E hemisphere. A substantial cross-field cometary current density can only be found in this hemisphere. As suggested by Richter et al. (2015), this third type of wave is triggered by an instability of the cometary ion cross-field current. We demonstrated that indeed a constant cross-field current density is not able to produce the fan-like structures. Our simulation study shows that the Mach cone and the Whistler wing are quasistationary structures, while the third type of wave is caused by an active process in the plasma flow. Furthermore, our simulations reveal that the plasma environment is dominated by the initial pick up of cometary ions and the superposition of these three wave structures.