The code used in these simulations is described by Swift [1996]. It is based on the hybrid code originally developed by Harned [1981]. In the hybrid code, the magnetic field is advanced in time from Faraday?s law.

  (1)

  (2)

  (3)

  (4)

  (5)

  (6)

The code used here has several important modifications from the original code described by Harned. One is that the code is written for a generalized curvilinear coordinate system. The particular coordinate system used is shown in Figure 1. One of the coordinate boundaries is the Earth?s surface at 1 R_E. Two other boundaries are the polar axes, which extend out to ± 11 R_E, while the outer "radial" coordinate is extended to x = -36 R_E to include much of the plasma sheet and magnetotail lobes. This particular coordinate system has one degree of latitude resolution near the Earth and has a reasonably high density of grid points across the plasma sheet. The coordinate system is specified as a table giving the GSM [x,z] positions of the coordinate points. The ignorable, y, coordinate points from dawn to dusk.

Another important feature of the code is provision for use of the fluid approximation in regions where the plasma is not expected to have kinetic behavior. This is generally true where plasma velocities are small, such as in the ionosphere and inside the plasmapause. The plasma in this region is dense, so use of the fluid approximation results in an enormous computational saving over having to follow the motion of a large number of individual particles. The code also allows the particle and fluid descriptions to occupy the same spatial volume. This is necessary because there will be energetic particles entering regions near the Earth along the high-latitude field lines.

Another feature of the code is the provision to subcycle the update of the fields to the particle push. The hybrid code has a Courant condition with respect to propagation of the whistler or Alfvén mode, depending on whether the gyroperiod is longer or shorter than the time step. In regions near the Earth, the Alfvén and whistler mode velocities are large and the distances between grid points are small. This necessitates use of a time step much smaller than that required to follow the motion of the ions. To make the code reasonably efficient, we subcycle the field update to the particle update. In the runs to be described below, thirty time steps were used to update the fields for each particle time step. On a vector computer, computer time for one particle time step is about one hundred times the time used for a field update, so the subcycling results in an enormous computational savings.

Other features included in the code are provisions for inclusion of collisional resistivity and electron inertia. The resistivity is applied to control slowly growing numerical instabilities that tend to arise where the grid is most non-uniform and where there is a large departure from orthogonality. Electron inertia, which also adds a stabilizing effect, is included.

One more feature that was added subsequent to the publication of the Swift [1996] paper is the use of the Villasenor-Bunemann [1992] exactly conservative current algorithm to compute the ion flux and density in (5). This algorithm was originally developed for electromagnetic codes where exact charge conservation is a requirement. Exact charge conservation is not as much a requirement in a hybrid code. However the algorithm has an advantage in non-orthogonal curvilinear coordinates in that the components of the ion flux are computed in cell-face normal components centered on cell faces, exactly the form needed in (3). Since the particle velocities are represented in Cartesian components, this saves conversion from Cartesian to curvilinear components and interpolation from cell center to cell faces.

Contents

Section 4

Section 5