<figure><img src="image_1.png" /></figure>

CrossMark

# A Fully Kinetic Perspective of Electron Acceleration around a Weakly Outgassing Comet

Andrey Divin<fnref n="1" /> ID, Jan Deca<fnref n="2" />,<fnref n="3" /> ID, Anders Eriksson<fnref n="4" />, Pierre Henri<fnref n="5" />,<fnref n="6" />, Giovanni Lapenta<fnref n="7" /> ID, Vyacheslav Olshevsky<fnref n="8" /> ID, and Stefano Markidis<fnref n="8" />

Received 2019 October 13; revised 2019 December 17; accepted 2019 December 30; published 2020 January 29

## Abstract

The cometary mission *Rosetta* has shown the presence of higher-than-expected suprathermal electron fluxes. In this study, using 3D fully kinetic electromagnetic simulations of the interaction of the solar wind with a comet, we constrain the kinetic mechanism that is responsible for the bulk electron energization that creates the suprathermal distribution from the warm background of solar wind electrons. We identify and characterize the magnetic field-aligned ambipolar electric field that ensures quasi-neutrality and traps warm electrons. Solar wind electrons are accelerated to energies as high as 50–70 eV close to the comet nucleus without the need for wave–particle or turbulent heating mechanisms. We find that the accelerating potential controls the parallel electron temperature, total density, and (to a lesser degree) the perpendicular electron temperature and the magnetic field magnitude. Our self-consistent approach enables us to better understand the underlying plasma processes that govern the near-comet plasma environment.

*Unified Astronomy Thesaurus concepts: Space plasmas* (1544); *Astronomical simulations* (1857); *Comets* (280)

## 1. Introduction

Cometary nuclei are small solar system bodies composed of dust, rocks, volatiles, and ices. As a cometary nucleus travels closer to the Sun, its surface warms up and the environment gets embedded with dust and gas, producing a large-scale comet. A neutral cold atmosphere is ionized by sunlight or via collisions with solar wind particles (Cravens & Gombosi 2004) or by suprathermal electrons (Galand et al. 2016) as the gas escapes slowly to space.

Over the course of the past several decades studies have focused mostly on quite active comets (Cravens & Gombosi 2004), such as comet Giacobini–Zinner (observed by the *International Cometary Explorer* at ~1.03 au with an estimated gas production rate of $2–4 \times 10^{28}$ s⁻<fnref n="1" />; Mendis et al. 1985), or comet Halley (observed by the VEGA 1 spacecraft at 0.79 au with a gas production rate of $1.3 \times 10^{30}$ s⁻<fnref n="1" />; Rème 1991). Most recently the ESA/*Rosetta* spacecraft accompanied comet 67P/Churyumov–Gerasimenko (hereafter 67P), providing a unique and unprecedented perspective into the evolution of a comet. *Rosetta* escorted 67P for more than two years (Glassmeier et al. 2007; Taylor et al. 2017). Compared to the targets of previous cometary missions, 67P is a weakly active object. At a heliocentric distance of 3.0 au a gas production rate of $10^{26}$ s⁻<fnref n="1" /> was found (Bieler et al. 2015), indicating a solar wind–coma interaction in a highly kinetic regime (Rubin et al. 2014).

One of the surprising findings of the *Rosetta* mission is the presence of suprathermal electrons in the close cometary plasma environment with energies up to about 100 eV. The population was present already during the weakly outgassing phases of 67P’s orbit around the Sun (Clark et al. 2015). Understanding the suprathermal electron population is

important, since increased fluxes of the latter have been shown to strongly affect also the cometary ionosphere via electron-impact ionization (Galand et al. 2016), charge exchange (Wedlund et al. 2017; Heritier et al. 2018), and is thought to affect dust grain charging processes (Gombosi et al. 2015).

A suggested mechanism to generate a suprathermal electron population is by interaction with lower hybrid waves (Broiles et al. 2016; André et al. 2017). A complementary mechanism discussed by Madanian et al. (2016) proposes that an ambipolar electric field established by the electron pressure gradient close to the nucleus provides sufficient nonturbulent acceleration to explain the observed electron velocity distributions. The latter was confirmed by Deca et al. (2017) using 3D fully kinetic particle-in-cell simulations, showing as well that the suprathermal population is primarily of solar wind origin. In this paper we investigate in detail the ambipolar acceleration mechanism that produces the observed suprathermal electron populations for comet 67P during its weakly outgassing phases.

## 2. Simulation Setup

This work makes use of the fully kinetic and electromagnetic particle-in-cell code iPIC3D (Markidis et al. 2010). The code solves the Vlasov–Maxwell coupled system of equations for multiple ion and electron species, and implements the implicit moment method that ensures numerical stability without the need for the computational grid and time step to resolve the Debye length and plasma frequencies (Brackbill & Forlund 1982; Lapenta et al. 2006). iPIC3D has been successfully applied to study, e.g., collisionless magnetic reconnection (Lapenta et al. 2010; Divin et al. 2016), magnetic null points (Olshevsky et al. 2015, 2016), and the solar wind

<footnote n="1">¹ St, Petersburg State University, Ulianovskaya, 1, 198504, St. Petersburg, Russia; andrey.div@gmail.com, a.divin@spbu.ru</footnote>

<footnote n="2">² Laboratory for Atmospheric and Space Physics (LASP), University of Colorado Boulder, Boulder, CO 80303, USA</footnote>

<footnote n="3">³ Institute for Modeling Plasma, Atmospheres and Cosmic Dust, NASA/SSERVI, Moffett Field, CA 94035, USA</footnote>

<footnote n="4">⁴ Swedish Institute of Space Physics, Box 537 SE-751 21, Uppsala, Sweden</footnote>

<footnote n="5">⁵ LPC2E, CNRS, 3 Avenue de la Recherche Scientifique, F-45071, Orléans, France</footnote>

<footnote n="6">⁶ Laboratoire Lagrange, OCA, UCA, CNRS, Nice, France</footnote>

<footnote n="7">⁷ Centre for mathematical Plasma Astrophysics (CmPA), Department of Mathematics, KU Leuven, Celestijnenlaan 200B, bus 2400 B-3001 Leuven, Belgium</footnote>

<footnote n="8">⁸ KTH Royal Institute of Technology, SE-100 44, Stockholm, Sweden</footnote>
<table><caption>**Table 1**
Simulation Parameters and Setup</caption><thead><tr><th></th><th></th><th></th><th></th></tr></thead><tbody><tr><td>T<sub>e,sw</sub></td><td>10 eV</td><td>n<sub>e,sw</sub>, n<sub>i,sw</sub></td><td>1 cm<sup>-3</sup></td></tr><tr><td>T<sub>i,sw</sub></td><td>7 eV</td><td>ρ<sub>e,sw</sub></td><td>8 km</td></tr><tr><td>v<sub>sw</sub></td><td>400 km s<sup>-1</sup></td><td>B<sub>sw</sub></td><td>6 nT</td></tr><tr><td>m<sub>i,sw</sub>/m<sub>e</sub></td><td>100</td><td>m<sub>i,c</sub>/m<sub>i,sw</sub></td><td>20</td></tr><tr><td>T<sub>e,c</sub></td><td>10 eV</td><td>Q</td><td>10<sup>25</sup> s<sup>-1</sup></td></tr><tr><td>T<sub>i,c</sub></td><td>0.026 eV</td><td>d<sub>i</sub></td><td>220 km</td></tr><tr><td colspan="2">Domain Size</td><td colspan="2">3200 × 2200 × 2200 km<sup>3</sup></td></tr><tr><td colspan="2">Resolution</td><td colspan="2">10 × 10 × 10 km<sup>3</sup></td></tr></tbody></table>

**Note.** The subscripts “e,” “i” represent electron and ion quantities, respectively, and “sw,” “c” represent solar wind and cometary parameters, respectively. Here, T, n, v, B, Q, dᵢ, ρₑ are temperature, density, bulk flow velocity, magnetic field, outgassing rate, solar wind ion inertial length, and electron gyroradius, respectively. Artificial solar wind ion (proton) to electron mass ratio is used, mᵢ,sw/mₑ = 100. Cometary ion to solar wind proton mass ratio is mᵢ,c/mᵢ,sw = 20.

interaction with lunar magnetic anomalies (Deca et al. 2018) and comets (Deca et al. 2017, 2019; Sishtla et al. 2019).

In this work, we use a similar setup as Deca et al. (2017), but assume a slightly lower outgassing rate that represents comet 67P/Churyumov-Gerasimenko at a heliocentric distance of approximately 4–4.5 au, a distance similar to when Rosetta first encountered the plasma signatures of the comet. We define the X-axis of the simulation domain along the solar wind flow. The interplanetary magnetic field points along Y and, hence, the solar wind convective electric field upstream of the comet is directed along the −Z direction. The comet location is defined at the origin of the computational domain. The mass loading process is modeled by injecting a Maxwellian distribution of cold cometary water ions and warm cometary electrons following an analytical profile that models the ionization of an expanding neutral gas with a $1/r^2$ radial density profile centered at the comet location. We assume a gas production rate of $Q = 10^{25}$ s$^{-1}$ (Deca et al. 2019). In addition, we lower the cometary ion temperature to $T_{i,c} = 300$ K to better reflect the observations (Nilsson et al. 2015). With a resolution of 10 km we do not resolve the cometary nucleus (~3 km) and therefore exclude it here to avoid numerical subgrid effects. The solar wind ions and electrons are injected into the computational domain with a Maxwellian distribution following the procedure outlined in Deca et al. (2015).

In our model, solar wind and cometary electrons are injected with equal temperatures (=10 eV). We denote electrons with energies ≤12.5 eV as “warm.” Electrons with energies between 12.5 and 25 eV are labeled “hot.” Those with energies exceeding 25 eV are marked as part of the “suprathelial” electron distribution. Note, because we do not include collisions in our model, the observed population of collision-cooled electrons (typically having energies as low as 0.1 eV; Eriksson et al. 2017; Engelhardt et al. 2018) is not present. An overview of all simulation parameters is given in Table 1.

### 3. Large-scale Interaction

The global morphology of the solar wind–comet interaction is presented in Figure 1. Using color-coded magnetic field lines we highlight the 3D structure of the magnetic pile-up region. The maximum field magnitude in this region reaches up to 2.6 times the interplanetary magnetic field magnitude, $B_{sw}$. The

magnetic field drapes around and passes through the region of the highest density near the nucleus (Koenders et al. 2016).

Cometary ions are accelerated (gray lines in Figure 1) predominantly in the direction of the solar wind convective electric field, $E_{sw}$ (corresponding to the −Z direction in Deca et al. 2017, 2019). The red-shaded isosurfaces indicate the cometary ion number density that is equal to the far upstream solar wind density ($n_{i,c} = 1$ cm$^{-3}$ contour). The density plume is located in the $x > 0, z < 0$ quadrant. The solar wind convective electric field accelerates the cometary ions up to 0.1–0.4 times the upstream solar wind speed, $v_{sw}$, before the particles reach the outflow boundary of the computational domain. Extending the latter would eventually result in pick-up acceleration up to $v_{sw}$. With the parameters considered in this simulation, the cometary ion density peaks at $14 n_{i,sw}$.

The blue isosurface in Figure 1 shows the 3D cometary electron profile for the $n_{e,c} = 0.2$ cm$^{-3}$ contour. Note that the cometary electrons are not coupled to the cometary ions and leave the source region along the magnetic field lines that intersect the population in the vicinity of the nucleus. A spatial separation of the particles of cometary origin can be seen in the density profiles of $n_{e,c}$ (Figure 2(a)) and $n_{i,c}$ (Figure 2(b)) in the $y = 0$ plane. The convective electric field carried by the solar wind is nonnegligible close to nucleus (Deca et al. 2019) and picks up the cometary electrons on a spatial scale comparable to the electron Larmor radius. After magnetization, they carry a velocity composed of an $E \times B$ drift component and a component parallel to their respective magnetic field line. Such a localized cometary electron current explains the cometary electron density profile in Figure 1, i.e., the two channels of cometary electrons emanating from the high-density part of the coma that move along the magnetic field lines in the solar wind frame. With the parameters considered in this simulation, the total electron density peaks at $14$ cm$^{-3}$ near the nucleus, combining the cometary (9 cm$^{-3}$) and solar wind (5 cm$^{-3}$) electron densities.

A large ambipolar electric field is established that is associated with the local electron density gradient. It temporarily traps the solar wind electron population (Madanian et al. 2016; Deca et al. 2017, 2019; Sishtla et al. 2019) that is visualized by blue “⊙” markers in Figure 1. Since cometary ions are much more inert (visualized by the magenta “⊕” markers in Figure 1), they need a considerably longer time to leave the region where they were ionized. The high solar wind electron density near the comet compensates for the charge imbalance (Figures 2(b), (e)), briefly trapping fast-moving electrons and confirming the coupled four-fluid scenario (Deca et al. 2017) also for a lower outgassing rate. Unlike the electrons, the solar wind ions are deflected in the +Z direction (Figure 2(f)), forming a weakly overdense region in the $x > 0, z > 0$ quadrant (Behar et al. 2018) and the solar wind ion cavity in the $x > 0, z < 0$ quadrant.

### 4. Electron Pitch-angle Distributions (PADs)

The presence of the ambipolar electric field parallel to the magnetic field provides favorable conditions to energize electrons above the solar wind temperature (Madanian et al. 2016). In order to disentangle the electron acceleration mechanisms, we plot an overview of the electron PADs along the Sun–comet line (X-direction) and through the center of the computation domain (Figure 3). The PADs are constructed by combining the electron velocity distribution in 12 uniformly
<figure><img src="image_2.png" /><figcaption>**Figure 1.** 3D overview of the simulation. Magnetic field lines are shown in color; cometary ion streamlines are shown in gray. A red isosurface displays the cometary ion density for $n_{i,c} = 1n_{i,sw} = 1 \text{ cm}^{-3}$. The blue isosurface shows the cometary electron density $n_{e,c} = 0.2n_{i,sw} = 0.2 \text{ cm}^{-3}$. The yellow surfaces (right panel) display the region of enhanced magnetic field. The red, green, and blue arrows in the right panel display the direction of $v_{sw}$, $B_{sw}$, and $E_{sw}$, respectively. Magenta “⊕” markers visualize electric charges created by the cometary ions, which are balanced the solar wind electrons (blue “⊙” markers).</figcaption></figure>

distributed pitch-angle bins, where the pitch angle $\theta$ is defined as the angle between the magnetic field and the velocity of the electron. The resulting distribution is then split in three energy channels and plotted separately for electrons of solar wind and cometary origin as follows: warm electrons (energies below 12.5 eV; Figure 3(a)), hot electrons (energies between 12.5 and 25 eV; Figure 3(b)), and suprathermal electrons (energies above 25 eV; Figure 3(c)). Note that the populations modeled here result from an initial Maxwellian velocity distribution; hence, we do not model and therefore do not have access to the halo and Strahl components of the solar wind. Observations have shown that the faint power-law tail may have a nonnegligible influence on the measured suprathermal electron fluxes (Broiles et al. 2016; Myllys et al. 2019). These contributions, however, cannot explain the total flux of suprathermal electrons, making the mechanism proposed here necessary to fully explain and understand the measurements by the Rosetta plasma instruments.

The distribution of the warm solar wind electron component is spatially uniform, except for a weak parallel enhancement near ~0° and ~180° close to the nucleus. We attribute this enhancement to a population of warm solar wind electrons that enter close to the nucleus and thus experience only minor acceleration by the ambipolar electric field. The hot and suprathermal electron populations close to the comet are generally of solar wind origin, as seen also by Deca et al. (2017). A strong day-night asymmetry exists with more energetic particles on the night side (Figure 3(d)).

We identify three regions in panels (a)–(d) of Figure 3: region A from 50 to 250 km upstream of the nucleus, region B surrounding the nucleus, and region C from 50 to 500 km downstream of the nucleus. We plot the omnidirectional electron phase space density (PSD) as a function of energy for the three regions in Figure 3(e) (at the point $(x, y, z) = (-100,0,0)$ km in region A), Figure 3(f) (at the point $(x, y, z) = (0,0,0)$ km in region B), and Figure 3(g) (at the point $(x, y, z) = (200,0,0)$ km in region C). Figures 3(h)–(j) compute the electron PSDs using the parallel ($F_\|$) and perpendicular ($F_\perp$)

component of the electron flux separately as illustrated in panel (k).

A different flavor of electron anisotropy is present in each of the identified regions. Region A corresponds to the pile-up region, where the magnetic field magnitude reaches its maximum. The solar wind population is compressed due to the magnetic flux pile-up there, causing betatron acceleration, which enhances the $\theta \sim 90^\circ$ electron flux (Figures 3(b) and 2(h)). Above 20 eV the perpendicular solar wind flux surpasses the parallel component (Figure 3(h)). The cometary electron flux is noisy and negligible there (Figure 3(e)).

In region B we observe the highest electron densities. The electron PSD is non-Maxwellian (Figure 3(f)) as it contains a main peak formed by the cometary electrons at energies below 10 eV and a peak formed by accelerated solar wind electrons with energies up to 50–70 eV. A knee is formed as a result of the superposition at about 15 eV.

Focusing on Figure 3(i), the warm cometary electron background shows a perpendicular flux greater than the parallel one at energies below 15 eV. For energies above 15 eV, however, the perpendicular flux drops sharply. The solar wind electrons on the other hand peak in the energy range of 20–50 eV. Here the parallel component dominates, as the perpendicular solar wind component drops above 20 eV. Thus, we conclude that the flux of energetic electrons seen near the nucleus is due to the parallel solar wind component. Note that for energies below 20 eV the solar wind is isotropic and behaves almost Maxwellian.

Region C centers on the wake close to the nucleus, where $|B| \le B_{sw}$. The warm solar wind population is nearly isotropic here, whereas the hot component is composed of solar wind electrons that are accelerated parallel to the magnetic field with $F_\|,sw > F_\perp,sw$. Cometary electrons are nearly absent here. Finally, note that the PSDs here display flat-tops and distributions with a clear knee, rather than monoenergetic peaks at energies comparable to the parallel potential energy, $\Phi_\|$ (Clark et al. 2015).
<figure><img src="image_3.png" /><figcaption>Figure 2. Profiles in the y = 0 plane. (a) Cometary electron density log₁₀(nₑ,c/cm⁻³). (b) Cometary ion density log₁₀(nᵢ,c/cm⁻³). (c) Parallel potential Φₖ (eV). (d) B_y/B_y,sw. (e) Solar wind electron density log₁₀(nₑ,sw/cm⁻³). (f) Solar wind ion density log₁₀(nᵢ,sw/cm⁻³). (g) Tₑ,ₖ (eV). (h) Tₑ,⊥ (eV). Quantities as functions of the electrostatic potential: total electron density (i), B_y/B_y,sw (j), Tₑ,ₖ (k), Tₑ,⊥ (l). Each dot represents a cell in the y = 0 plane. Blue lines are linear regression fits (see the respective titles of panels (i)–(l)).</figcaption></figure>

5. Discussion

In the case of a weakly outgassing comet the solar wind electrons move to balance the positively charged cometary ions (Deca et al. 2017). In order to maintain quasi-neutrality, i.e., to increase the solar wind electron density in regions where the cometary ion density is substantially higher, a strong parallel electric field needs to be generated to locally and temporarily trap electrons. We have shown above that the production of suprathermal (solar wind) electrons is co-located with the region where the parallel electric field component dominates. To constrain the mechanism responsible for the anisotropic electron heating present in the latter region, we introduce an effective acceleration (ambipolar) potential (Egedal et al. 2008, 2010), Φₖ, defined by the integral

$$
\Phi_{\parallel}(\mathbf{x}) = \oint \mathbf{E} \cdot d\mathbf{l}. \qquad (1)
$$

Here the integration is performed along a field line starting
from a point **x** at y = 0 up to the ambient solar wind,
where **E** ⋅ **B** ~ 0. Φ<sub>‖</sub> is a measure for the work performed by
the electric field. Electrons that have a parallel kinetic energy
greater than Φ<sub>‖</sub> that pass through the region gain energy while
traversing into the potential well, and lose it again on their way
out. Those that enter the potential well with a parallel kinetic
energy below the local Φ<sub>‖</sub> remain trapped and bounce along the
magnetic field line while their motion is simultaneously altered
by the local **E** × **B** drift component.

Typically Φₖ is found to scale with Tₑ (Clark et al. 2015; Madanian et al. 2016) and values as high as 7 Tₑ ~ 70 eV are observed in our simulation. Consistent with the spatial profile of Φₖ in the y = 0 plane (Figure 2(c)), suprathermal electrons up to 50–70 eV (depending on the distance) are present (Figure 2(g), also distribution functions at Figures 3(b), (c)) that are accelerated parallel to the **B** field (Figure 2(d)). Both Tₖ,ₐ and Φₖ peak inside the cometary ion-rich plasma, indicative of acceleration by the parallel electric field Eₖ, which is visualized by dashed magenta arrows in Figure 1. The distributions presented in Figure 3 are characteristic for the cometary ion-filled wake and not only for the source region (where most of the injection occurs). It is important to note that the acceleration is not spatially localized and thus no electron prominent beam-like features are formed. First, the ambipolar electric field is roughly governed by the spatial distribution of implanted comet plasma, which in turn has a gradient scale comparable to the distance to the comet. Second, particles with energies less than Φₖ gradually fill the phase space as flux tubes pass slowly (defined by **E** × **B** drift transit time) through the high-Φₖ region close to nucleus. Third, such acceleration is nonturbulent and adiabatic (Egedal et al. 2010), but we do not exclude the possibility of electron beam modes development in a more realistic (but more computationally expensive) simulation with higher resolution and a higher mᵢ,sw/mₑ mass ratio.

Cometary electrons are picked up and follow the magnetic
field lines (i.e., the electron channels seen in Figure 1, left
panel). Downstream of the nucleus, there is a region in which
<figure><img src="image_4.png" /><figcaption>**Figure 3.** Electrons PADs showing the spatial variation of (a) warm ($E < 12.5$ eV), (b) hot ($12.5$ eV $\le E < 25$ eV), and (c) suprathermal ($E \ge 25$ eV) electrons. Differential energy flux is color-coded in arbitrary units. (d) Differential energy flux. Solar wind electrons and cometary electrons are counted separately (left and right panels, respectively). PADs are sampled along a line in the X-direction through the nucleus. Electron PSD at (x, y, z) = (−100,0,0) km (e), (h); (x, y, z) = (0,0,0) km (f), (i); (x, y, z) = (200,0,0) km (g), (j). Panels (e)–(g) display the omnidirectional electron flux. Panels (h)–(j) display the parallel $F_\|$ and the perpendicular $F_\perp$ electron fluxes with respect to the local magnetic field. Panel (k) shows the particle selection criteria for the computation of the fluxes $F_\|$: $\{\theta \le 25^\circ, \theta \ge 155^\circ\}$ and $F_\perp$: $\{|\theta - 90^\circ| \le 10^\circ\}$. $\theta$ is the pitch angle.</figcaption></figure>

the excess of cometary electrons produces an ambipolar potential of the “wrong” sign (in contrast to the more usual excess of ions), visible at right in Figure 2(c) as a region reaching $\Phi_\| \sim -7$ eV.

This region appears to be a potential barrier for solar wind electrons, leading to a reduction of $T_{e,\|}$ (Figure 2(g)). This effect indicates that the ambipolar electric field is not spherically (or cylindrically) symmetric as is often assumed in theoretical considerations of the kinetic expansion of a cometary electron-ion plasma cloud (Nilsson et al. 2018).

It can be seen in Figures 2(a)–(h) that the mass loading by heavy cometary ions produces perturbations in plasma parameters. Figures 2(i)–(l) display the variations of $n_e$, $T_{e,\|}$, $T_{e,\perp}$, and $B_y$ as a function of $\Phi_\|$. The density and parallel electron temperature exhibit a quasi-linear behavior for $-7$ eV $< \Phi_\| < 25$ eV. Potentials above 25 eV are present only very close to the nucleus ($<100$ km) where a large flux of warm (unaccelerated) cometary electrons reduces $T_{e,\|}$. The slope of $T_{e,\perp}(\Phi_\|)$ is considerably flatter than $T_{e,\|}(\Phi_\|)$, confirming that parallel acceleration by $\Phi_\|$ is preferred by the system. Hence, the accelerating potential $\Phi_\|$ is indeed a critical

quantity that controls the trapping of electrons, and by extension the nature of the electron closure relation (equation of state; Le et al. 2010; Deca et al. 2019; Sishtla et al. 2019).

## 6. Conclusions

To conclude, in this study we have presented a three-dimensional fully kinetic simulation of the solar wind interaction with a weakly outgassing comet. Our main results are:

1. The solar wind electron population persists throughout the coma alongside electrons of cometary origin. A hot population with energies up to 50–70 eV (equivalent to several times the average energy of the ambient cometary electrons) near the nucleus constitutes mostly solar wind electrons.

2. We show that the solar wind electrons are accelerated by $\Phi_\|$ and that $T_{e,\|}$ depends linearly on the latter. The total electron density correlates well with the acceleration potential for $\Phi_\|/T_e \le 25$ eV but shows a steeper slope with more scatter for $\Phi_\|/T_e > 25$ eV.

3. In addition to the accelerated electrons that stream parallel to $B$, there exists a population of compressed solar wind
electrons with pitch angles ~90° (for energies from 20 eV up to 100 eV at x = -100 km). They are formed as a response to the magnetic pile-up region upstream of the comet when looking along the Sun-comet line. The presented simulation features a relatively weak outgassing rate, resulting in a peak magnetic field 2.6 times the solar wind magnetic field. We expect that for a higher outgassing rate this number will become greater (Koenders et al. 2013; Rubin et al. 2014) and form more electrons with θ ~ 90°, eventually leading to a (infant) bow shock (Gunell et al. 2018).

Simulations were conducted using resources provided by the Swedish National Infrastructure for Computing (SNIC) at KTH, Stockholm, Sweden, grants m.2016-1-552, m.2016-1-457, and m.2019-3-128. We acknowledge PRACE for awarding us access to Curie at GENCI@CEA, France. Test simulations were performed at the Lomonosov supercomputing facility (Moscow State University) under project Nos. 1576 and 1658. Work of A.D. was supported by RFBR grant 19-02-00993. V.O. was funded by Swedish Research (VR) Council project 2017-04508. Work at LPC2E/CNRS was supported by CNES and by ANR under the financial agreement ANR-15-CE31-0009-01. Access to simulation data can be provided upon request. This work was supported in part by NASA Solar System Exploration Research Virtual Institute (SSERVI): Institute for Modeling Plasmas, Atmosphere, and Cosmic Dust (IMPACT), and the NASA High-End Computing (HEC) Program through the NASA Advanced Supercomputing (NAS) Division at Ames Research Center. J.D. acknowledges partial support from Contract No. JPL-1502225 at the University of Colorado Boulder. Part of this work was inspired by discussions within International Team 402: "Plasma Environment of Comet 67P after Rosetta" at the International Space Science Institute, Bern, Switzerland.

ORCID iDs

Andrey Divin https://orcid.org/0000-0002-5579-3066

Jan Deca https://orcid.org/0000-0002-2272-9818

Giovanni Lapenta https://orcid.org/0000-0002-3123-4024

Vyacheslav Olshevsky https://orcid.org/0000-0001-8467-0041

References

André, M., Odelstad, E., Graham, D. B., et al. 2017, MNRAS, 469, S29

Behar, E., Tabone, B., Saillenfest, M., et al. 2018, A&A, 620, A35

Bieler, A., Altwegg, K., Balsiger, H., et al. 2015, A&A, 583, A7

Brackbill, J. U., & Forslund, D. W. 1982, JCoPh, 46, 271

Broiles, T. W., Burch, J., Chae, K., et al. 2016, MNRAS, 462, S312

Clark, G., Broiles, T. W., Burch, J. L., et al. 2015, A&A, 583, A24

Cravens, T., & Gombosi, T. 2004, AdSpR, 33, 1968

Deca, J., Divin, A., Henri, P., et al. 2017, PhRvL, 118, 205101

Deca, J., Divin, A., Lembège, B., et al. 2015, JGRA, 120, 6443

Deca, J., Divin, A., Lue, C., Ahmadi, T., & Horányi, M. 2018, CmPhy, 1, 12

Deca, J., Henri, P., Divin, A., et al. 2019, PhRvL, 123, 055101

Divin, A., Semenov, V., Korovinskiy, D., et al. 2016, GeoRL, 43, 10565

Egedal, J., Fox, W., Katz, N., et al. 2008, JGRA, 113, A12207

Egedal, J., Lê, A., Zhu, Y., et al. 2010, GeoRL, 37, L10102

Engelhardt, I. A. D., Eriksson, A. I., Vigren, E., et al. 2018, A&A, 616, A51

Eriksson, A. I., Engelhardt, I. A., André, M., et al. 2017, A&A, 605, A15

Galand, M., Héritier, K., Odelstad, E., et al. 2016, MNRAS, 462, S331

Glassmeier, K.-H., Boehnhardt, H., Koschny, D., Kührt, E., & Richter, I. 2007, SSRv, 128, 1

Gombosi, T. I., Burch, J. L., & Horányi, M. 2015, A&A, 583, A23

Gunell, H., Goetz, C., Wedlund, C. S., et al. 2018, A&A, 619, L2

Heritier, K., Galand, M., Henri, P., et al. 2018, A&A, 618, A77

Koenders, C., Glassmeier, K.-H., Richter, I., Motschmann, U., & Rubin, M. 2013, P&SS, 87, 85

Koenders, C., Goetz, C., Richter, I., Motschmann, U., & Glassmeier, K.-H. 2016, MNRAS, 462, S235

Lapenta, G., Brackbill, J. U., & Ricci, P. 2006, PhPl, 13, 055904

Lapenta, G., Markidis, S., Divin, A., Goldman, M., & Newman, D. 2010, PhPl, 17, 082106

Le, A., Egedal, J., Fox, W., et al. 2010, PhPl, 17, 055703

Madanian, H., Cravens, T. E., Rahmati, A., et al. 2016, JGRA, 121, 5815

Markidis, S., Lapenta, G., & Rizwan-uddin 2010, Math. Comput. Simul., 80, 1509

Mendis, D. A., Houpis, H. L. F., & Marconi, M. L. 1985, FCPh, 10, 2

Myllys, M., Henri, P., Galand, M., et al. 2019, A&A, 630, A42

Nilsson, H., Gunell, H., Karlsson, T., et al. 2018, A&A, 616, A50

Nilsson, H., Stenberg Wieser, G., Behar, E., et al. 2015, Sci, 347, aaa0571

Olshevsky, V., Deca, J., Divin, A., et al. 2016, ApJ, 819, 52

Olshevsky, V., Divin, A., Eriksson, E., Markidis, S., & Lapenta, G. 2015, ApJ, 807, 155

Rème, H. 1991, in Cometary Plasma Processes, ed. A. Johnstone, Vol. 61 (Washington, DC: AGU), 87

Rubin, M., Koenders, C., Altwegg, K., et al. 2014, Icar, 242, 38

Sishtla, C. P., Divin, A., Deca, J., Olshevsky, V., & Markidis, S. 2019, PhPl, 26, 102904

Taylor, M. G. G. T., Altobelli, N., Buratti, B. J., & Choukroun, M. 2017, RSPTA, 375, 20160262

Wedlund, C. S., Alho, M., Gronoff, G., et al. 2017, A&A, 604, A73