Self-consistent electrostatic simulations of reforming double layers in the downward current region of the aurora

Abstract

The plasma on a magnetic field line in the downward current region of the aurora is simulated using a Vlasov model. It is found that an electric field parallel to the magnetic fields is supported by a double layer moving toward higher altitude. The double layer accelerates electrons upward, and these electrons give rise to plasma waves and electron phase-space holes through beam–plasma interaction. The double layer is disrupted when reaching altitudes of 1–2 Earth radii where the Langmuir condition no longer can be satisfied due to the diminishing density of electrons coming up from the ionosphere. During the disruption the potential drop is in part carried by the electron holes. The disruption creates favourable conditions for double layer formation near the ionosphere and double layers form anew in that region. The process repeats itself with a period of approximately 1 min. This period is determined by how far the double layer can reach before being disrupted: a higher disruption altitude corresponds to a longer repetition period. The disruption altitude is, in turn, found to increase with ionospheric density and to decrease with total voltage. The current displays oscillations around a mean value. The period of the oscillations is the same as the recurrence period of the double layer formations. The oscillation amplitude increases with increasing voltage, whereas the mean value of the current is independent of voltage in the 100 to 800 V range covered by our simulations. Instead, the mean value of the current is determined by the electron density at the ionospheric boundary.

Introduction

In the downward current region of the aurora, electrons are accelerated upward into space by electric fields that are parallel to the magnetic field. Measurements performed using instruments on the Freja spacecraft showed that beams of upgoing electrons are present within structures of diverging perpendicular electric fields at altitudes down to 800 $\mathrm{km}$, implying that parallel downward electric fields at times exist at even lower altitudes . This was interpreted as a U-shaped potential structure with a net positive charge, giving rise to the diverging fields, which is similar to the upward current region where a net negative charge is the source of converging electric fields. Multi-spacecraft measurements by the Cluster satellites showed that the structure with the diverging fields extends up to about $\mathrm{2}\times {10}^{\mathrm{4}}\mathrm{km}$ in altitude and persists on timescales of several minutes . Data from the Fast Auroral Snapshot (FAST) spacecraft supported a pressure-cooker model of the downward current region, wherein the downward electric field keeps ions confined to low altitudes, where they undergo intense perpendicular heating until they have gained enough energy to be able to escape and contribute to the ion outflow . proposed a modification to the pressure-cooker model by replacing the static electric field by a moving double layer.

observed an electric double layer in the downward current region using the FAST satellite. The double layer was moving upward along the magnetic field; there was a region with wave activity and electron phase-space holes on the high potential side, and there was a gap without significant electric fields or waves between the wave-dominated region and the double layer itself. Similar properties have been reported for double layers in laboratory experiments . Analysis of several double layer encounters by the FAST spacecraft showed that a well-defined gap region is favoured by the presence of a suprathermal electron population on the high potential side of the double layer . This conclusion was supported by Vlasov simulations .

showed that a stationary double layer can exist in a magnetic mirror field configuration, if the polarity is such that the electrons are accelerated in the direction where the magnetic field is stronger. This is the case in the upward current region of the aurora. In the downward current region no stable equilibrium position exists, and therefore double layers will always be in motion in this region. For moving double layers, steady-state theory can sometimes be applied in a moving frame of reference .

In this article we use Vlasov simulations to model the plasma in the downward current region. We study repeatedly reforming double layers and their relations with circuit parameters and phenomena such as electron phase-space holes. A brief description of the simulation model is given in Sect. , and the initial and boundary conditions are treated in Sect. . In Sect.  we report simulations of what a stationary observer on a spacecraft would see when a double layer passes by. Double layer motion and formation are explored in Sect. . The parameter dependence is studied in Sect. , and the conclusions are discussed in Sect. .

Simulation model

In order to model the plasma in the return current region we use the Vlasov simulation code published by . The model is one-dimensional in configuration space and two-dimensional in velocity space. The distribution function is written $f(z,v_z,\mathit{\mu },t)$, where $z$ is the coordinate along the magnetic field direction, $v_z$ is the velocity component in that direction, $\mathit{\mu }$ is the magnetic moment, and $t$ is time. The magnetic moment is an adiabatic invariant; therefore $\dot{\mathit{\mu }}=\mathrm{0}$. The forces on the plasma come from the parallel electric field, the magnetic mirror field and the gravitational field. The Vlasov equation is coupled to a Poisson type equation adapted to the converging magnetic field geometry. Thus, the system of equations that we solve is $$\begin{array}{}& & \frac{\partial f}{\partial t}+v_z\frac{\partial f}{\partial z}+\frac{\mathrm{1}}{m}\left(q,E,-,\mathit{\mu },\frac{\mathrm{d}B}{\mathrm{d}z},+,m,a_{\mathrm{g}}\right)\frac{\partial f}{\partial v_z}=\mathrm{0},& & \frac{\mathrm{d}}{\mathrm{d}z}\left(\frac{B_{\mathrm{S}}}{B},E\right)=\frac{{\mathit{\rho }}_l}{S{\mathit{ϵ}}_{\mathrm{r}}{\mathit{ϵ}}_{\mathrm{0}}},\end{array}$$ where $B_{\mathrm{S}}$ is the magnetic flux density at the reference point, and $S$ is the flux tube cross section at that point. The charge per unit length of the flux tube is denoted ${\mathit{\rho }}_l$ in Eq. (), and it is computed as a sum of the integrals of the distribution function over all of velocity space for all species $s$: $${\mathit{\rho }}_l=\sum _s{q}_s\int f_s(v_z,\mathit{\mu })\mathrm{d}\mathit{\mu }\mathrm{d}{v}_z.$$ The constant ${\mathit{ϵ}}_{\mathrm{r}}$ that appears on the right-hand side of Eq. () is artificial, and it was introduced to reduce the computational effort. Since ${\mathit{\lambda }}_{\mathrm{D}}\sim \sqrt{{\mathit{ϵ}}_{\mathrm{r}}}$ and ${\mathit{\omega }}_{\mathrm{p}}\sim \mathrm{1}/\sqrt{{\mathit{ϵ}}_{\mathrm{r}}}$, larger grid cells and longer time steps may be used. This means that the widths of double layers and phase space holes will be overestimated by a factor of $\sqrt{{\mathit{ϵ}}_{\mathrm{r}}}$, but as long as these widths are much smaller than the typical length scales of the overall changes of the plasma properties the exact values of the widths are not important for the results of the simulation. In order to test this method performed a series of four simulation runs, successively decreasing ${\mathit{ϵ}}_{\mathrm{r}}$ from $4.98\times {10}^{\mathrm{8}}$ down to $4.98\times {10}^{\mathrm{4}}$, while observing how the gradients got sharper as this series of runs converged to a solution (see Fig. 3c of that paper). In the simulations reported here, ${\mathit{ϵ}}_{\mathrm{r}}=8100$ is used; the time step is $\mathrm{\Delta }t=1.0\times {10}^{-\mathrm{5}}\mathrm{s}$; the grid is non-uniform with the smallest grid cell size, $\mathrm{\Delta }z=622\mathrm{m}$, at the ionosphere and the largest, $\mathrm{\Delta }z=1.03\times {10}^{\mathrm{4}}\mathrm{m}$, at the magnetospheric equator.

A stationary magnetic field is prescribed. We use a magnetic dipole model, where we approximate the $L=\mathrm{7}$ shell $$\begin{array}{ccc}{\displaystyle }{B}_{\mathrm{z}}\left(z\right)& {\displaystyle }=& {\displaystyle }{B}_{\mathrm{M}}\exp \left[{\left({\displaystyle {\displaystyle \frac{z}{L_{\mathrm{z}}}}}\right)}^{\mathrm{2}},\left(\ln ,\left({\displaystyle {\displaystyle \frac{B_{\mathrm{I}}}{B_{\mathrm{M}}}}}\right),-,0.6\right.\right.\\ & {\displaystyle }& {\displaystyle }& {\displaystyle }\left..\left..-,1.8,{\left({\displaystyle {\displaystyle \frac{z}{L_{\mathrm{z}}}}}\right)}^{\mathrm{2}},+,2.4,{\left({\displaystyle {\displaystyle \frac{z}{L_{\mathrm{z}}}}}\right)}^{\mathrm{6}}\right)\right],\end{array}$$ where $B_{\mathrm{I}}=56$ $\mathrm{\mu }\mathrm{T}$ and $B_{\mathrm{M}}=0.0864$ $\mathrm{\mu }\mathrm{T}$ are the magnetic flux densities at the ionospheric and the magnetospheric ends of the system respectively, and $L_{\mathrm{z}}=5.5\times {10}^{\mathrm{7}}$ $\mathrm{m}$ is the length of the system. The $z$ axis is defined so that $z=\mathrm{0}$ at the magnetospheric end of the system, and $z=L_{\mathrm{z}}$ at the ionosphere, which, in this model, is at an altitude of 120 $\mathrm{km}$.

Alfvén waves are often observed in the auroral zone e.g., and their contribution to auroral acceleration differs between different regions . Observations by the Cluster spacecraft have shown both Alfvénic and electrostatic aurora and combinations of the two along a single spacecraft path . The simulation model used in the present paper is electrostatic ($\mathrm{\nabla }\times \mathbf{E}=-\partial \mathbf{B}/\partial t=\mathrm{0}$), and therefore it cannot be used to study Alfvén waves. We are thus limited to study electrostatic aurora and the electrostatic aspects of aurora that have both electrostatic and Alfvénic components.

For more information about the simulation model, see the paper by , which also includes the Fortran code itself. The same code has also been used to study trapping and loss of electrons in the upward current region of the aurora and to assess the possibility of constructing a laboratory setup to model auroral acceleration .

Initial and boundary conditions

In the simulations reported in the present paper, we have modelled the plasma in the return current region using the parameters that are shown in Table  as boundary conditions.

Parameters used in the simulations at the magnetospheric and ionospheric boundaries. The multiple values given for $V_{\mathrm{p}}$ and $n$ at the ionosphere correspond to different simulation runs.

	Magnetosphere	Ionosphere
$z$	0	$5.5\times {10}^{\mathrm{7}}$ $\mathrm{m}$
$B$	0.086 $\mathrm{\mu }\mathrm{T}$	56 $\mathrm{\mu }\mathrm{T}$
$k_{\mathrm{B}}{T}_{\mathrm{e}}$	500 $\mathrm{eV}$	50 $\mathrm{eV}$
$k_{\mathrm{B}}{T}_{{\mathrm{H}}^{+}}$	2500 $\mathrm{eV}$	50 $\mathrm{eV}$
$V_{\mathrm{p}}$	$\mathrm{0}\mathrm{V}$	($-$100, $-$400, $-$800) $\mathrm{V}$
$n$	$\mathrm{3}\times {10}^{\mathrm{5}}$ ${\mathrm{m}}^{-\mathrm{3}}$	$(\mathrm{1},\mathrm{2},\mathrm{4})\times {10}^{\mathrm{7}}$ ${\mathrm{m}}^{-\mathrm{3}}$

The plasma potential, $V_{\mathrm{p}}$, at the magnetospheric end of the system is always 0 $\mathrm{V}$. The total voltage is controlled by changing the potential of the ionosphere. Three simulation runs were performed with $V_{\mathrm{p},\mathrm{I}}=-100$ $\mathrm{V}$, $-$400, and $-$800 $\mathrm{V}$ respectively. In these runs the density at the ionospheric boundary was $n_{\mathrm{I}}=\mathrm{2}\times {10}^{\mathrm{7}}$ ${\mathrm{m}}^{-\mathrm{3}}$. Another series of three runs was conducted for $n_{\mathrm{I}}=\mathrm{1},\mathrm{2}$, and $\mathrm{4}\times {10}^{\mathrm{7}}$ ${\mathrm{m}}^{-\mathrm{3}}$, keeping $V_{\mathrm{p},\mathrm{I}}$ at $-$400 $\mathrm{V}$. The run with $V_{\mathrm{p},\mathrm{I}}=-400$ $\mathrm{V}$ and $\mathrm{2}\times {10}^{\mathrm{7}}$ ${\mathrm{m}}^{-\mathrm{3}}$ belongs to both series. The parameters in Table  are chosen to reproduce the conditions of the double layer observation by . The run that comes closest is the one with $V_{\mathrm{p},\mathrm{I}}=-400$ $\mathrm{V}$ and $n_{\mathrm{I}}=\mathrm{4}\times {10}^{\mathrm{7}}$ ${\mathrm{m}}^{-\mathrm{3}}$. We present that run in detail in the following sections, and then we use the complete series to compare the results for different voltages and ionospheric densities.

Electrons and ions are included in the simulations, and all ions that we consider are protons. This means that we do not include phenomena that are caused by the presence of different ion species such as the hydrogen–oxygen ion instabilities reported by in simulations of the upward current region. Particle populations entering the simulated region from the two boundaries are treated as different species. Thus we have four species, namely, electrons from the magnetosphere; ions from the magnetosphere; electrons from the ionosphere, and ions from the ionosphere. Phase space densities for these species are shown in Fig. c–f.

Simulation state at $t=130\mathrm{s}$ in the run with $V_{\mathrm{p},\mathrm{I}}=-400$ $\mathrm{V}$ and $\mathrm{4}\times {10}^{\mathrm{7}}$ ${\mathrm{m}}^{-\mathrm{3}}$. (a) Plasma potential as a function of $z$. (b) Densities as functions of $z$. The thick blue curve shows the plasma density. The thin solid curves show protons (red) and electrons (blue) from the magnetospheric end of the system. The dashed curves show protons (red) and electrons (blue) originating from the ionosphere. (c)–(f) Phase space densities for (c) magnetospheric electrons; (d) magnetospheric protons; (e) ionospheric electrons; and (f) ionospheric protons. The colour scales have been normalised so that integrals over all $v_z$ yield $n_s/B$. The unit for $f(z,v_z)$ is ${\mathrm{m}}^{-\mathrm{4}}{\mathrm{T}}^{-\mathrm{1}}\mathrm{s}$.

[Figure omitted. See PDF]

The figure shows the simulation state at $t=130\mathrm{s}$. At $t=\mathrm{0}$ the system was filled with a constant density ($n=\mathrm{3}\times {10}^{\mathrm{5}}$ ${\mathrm{m}}^{-\mathrm{3}}$) of the magnetospheric species (the species that have their source at $z=\mathrm{0}$). The initial density of the ionospheric species was zero throughout the system. For $t\mathit{>}\mathrm{0}$ both ions and electrons from the ionosphere are allowed to enter the system, which is populated as is seen in Fig. e and f. Here and throughout this paper, the colour scales showing distributions $f_s(z,v_z)$ have been normalised so that integrals over all $v_z$ yield $n_s/B$. The unit for $f(z,v_z)$ is ${\mathrm{m}}^{-\mathrm{4}}{\mathrm{T}}^{-\mathrm{1}}\mathrm{s}$.

reported parallel electron temperatures in the 20–100 $\mathrm{eV}$ range when the FAST spacecraft was between the double layer and the ionosphere, while the perpendicular temperature was below 1 $\mathrm{eV}$. In the ionosphere, both the perpendicular and the parallel temperatures are below 1 $\mathrm{eV}$. The observations indicate that waves heat the plasma at low altitudes and that there is transfer from perpendicular to parallel energy as the plasma moves upward into a weaker magnetic field. The heating involves processes that we cannot model self-consistently with an electrostatic code. We therefore set the temperature of the ionosphere to 50 $\mathrm{eV}$ as a boundary condition, and that produces temperatures that are within the range of the observations at the altitudes where these were made.

Although the heating that occurs in the simulation is insufficient to heat a cold ionosphere to 50 $\mathrm{eV}$, some heating of the ionospheric species can be seen in Fig. e and f. Also, some of the magnetospheric plasma is able to reach the ionosphere. This means that there is a mismatch between the distribution that is used as a boundary condition at the ionosphere and that of the plasma just above it. As a result a sheath of a few kilometres in width develops at the boundary. It is too thin to be completely resolved on the scale of Fig. , but it can be seen in panel b as the increase of the plasma density (thick blue line) from $1.3\times {10}^{\mathrm{7}}$ to $\mathrm{2}\times {10}^{\mathrm{7}}$ ${\mathrm{m}}^{-\mathrm{3}}$ at the right-hand side of the figure. In space, there must be a transition region between the cold ionosphere and the hotter plasma just above it. In the simulation, this transition region is represented by the sheath. The $z$ dependence of the simulated quantities in the sheath is unlikely to be a good representation of this dependence in space, but most of the plasma, from $z=\mathrm{0}$ to $z=5.499\times {10}^{\mathrm{7}}\mathrm{m}$ (altitude 130 $\mathrm{km}$), is outside the sheath and can be represented by this model. The effect the sheath has on the rest of the system is to reduce the density that is available at the ionospheric boundary, and that aspect is discussed in Sect. . At the start of the simulation the density gradient at the ionosphere is particularly large. The density changes two orders of magnitude over one grid cell. This causes the sheath to carry a larger voltage of approximately 200 $\mathrm{V}$, but it is a transient that lasts less than a second.

Double layers seen by a stationary observer

It is seen in Fig.  that there is a double layer at $z\approx 5.3\times {10}^{\mathrm{7}}\mathrm{m}$ (altitude $2.1\times {10}^{\mathrm{3}}\mathrm{km}$), where there is a sharp drop in the plasma potential. Double layers are space charge structures embedded in the plasma that can carry a large potential difference. They are named double layers because they consist of at least two layers of different net charge. See for a review of double layer physics. The double layer voltage is about 700 $\mathrm{V}$, which is more than the total voltage over the system. There is a potential minimum on the low potential side, and there is a region of positive potentials on the high potential side. This region corresponds to the bump in the plasma density curve in Fig. b between $z\approx \mathrm{4}\times {10}^{\mathrm{7}}\mathrm{m}$ and $z\approx 5.3\times {10}^{\mathrm{7}}\mathrm{m}$. The presence of waves in this region can be seen both in the plasma potential (Fig. a) and in the phase space diagrams (Fig. c–f).

The double layer separates the cold ionospheric from the hot magnetospheric plasma. Electrons from the ionosphere are accelerated upward in the double layer (Fig. e). The double layer also affects the distribution of electrons that come from the magnetosphere, which is seen by the sharp cutoff at $z\approx 5.3\times {10}^{\mathrm{7}}\mathrm{m}$ in Fig. c. The influence of the double layer on the magnetospheric ions is insignificant, as is seen in Fig. d. The reason for this is that for the hot magnetospheric ions the voltage equivalent of the temperature is $k_{\mathrm{B}}{T}_{\mathrm{i},\mathrm{M}}/e=2500\mathrm{V}$, which is greater than the double layer voltage. In panel f ions of ionospheric origin are seen also on the high potential side of the double layer. At $t=130\mathrm{s}$ the part of this population that is in contact with the double layer is being accelerated downward through it, while those at higher altitudes (lower $z$) continue to move upward.

Electric field, current, and density as functions of time in the run with $V_{\mathrm{p},\mathrm{I}}=-400$ $\mathrm{V}$ and $\mathrm{4}\times {10}^{\mathrm{7}}$ ${\mathrm{m}}^{-\mathrm{3}}$. (a) The electric field in a $z-t$ diagram for $z\ge 3.5\times {10}^{\mathrm{7}}\mathrm{m}$. (b) Current density scaled to the ionospheric side of the system. (c) Electric field at $z=5.265\times {10}^{\mathrm{7}}\mathrm{m}$. (d) Plasma density at $z=5.265\times {10}^{\mathrm{7}}\mathrm{m}$. The line at $z=5.265\times {10}^{\mathrm{7}}\mathrm{m}$ in panel (a), corresponding to an altitude of 2470 $\mathrm{km}$, indicates the position for which the electric field and density are shown in panels (c) and (d).

[Figure omitted. See PDF]

A $z-t$ diagram of the electric field is shown in Fig. a. Only the region $z\ge 3.5\times {10}^{\mathrm{7}}\mathrm{m}$ is shown, since the electric field is close to zero for $z\mathit{<}3.5\times {10}^{\mathrm{7}}\mathrm{m}$. Double layers can be seen as thin lines extending diagonally to $z\approx 4.5\times {10}^{\mathrm{7}}\mathrm{m}$. This means that they move toward higher altitude until they no longer can be maintained. Then a new double layer forms close to the ionosphere, and that newly created double layer moves toward higher altitudes and the process repeats itself. A region of fast varying electric fields is seen on the high potential side of the double layer, for example around $t=130\mathrm{s}$, which is the time shown in Fig. . These fields correspond to waves that are generated by electron–beam plasma interaction. The wave field is not resolved on the spatial and temporal scales of Fig. a, which is designed to show the large-scale motion of the double layer. In a disruption, waves first appear on both sides of the double layer, and then there is period when no distinct double layer is seen, before the new double layer is formed. We identify the double layer that is seen as a distinct thin line with what has been called a laminar double layer, and the indistinct, wave-dominated, structure with a turbulent double layer .

Figure b shows the current density, scaled to the ionosphere. The current increases during periods where there is one laminar double layer in the system, and during the disruptions and turbulent periods the current is disrupted too, and oscillations occur.

The horizontal black line at $z=z_{\mathrm{1}}=5.265\times {10}^{\mathrm{7}}\mathrm{m}$ in Fig. a indicates the $z$ coordinate for which the electric field $E_{\mathrm{1}}$ and plasma density $n_{\mathrm{1}}$ are shown in panels c and d respectively. $z_{\mathrm{1}}=5.265\times {10}^{\mathrm{7}}\mathrm{m}$ corresponds to an altitude of 2470 $\mathrm{km}$, which is the altitude of the FAST double layer observation by . A double layer passes $z_{\mathrm{1}}$ at $t=130.4\mathrm{s}$, just after the time shown in Fig. . In the interval $110\mathrm{s}\mathit{\lesssim }t\mathit{\lesssim }130\mathrm{s}$ waves are seen in both the electric field trace in panel c and in the density in panel d. The passage of the double layer brings $z_{\mathrm{1}}$ to the low potential side, where the density is higher as is seen in Fig. b. This explains the sharp rise in density at $t=130\mathrm{s}$. A similar event happened at $t=50\mathrm{s}$, and one is about to start at the end of the simulation at $t=180\mathrm{s}$.

Figure  shows a closeup of the double layer passage at $z=z_{\mathrm{1}}=5.265\times {10}^{\mathrm{7}}\mathrm{m}$ and $t=130\mathrm{s}$ (altitude 2470 $\mathrm{km}$). In the left-hand column, panels a–c, the distribution function for electrons of ionospheric origin, the electric field, and the plasma potential are shown as functions of space for $t=130\mathrm{s}$, and in the right-hand column, panels d–f, the density, electric field, and plasma potential are shown as functions of time for $z=z_{\mathrm{1}}$. The double layer is seen as the potential drop of about 700 $\mathrm{V}$ at the same location as the large positive electric field peaking at 11 $\mathrm{mV}/\mathrm{m}$. In Fig. a the electron beam that forms by acceleration in the double layer is seen. There is an electron phase space hole at the double layer and more are being formed further downstream, where $z\mathit{<}{z}_{\mathrm{1}}$. The width of the double layer and the widths of the phase space holes is approximately 100 $\mathrm{km}$. The simulation was run with ${\mathit{ϵ}}_{\mathrm{r}}=8100$, and since the width of double layers and phase space holes scale with the Debye length, which is proportional to $\sqrt{{\mathit{ϵ}}_{\mathrm{r}}}=90$, the predicted widths in space are about 1 $\mathrm{km}$. The electric field scales in the same way, and thus the electric field in the double layer should be about 1 V m${}^{-\mathrm{1}}$, which is similar to the 0.85 V m${}^{-\mathrm{1}}$ reported by . Although the absolute values of the density in Fig. d are lower than those seen by by a factor of 2–3, the behaviour of the density as a function of time is similar to those observations. A density minimum was detected just before the double layer passed the spacecraft, and an increase to higher values occurred afterwards.

Figure e shows the electric field seen by a stationary observer at $z=z_{\mathrm{1}}$. First, there are waves with frequencies on the order of the local plasma frequency. Then the double layer passage can be seen as a positive peak around $t=130.4\mathrm{s}$. Comparing panels b and e, we may see the double layer as a 100 $\mathrm{km}$ wide structure that passes a stationary observer at a speed of approximately 200 km s${}^{-\mathrm{1}}$, giving rise to an electric field pulse lasting half a second. The subsequent increase in density that is seen in Fig. d happens on a similar 1 $\mathrm{s}$ timescale. The waves on the high potential side are much faster. Only the double layer itself can be seen as a quasi-stationary structure in a moving frame of reference. Including the waves on the high potential side, we have an entity that is dynamic in nature.

The double layer passing $z=z_{\mathrm{1}}=5.265\times {10}^{\mathrm{7}}\mathrm{m}$ (altitude 2470 $\mathrm{km}$) in the run with $V_{\mathrm{p},\mathrm{I}}=-400$ $\mathrm{V}$ and $\mathrm{4}\times {10}^{\mathrm{7}}$ ${\mathrm{m}}^{-\mathrm{3}}$. The left column shows plasma properties at $t=130\mathrm{s}$ as functions of $z$: (a) distribution function $f(z,v_z)$ in units of ${\mathrm{m}}^{-\mathrm{4}}{\mathrm{T}}^{-\mathrm{1}}\mathrm{s}$ for electrons of ionospheric origin. (b) Electric field; and (c) plasma potential. The right column shows plasma properties at $z=z_{\mathrm{1}}$ as functions of time: (d) plasma density; (e) electric field; and (f) plasma potential.

[Figure omitted. See PDF]

Double layer motion and formation

Our simulations suggest that double layers form in the plasma above the ionosphere, they move toward higher altitudes, where they subsequently disrupt, and a new double layer forms close to the ionosphere. In Fig. , a sequence of images of this development is shown. Each row shows (a) the plasma potential $V_{\mathrm{p}}\left(z\right)$, (b) the phase space density of the ionospheric electrons, and (c) the phase space density of the ionospheric ions. There is one row for every 10 s from $t=50\mathrm{s}$ to $t=120\mathrm{s}$. The subscript of the panel label indicates the number of seconds for which the distribution function is shown in that panel.

Plasma potential and phase space densities for the ionospheric species in the run with $V_{\mathrm{p},\mathrm{I}}=-400$ $\mathrm{V}$ and $\mathrm{4}\times {10}^{\mathrm{7}}$ ${\mathrm{m}}^{-\mathrm{3}}$. (a) Plasma potential $V_{\mathrm{p}}\left(z\right)$. (b) Phase space density $f_{\mathrm{e},\mathrm{I}}(z,v_z)$ of ionospheric electrons. (c) Phase space $f_{\mathrm{i},\mathrm{I}}(z,v_z)$ of ionospheric ions. The colour scales have been normalised so that integrals over all $v_z$ are equal to $n_{\mathrm{s}}/B$, where $n_{\mathrm{s}}$ denotes the density of species in question, and the unit for $f(z,v_z)$ is ${\mathrm{m}}^{-\mathrm{4}}{\mathrm{T}}^{-\mathrm{1}}\mathrm{s}$. The subscript of the panel label indicates the number of seconds for which the distribution function is shown in that panel.

[Figure omitted. See PDF]

Figure  follows a newly formed double layer at $t=50\mathrm{s}$ through its motion toward higher altitudes to its eventual breakup at $t=100\mathrm{s}$, the subsequent reformation of a double layer at low altitude around $t=110\mathrm{s}$, and the early stage of the upward motion of that double layer at $t=120\mathrm{s}$. From $t=70\mathrm{s}$ to $t=90\mathrm{s}$ the double layer moves at a constant speed of $v_{\mathrm{DL}}=21$4 km s${}^{-\mathrm{1}}$ as estimated in Fig. f.

Comparing $t-z$ diagrams of $E$ in different simulation runs. (a–c) Constant ionospheric density, $\mathrm{2}\times {10}^{\mathrm{7}}{\mathrm{m}}^{-\mathrm{3}}$, and three different voltages, namely 100, 400, and 800 $\mathrm{V}$. (d–f) Constant voltage, 400 $\mathrm{V}$, and three different densities at the ionospheric boundary, 1, 2, and $\mathrm{4}\times {10}^{\mathrm{7}}{\mathrm{m}}^{-\mathrm{3}}$. Panels (b) and (e) show the same run. The standard run presented in the other figures corresponds to panel (f). Examples of the period between double layers passing $z=\mathrm{5}\times {10}^{\mathrm{7}}\mathrm{m}$ (altitude $5.1\times {10}^{\mathrm{3}}$ $\mathrm{km}$) and of double layer velocities are shown in each panel.

[Figure omitted. See PDF]

In the electron phase space diagrams in Fig. b${}_{80}$ and b${}_{90}$ there is a laminar flow of electrons through the double layer, and these electrons move as a beam through a short region on the high potential side before the beam–plasma instability has caused a large enough wave growth to create electron phase-space holes. The region where there is a distinct beam is known as the gap region due to the relatively small electric fields found there. The gap is also seen in Fig. a${}_{90}$, where the plasma potential curve is flat immediately to the left of the double layer.

A potential minimum is present on the low potential side of the double layer, and it is followed by an expanding ionospheric plasma. The explanation, provided by , for this behaviour is that the ions on the low potential side move more slowly than the double layer. The expansion of the ionospheric ion population would proceed at approximately the ion acoustic speed, which in our case is $c_{\mathrm{s}}=\sqrt{k_{\mathrm{B}}(T_{\mathrm{e}}+\mathrm{3}{T}_{\mathrm{i}})/m_{\mathrm{i}}}\approx 140$ km s${}^{-\mathrm{1}}$, and less than $v_{\mathrm{DL}}=214$ km s${}^{-\mathrm{1}}$ estimated in Fig. f. While the ions are slower than the double layer the electrons are faster, and they can keep up with its speed. This causes a net excess of electrons at the low potential side of the double layer creating a potential dip at the foot of the double layer. The potential of the dip is lower than the potential of the ionosphere, and therefore there is an ambipolar electric field between the double layer and the ionosphere. This ambipolar field pulls the ions upward, as can be seen in Fig. c, where $v_z$ for the ions gets more negative as they come out of the right-hand boundary of the system toward the double layer.

When the double layer is disrupted and a new one forms at low altitude, the bulk of this ion population is unaffected and continues to move upwards. Thus, the ions that are present at $5.2\times {10}^{\mathrm{7}}\mathrm{m}\le z\le 5.5\times {10}^{\mathrm{7}}\mathrm{m}$ in Fig. c${}_{50}$ can be found in the range $\mathrm{4}\times {10}^{\mathrm{7}}\mathrm{m}\le z\le 5.3\times {10}^{\mathrm{7}}\mathrm{m}$ in Fig. c${}_{120}$, and those that are found at $z\mathit{<}5.2\times {10}^{\mathrm{7}}\mathrm{m}$ in Fig. c${}_{50}$ are spread out over a large $z$ range, $z\mathit{<}\mathrm{4}\times {10}^{\mathrm{7}}\mathrm{m}$, in Fig. c${}_{120}$. Although the total voltage in the downward current region acts to accelerate ions downward, the motion of the double layer and its repeated reformation move a fraction of the ions in the opposite direction. However, some of these ions are caught by the new double layer and accelerated toward the ionosphere forming an ion beam that moves through the low potential plasma. This acceleration is seen in all panels of Fig. c with the most obvious beam populations in Fig. c${}_{90}$ and c${}_{100}$. The outflowing ions and those accelerated back downward are also seen in Fig. f.

The ion beam on the low potential side causes waves that are seen in the ion phase space in Fig. c and may contribute to the disruption of the double layer . Two competing double layers are present between $t=55\mathrm{s}$ and $t=65\mathrm{s}$ as the electric field $z-t$ diagram in Fig. a shows. At $t=60\mathrm{s}$, which is the time shown in Fig. a${}_{60}$ the higher $z$ (lower altitude) double layer carries only little voltage. It can still be seen in the ion population in Fig. c${}_{60}$ that ions have been accelerated at $z\approx 5.3\times {10}^{\mathrm{7}}\mathrm{m}$ a few seconds before that time, and also the phase space density of the electrons (Fig. b${}_{60}$) changes at the same position. Finally, it was the double layer at higher $z$ (lower altitude) that survived, and at $t=70\mathrm{s}$ it is the only one that remains. This phenomenon was also observed in simulations by , who saw double layers first being destroyed by Buneman instabilities and then form anew on the low potential side.

showed that a condition for the existence of a double layer is that $${\displaystyle \frac{i_{\mathrm{e}}}{i_{\mathrm{i}}}}=\sqrt{{\displaystyle \frac{m_{\mathrm{i}}}{m_{\mathrm{e}}}}},$$ where $i_{\mathrm{e}}$ and $i_{\mathrm{i}}$ are the electron and ion currents through the double layer respectively. This condition only applies to strong double layers, where the double layer voltage is much greater than the voltage equivalent of the electron and ion temperatures. When the temperatures are not negligible a generalised Langmuir condition can be derived, requiring that the total pressure – that is to say, dynamic plus thermal – is constant across the double layer . Applying this condition in a frame of reference that moves with the same velocity, $v_{\mathrm{DL}}$, as the double layer we have $$p_z(z_{\mathrm{1}},v_{\mathrm{DL}})=p_z(z_{\mathrm{2}},v_{\mathrm{DL}}),$$ where $z_{\mathrm{1}}$ and $z_{\mathrm{2}}$ are coordinates on the two respective sides of the double layer and $$p_z(z,v_{\mathrm{DL}})=\sum _s{m}_s\int (v_z-v_{\mathrm{DL}}{)}^{\mathrm{2}}{f}_s(z,v_z)\mathrm{d}{v}_z$$ is the field-aligned pressure. The application of the generalised Langmuir condition to the double layer at $t=\mathrm{3}\mathrm{s}$ is illustrated in Fig. .

Application of the generalised Langmuir condition at $t=\mathrm{3}\mathrm{s}$ in the run with $V_{\mathrm{p},\mathrm{I}}=-400$ $\mathrm{V}$ and $\mathrm{4}\times {10}^{\mathrm{7}}$ ${\mathrm{m}}^{-\mathrm{3}}$. The panels show (a) the plasma potential; (b) the field-aligned pressure $p_z$ in the moving frame (solid line) and in the stationary frame (dashed line); and (c) $p_z\left(v_{\mathrm{DL}}\right)$ for $z_{\mathrm{1}}=5.34\times {10}^{\mathrm{7}}\mathrm{m}$ (red) and $z_{\mathrm{2}}=5.41\times {10}^{\mathrm{7}}\mathrm{m}$ (blue). These positions are also marked in panels (a) and (b). The vertical dashed line in panel (c) marks $v_{\mathrm{DL}}=331$ km s${}^{-\mathrm{1}}$, which is the velocity of the moving frame used to plot $p_z\left(z\right)$ in panel (b).

[Figure omitted. See PDF]

Panel (a) shows the plasma potential as a function of $z$, and the points $z=z_{\mathrm{1}}=5.34\times {10}^{\mathrm{7}}\mathrm{m}$ (altitude 1720 $\mathrm{km}$) and $z=z_{\mathrm{2}}=5.41\times {10}^{\mathrm{7}}\mathrm{m}$ (altitude 1020 $\mathrm{km}$) on each side of the double layer have been marked in red and blue respectively. The velocity of the double layer was $v_z=v_{\mathrm{DL}}=-331$ km s${}^{-\mathrm{1}}$ as determined from Fig. f. The solid line in Fig. b shows the field-aligned pressure $p_z$ in the moving frame given by Eq. () for $v_{\mathrm{DL}}=-331$ km s${}^{-\mathrm{1}}$. For comparison the dashed line shows $p_z$ computed for $v_{\mathrm{DL}}=\mathrm{0}$. We see that the pressure can be balanced across the double layer to meet the generalised Langmuir condition in the moving frame of reference. The pressures on either side of the double layer are marked with red and blue crosses in Fig. b. Fig. c shows $p_z(z_{\mathrm{1}},v_{\mathrm{DL}})$ (red line) and $p_z(z_{\mathrm{2}},v_{\mathrm{DL}})$ (blue line) as functions of $v_{\mathrm{DL}}$. The two pressures are nearly equal over a wide range of velocities.

While there is a wide velocity range where the generalised Langmuir condition can be satisfied when the double layer is close to the ionosphere, this becomes increasingly difficult as it moves to higher altitudes. The densities of the species that originate at the magnetospheric end of the system is constant, whereas the ionospheric species become less dense as the flux tube expands, as is seen in Fig. b. Thus, there is an altitude where the density of the ionospheric electrons is too low to counterbalance the pressure of the magnetospheric ions that have been accelerated through the double layer. When this happens the double layer can no longer exist. In Fig.  this is seen at $t=100\mathrm{s}$, where instead of a double layer the voltage is assumed by a series of electron phase-space holes. At the same time the depth of the dip on the low potential side decreases; the electrons near the ionosphere are no longer held back by the electric field, and that leads to a faster relative drift between electrons and ions, creating conditions favourable for double layer formation in the plasma above the ionosphere. The altitude at which the new double layers form can be determined by enlarging Fig. . While there are instances when the formation is at or almost at the ionospheric boundary, in most cases formation takes place in the range $5.4\times {10}^{\mathrm{7}}\mathrm{m}\mathit{\lesssim }z\mathit{\lesssim }5.48\times {10}^{\mathrm{7}}\mathrm{m}$, which corresponds to altitudes between 300 and 1100 $\mathrm{km}$. Occasionally, reformation occurs at even higher altitudes, such as at $t\approx 60\mathrm{s}$ in Fig. a, where the altitude of the newly formed double layer is about 2000 $\mathrm{km}$.

To further illustrate the double layer motion, videos showing the development of the plasma potential and distributions of the two ionospheric species for each case in Fig.  have been deposited with this article as Supplement.

Parameter dependence

In order to study how the double layer motion depends on total voltage and ionospheric density we have made five different simulation runs with different parameters. In Fig. , $t-z$ diagrams of the electric field are shown for these runs. In Fig. a, b, and c the voltage is 100, 400, and 800 $\mathrm{V}$ respectively, while the density at the ionospheric boundary is kept at $n_{\mathrm{I}}=\mathrm{2}\times {10}^{\mathrm{7}}{\mathrm{m}}^{-\mathrm{3}}$. In Fig. d, e, and f the voltage is constant at 400 $\mathrm{V}$, and the ionospheric density is $n_{\mathrm{I}}=\mathrm{1},\mathrm{2}$, and $\mathrm{4}\times {10}^{\mathrm{7}}{\mathrm{m}}^{-\mathrm{3}}$ in the three panels respectively. Panels (b) and (e) show the same run, it being part of both series. Examples of double layer velocities and periods between double layers passing $z=\mathrm{5}\times {10}^{\mathrm{7}}\mathrm{m}$ (altitude $5.1\times {10}^{\mathrm{3}}\mathrm{km})$ are shown in each panel.

The recurrence period is seen to decrease with increasing voltage in Fig. a–c and it increases with increasing ionospheric density in Fig. d–f. The double layer velocity is approximately $-$200 km s${}^{-\mathrm{1}}$, but it varies both between the different simulation runs and between different intervals within the same run, and no clear trend in parameter dependence can be seen. With approximately the same double layer velocity the different recurrence frequencies can be explained by the difference in the distance the double layer moves before it is disrupted. In the experiments with long recurrence periods shown in Fig. a and f, the double layers reach higher altitudes before they are disrupted than in the runs shown in Fig. c and d.

A higher density at the ionosphere leads to a higher density at higher altitude, which allows the double layer to move farther toward lower $z$ values before reaching the point where the Langmuir condition no longer can be satisfied. For quantitative accuracy one needs to use the generalised Langmuir condition in Eq. (), but qualitatively, reasoning based on the standard Langmuir condition (Eq. ) can be used to explain the density dependence.

Density of the ionospheric electrons in the 100 $\mathrm{V}$ simulation at $t=215.5\mathrm{s}$ (black line) and in the 800 $\mathrm{V}$ simulation at $t=138.5\mathrm{s}$ (red line).

[Figure omitted. See PDF]

For the voltage dependence it is necessary to base even a qualitative argument on the generalised Langmuir condition, because the standard Langmuir condition does not depend on the double layer voltage at all. The dominating species are ions from the magnetosphere that are accelerated toward higher $z$ values in the double layer and electrons from the ionosphere that are accelerated in the opposite direction. The double layer provides the same amount of energy to an ion as to an electron, but its effect is quite different on the density of the different species. Hence, pressure balance becomes a question of density only. The density of ions from the magnetosphere is almost unaffected by the double layer, as is seen by the solid red line in Fig. b. This is a result of their temperature ($k_{\mathrm{B}}{T}_{\mathrm{i},\mathrm{M}}=2500\mathrm{eV}$) being much higher than the double layer voltage in all the simulated cases. For the colder plasma from the ionosphere ($k_{\mathrm{B}}{T}_{\mathrm{e},\mathrm{I}}=50\mathrm{eV}$) this is not the case, and the density of the ionospheric electrons on the high potential side of the double layer is reduced more in the 800 $\mathrm{V}$ case than in the 100 $\mathrm{V}$ case. This is illustrated in Fig. , which shows the density of ionospheric electrons in the two cases for times when there was a double layer at $z=5.16\times {10}^{\mathrm{7}}\mathrm{m}$ (altitude 3520 $\mathrm{km}$). This behaviour is a result of flux continuity, which, for a cold population, makes density inversely proportional to velocity. Because the density is reduced more for higher voltages, the contribution of ionospheric electrons to $p_z(z_{\mathrm{1}},v_{\mathrm{DL}})$ in Eq. () is diminished more, and this renders the double layer unable to reach as far from the ionosphere as it can for lower voltages.

The field-aligned current densities, scaled to the ionosphere, are shown for the different runs in Fig. . Currents as functions of time for the simulation runs with different voltages are shown in Fig. a. In the first minute the current rises from zero to an equilibrium value around which there are oscillations. Figure c shows the mean value for $t\mathit{>}60\mathrm{s}$ as a function of total voltage. The error bars show the root mean square value of the oscillations. The amplitude of the oscillations increases with voltage, but the mean does not. The difference in the mean current between the runs with different voltages is smaller than the current fluctuations. Thus, we find no DC current–voltage relation for the return current region in the parameter range we have investigated. Outside this range there must be a transition from this behaviour to that of the upward current region where the DC current does depend on voltage. A model of the downward current region with a constant current and a transition region was used by . When a case with zero voltage was computed with the same Vlasov model as we use here, the resulting current was found to be close to zero . Thus we can confine the transition to the $-100\mathrm{V}\mathit{<}{V}_{\mathrm{p},\mathrm{I}}\mathit{<}\mathrm{0}$ range.

For the runs with different ionospheric densities, current densities are shown as functions of time in Fig. b. Here, both the mean value of the current and the fluctuation amplitude increase with increasing ionospheric density. This is also seen in Fig. d, which shows the mean current density ($\times$) and the root mean square fluctuation (error bars) as a function of the density at the ionospheric boundary. The black line in Fig. d shows the current density corresponding to the upward flux of electrons from a $k_{\mathrm{B}}{T}_{\mathrm{e}}=50\mathrm{eV}$ Maxwellian as a function of density $n_{\mathrm{I}}$. The electron density that is available to carry a current is that just above the sheath that is formed at the boundary, as was discussed in Section . The circles in Fig. d show the mean current density as a function of the density 3.5 $\mathrm{km}$ above the boundary, which is outside the sheath. There are fluctuations also in the density, and the values that are used are formed by taking the mean at that position over $t\mathit{>}60\mathrm{s}$. The agreement between the circles and the line is reasonable, and we conclude that the mean current is determined by the electron density available above the sheath in our simulations or above the corresponding transition region in space.

Currents in simulation runs with (a) different voltages and (b) different ionospheric densities; the current as a function of (c) voltage and (d) ionospheric density. The mean and the standard deviation, shown by the error bars, in panels (c) and (d) have been computed for $t\mathit{>}60\mathrm{s}$. The densities that were used for the circles in panel (d) were taken 3.5 $\mathrm{km}$ from the ionospheric boundary, and the black line shows the current corresponding to the electron flux of a Maxwellian with $k_{\mathrm{B}}{T}_{\mathrm{e}}=50\mathrm{eV}$.

[Figure omitted. See PDF]

Discussion

We have simulated the downward current region of the aurora and found double layers moving upward. In contrast to the upward current region, where there is a stable equilibrium double layer position, the upward motion is a necessity for the existence of double layers in the downward current region, as was shown by . The density of the ionospheric electrons decreases with altitude as the flux tube becomes wider due to the decreasing magnetic field. Thus, when the double layer moves upward it eventually reaches a point where the Langmuir condition no longer can be satisfied even in the moving frame of reference. Then the double layer is disrupted and a new one forms near the ionosphere.

Both the standard and the generalised Langmuir conditions are derived for steady-state conditions. In Sect.  we applied Eq. () in a moving frame of reference. This is possible for laminar double layers, which are quasi-stationary in that frame. Beyond the gap region, phase space holes and large wave fields lead to a non-stationary situation, which is seen in Fig. . For the turbulent double layers that do not have a gap region there is no stationary state at all, and numerical application of conditions as those in Eqs. () and () becomes unsuitable. Even in the example in Fig. , where the double layer is in a particularly laminar state, pressure balance is nearly achieved over a wide range of velocities, and therefore the use of the criterion to find a numerical value for the double layer velocity cannot be accurate. When pressure balance cannot be maintained the double layer must vanish. First, it undergoes a transition from a laminar to a turbulent state. Then, it disappears altogether. The demise of the double layer enables the rise of its successor. The conditions for double layer formation close to the ionosphere are made favourable by the removal of both the potential dip on the low potential side and the ambipolar electric field in the space between the ionosphere and the double layer.

For the upward current region stationary solutions to the Vlasov–Poisson system may be found e.g., and in time-dependent Vlasov simulations there is an equilibrium large-scale solution even if small-scale waves still are present . The particular solution that is found depends on the boundary conditions, for example the double layer position in the upward current region is at lower altitudes for higher acceleration voltages . Hence, if there are fluctuations in the boundary conditions, the solution will change between different equilibrium states, and this would cause double layers to move. Fluctuations could, for example, be caused by the removal of electrons from the ionosphere in the return current region, affecting the conductivity of the ionosphere, or by changes in the generator that drives the current circuit. The Cluster spacecraft have observed fluctuations on the timescale of minutes, or possibly faster as the temporal resolution is limited by the spacecraft separation . Simulations have shown that there are hysteresis effects in the position of the double layer which may move a few hundred kilometres to a new position, if the voltage makes an excursion and returns to its original value . Thus, in the upward current region, there may be fluctuations around a stable equilibrium, or the system may be in transition between different equilibria. In contrast, there is no such equilibrium in the downward current region over the parameter range that is covered in this work.

In the upward current region a relatively simple DC current–voltage relationship exists . Observational studies have confirmed a linear current–voltage relationship . combined ground-based and spacecraft observations and found good agreement between their observed conductance and that derived by . , on the other hand, found currents several times higher than predicted by Knight's relation during periods of changing voltage. Knight's relation is derived under the assumption of a steady state, and the reason for the difference between the observations could be that this assumption was correct for the cases studied by but not for those studied by . In the downward current region there is no steady-state solution to the Vlasov–Poisson system and no DC current–voltage relationship. Instead there are large fluctuations in the current density, and its mean value is limited by the density and temperature of the ionosphere rather than the voltage over the flux tube. Nevertheless, there may be ways that are outside the model used here by which the voltage can influence the current. We include the ionosphere in our model only through the boundary conditions, and no feedback between these two elements of the auroral current circuit have been taken into account. found a correlation between the current and the parallel electric field. It is possible that a higher voltage, and hence a larger electric field, leads to an increased heating of the ionospheric plasma through better confinement of ions to a region where they are heated by waves of kinds that are not included in our model. That is the basic principle of heating in a static pressure-cooker model. The higher temperature would then allow a higher upward electron flux. However, the correlation may have other causes. What we can say is that it is not a simple current–voltage relationship as that in the upward current region.

The oscillations of the current seen in Fig.  would be expected to perturb the magnetic field, and thus induced perturbations would propagate as Alfvén waves along the field lines. Such waves cannot be simulated by our electrostatic model, but we can use Ampère's law to make a simple estimate of what the amplitude would be. If we take the current channel to have a circular cross section with a radius of $r_{\mathrm{I}}=\mathrm{5}\mathrm{km}$ at the ionosphere and use the root mean square value of the fluctuations in Fig.  for the current density $J_{\mathrm{I}}=\mathrm{1}\mathrm{\mu }\mathrm{A}/{\mathrm{m}}^{\mathrm{2}}$, we obtain an azimuthal field $B_{\mathit{\varphi }}={\mathit{\mu }}_{\mathrm{0}}{J}_{\mathrm{I}}{r}_{\mathrm{I}}/\mathrm{2}\approx \mathrm{3}\mathrm{nT}$. This is much smaller than the unperturbed field $B_z=56\mathrm{\mu }\mathrm{T}$, but it could possibly be detected by spacecraft-based magnetometers. The fundamental period of the oscillations in Fig.  is about a minute, and that is longer than it takes a spacecraft to pass through the return current region. The faster oscillations that are also seen in Fig.  would be more likely candidates for detection.

simulated a moving double layer in a Vlasov simulation, where perpendicular heating was introduced by widening of the distribution in the $\mathit{\mu }$ dimension in proportion to the energy density of the parallel electric field. They found that the double layer speed increases with increased perpendicular heating. Perpendicular kinetic energy is converted to parallel kinetic energy through the action of the magnetic mirror force. In the simulations presented in the present paper, no mechanism for perpendicular heating has been included. Instead we fixed the temperature at the ionosphere to 50 $\mathrm{eV}$ as a boundary condition. This gives realistic temperatures at the double layer position although the feedback from the waves that provide the heating is not included in the model. A higher temperature allows the plasma on the low potential side to expand faster, and that also leads to a higher double layer speed. As we have seen above, the Langmuir condition can be satisfied over a wide range of velocities. The exact value of the double layer velocity is determined by the expansion of the low potential plasma, which also depends on double layer properties such as the potential dip on its low potential side and the associated ambipolar electric field.

proposed a modified pressure-cooker model, where the static electric field is replaced by a moving double layer, and they found that the bulk of the outflowing ions do not pass through this double layer. Our results agree with this model. The analogy with a pressure cooker is not as striking as before, since the constantly moving double layer does little to confine the ions. Nevertheless, double layers play an essential role in ion outflow in the downward current region. The double layer accelerates the ions that heat the plasma on the low potential side. This allows the plasma to expand, and the ions in that plasma continue to move upward when the double layer is disrupted and a new one forms. Some of these ions are accelerated downward again by the new double layer, and thus they take part in the heating process. The rest continue to move upward, contributing to the outflow.

To make predictions on timescales of several minutes one would need to model the complete circuit, including both the ionosphere and the magnetosphere as well as the upward and downward current regions. That would increase the computational cost, and knowledge of the state of the plasma in large regions of the ionosphere and the equatorial magnetosphere would be required in order to fix the initial and boundary conditions. While we do not model the feedback from a changing ionosphere, we can draw some conclusions from the results of simulations performed with different parameters. The time between two consecutive double layers passing a specific altitude spans from about 40 to about 70 $\mathrm{s}$ in Fig. . This is a 30 $\mathrm{s}$ interval. If one spacecraft observes a double layer passage and another follows along the same path 55 $\mathrm{s}$ later – that is to say, centred on that interval – we could make an order of magnitude estimate of the probability of the second spacecraft also observing a double layer passage. Assuming that the recurrence period is uniformly distributed between 40 and 70 $\mathrm{s}$ and that it takes the spacecraft 1 $\mathrm{s}$ to pass the return current region we arrive at a probability of $\mathrm{1}/30$ for a second observation given that the first one has occurred.

The Supplement related to this article is available online at doi:10.5194/angeo-33-1331-2015-supplement.

Acknowledgements

This work was supported by the Belgian Science Policy Office through the Solar–Terrestrial Centre of Excellence and by PRODEX/Cluster contract 13127/98/NL/VJ(IC)-PEA 90316. This research was conducted using the resources of the High Performance Computing Center North (HPC2N) at Umeå University in Sweden. The topical editor L. Blomberg thanks C. Watt and one anonymous referee for help in evaluating this paper.