The penetration of ions into the magnetosphere through the magnetopause turbulent current sheet

This paper reports the results of numerical modeling of magnetosheath ion motion in the magnetopause current sheet (MCS) in the presence of magnetic fluctuations. Our model of magnetic field turbulence has a power law spectrum in the wave vector space, reaches maximum intensity in the center of MCS, and decreases towards the magnetosheath and magnetosphere boundaries. We calculated the density profile across the MCS. We also calculated the number of particles entering the magnetosphere, reflected from the magnetopause and escaping from the flanks, as a function of the fluctuation level of the turbulence and magnetic field shear parameter. All of these quantities appeared to be strongly dependent on the fluctuation level, but not on the magnetic field shear parameter. For the highest fluctuation levels the number of particles entering the magnetosphere does not exceed 15% of the total number of particles launched from the magnetosheath side of the MCS; the modeling also reproduced the effective reflection of the magnetosheath flow from very high levels of magnetic fluctuations.


Introduction
The magnetopause, a thin (∼500-1000 km) boundary (Phan and Paschmann, 1996) where the magnetic field strength and orientation change from the magnetosheath (MSH) to magnetosphere (MSP) values, was formerly considered as an impenetrable boundary, separating cold and dense MSH plasma from hot and tenuous MSP plasma.This, however, is at variance with the recent spacecraft measurements which indicate the existence of solar wind-like ion composition in magnetosphere plasmas and thus provide the evidence for plasma entry from MSH into the MSP (Sibeck et al., 2000).
Correspondence to: A. Greco (greco@fis.unical.it)Our present study is devoted to the analysis of the magnetosheath plasma transport through the magnetopause, which plays the role of a magnetic shield against the penetration of solar particles to the near-Earth environment.This is one of the challenging questions of the magnetospheric physics.
Recently, Borovsky (2002) showed that there is a strong correlation between the turbulence in the solar wind (SW) and the geomagnetic activity, so that turbulence is an important mechanism to provide the coupling of the solar wind to the Earth's magnetosphere.INTERBALL and POLAR satellite measurements (Savin et al., 1998(Savin et al., , 2001(Savin et al., , 2002) ) indicated that magnetic field turbulence is a common feature of the high-latitude magnetopause and regions close to it.The observed magnetic fluctuation energy content could reach up to 10% of the total kinetic energy of the SW plasma particles.
Traversing the magnetopause, the magnetic field turns from the IMF direction to the geomagnetic direction, so the magnetopause can be roughly described as a rotational discontinuity.We model the MCS global magnetic field geometry as a superposition of the sign reversing magnetic field component B z and the shear field component B y (see Fig. 1).On this non-fluctuating (regular) field our model superimposes three-dimensional stochastic magnetic field fluctuations that have a maximum intensity in the center of MCS and decrease towards the MSH and MSP boundaries.This magnetic turbulence has a power law spectrum in wave vector space, which corresponds to the well-developed turbulence at the stage of nonlinear saturation.There is a kink in the spectrum corresponding to the POLAR and INTERBALL satellite observations (Savin et al., 2002) in a frequency space.
Below we present the results of the test particle simulation in such a configuration.The numerical technique is described in detail by Veltri et al. (1998).Ions launched from the MSH side have a shifted Maxwellian initial velocity distribution with a flow velocity vector located in the same plane as the MSH magnetic field.In our previous studies (Taktakishvili et al., 2002) we examined a randomly uniform spatial distribution (in a statistical sense) of the magnetic fluctuations in the MCS volume.In that case the results of calculations depend on the geometry of the model configuration, e.g. the number of penetrated particles would be inversely proportional to the thickness of the magnetopause/MCS.The present analysis, with the fluctuations concentrated in the center of MCS, makes the model much more robust to changes in the spatial scales of the simulations and is more relevant from the point of view of observations, which clearly demonstrate a non-uniform distribution of the fluctuations in the turbulent boundary layer between MSH and MSP (Savin et al., 2002).In the presented analysis we varied the fluctuation level of the turbulence and the magnetic field shear parameter while keeping the other parameters (such as flow velocity and temperature; the subject of future studies) unchanged.We counted the number of particles that penetrated the MSP side of the MCS, those that are reflected back to the MSH, and those that escaped from the flanks of the MCS, and calculated the particle density profile across the MCS.
Comparing our results with those of our previous study, we find that the number of particles crossing the MCS is only slightly reduced when the magnetic turbulence profile is peaked in the center, so we could say that the main controlling parameter for ion transport is the fluctuation level.

Model geometry and particle trajectories
To describe the most principal features of the real magnetopause configuration we use a magnetic field model geometry that consists of three parts: two constant, non-fluctuating field components: 1.A one-dimensional magnetic field reversal in the zdirection, changing sign in x = 0 plane; 2. A magnetic field component in the y-direction, which causes the magnetic field shear around the x axis (we call it the "shear component") and; 3. Superimposed magnetic fluctuations, which have a three-dimensional power law spectrum.The schematic of this model geometry is plotted in Fig. 1 and is represented by the following equations: The first term in Eq. ( 1), describing a sign reversing component, is taken as the Harris (1962) current sheet magnetic field, B z (x) = B 0z (x)e z ∼ B 0 tanh(x/L x )e z , slightly modified in such a way that the asymptotic value B 0 is reached at the edges of our simulation box, x = ±0.5 L (see Fig. 1).Here, L is the total thickness of our magnetic configuration in the cross layer (x) direction (for numerical estimations we take L = 500 km), L x = 0.25 L is the MCS half thickness.The simulation box is extending up to 60 L in both the y and z directions (−30 L ≤ y = 30 L, 0 ≤ z ≤ 60 L), and it is considered to be local, in the sense that the average macroscopic features do not change within the box.
In order to accurately take into account particle Larmor gyration in the vicinity of the simulation box borders x = ±0.5L, a particle is considered completely gone (exited) from the simulation box only when it crosses, respectively, the x = 0.5L + 0.25L = 0.75L and the x = −0.5L− 0.25L = −0.75Lsurfaces (the Larmor radius is of the order of 0.05L, but there could be particles with larger Larmor radii in the tails of the injection velocity distribution function).The average magnetic field assumes the value of B 0 from the planes x = ±0.5L to X = ±0.75L.The vector sum of the B z and B y fields describes the magnetic field shear around the x axis (see Fig. 1), and the ratio b y = B y /B 0 is defined as the shear parameter.The fluctuating field given by Eq. ( 2) is the sum of the Fourier static magnetic perturbations.Here, e σ (k) are the polarization unit vectors, φ σ k are random phases, simulating the irregular nature of the magnetic fluctuations, and the amplitude of perturbations has a power law dependence on the wave number: where k x , k y and k z are chosen on a discrete grid, l i are the turbulence correlation lengths in the x, y, and z directions and α is the spectral index.
In order to reproduce the concentration of magnetic fluctuations in the current sheet region, we assume the following parity conditions for the three components of the magnetic field with respect to the x direction: δB x is even in x, i.e. δB x (x, y, z) = δB x (−x, y, z); δB y , and δB z are odd in x, i.e. δB z (x, y, z) = −δB z (−x, y, z), δB y (x, y, z) = −δB y (−x, y, z).At the same time, the phases of the fluctuating magnetic field components in Eq. ( 2) are chosen randomly for different k z , k y , but do not depend on the wave vector component in the x-direction k x , that is, if we change k x only, the random phase is the same.As a consequence, all wave modes with the same k y and k z , but different k x , have the same phase and they sum up in a coherent way when one moves along the normal direction x only, leading to the formation of coherent structures, like magnetic islands.Together with the parity rules, this provides the concentration of the perturbations in the vicinity of the central plane x = 0 (see Veltri et al. (1998) for details).The correlation lengths l x , l y and l z are taken in such a way to reproduce the features of the actual geometry of the magnetosphere: the correlation length in the direction perpendicular to the current sheet, i.e. the x direction is smaller than the correlation lengths in the z and y directions, the last two being of the same order:  l z , l y , l z ≈ l y .From the condition that at least three modes in each direction are localized in the injection zone in k-space of the considered turbulence represented by Eq. (3), i.e. between k j = 0 and k j = 1/ l j and from the restriction on the maximum value of the tearing mode wave length (k x l x < 1), we finally obtain the estimation for l x ∼ 0.05L.For l y and l z we took correspondingly l y = l z = 10l x .
The spectrum of the fluctuations in the wave number space was taken in correspondence to an observational sample in the frequency space (Savin et al., 2002) and has a power law dependence on the wave number.In Fig. 2 our model spectrum is presented schematically in logarithmic scale, with the magnetic field normalized to B 0 and the wave vector to 1/L.The spectrum length in k-space in our model is determined by the ratio between its extreme values (and is restricted by the computer resources), k(max)/k(min) = 20, k(min) corresponding to the largest wave length, which was taken to be equal to L, while the minimum wave length in the x direction, L/20, is approximately equal to the ion thermal Larmor radius in the B 0 field.The spectrum changes its slope (power index), becoming steeper for higher wave numbers and it also reproduces the INTERBALL/POLAR slopes for the constant phase velocity V ∼ ω/k.The kink in the spectrum was placed at k kink = 10 k(min).We took α = 1.5 for k < k kink and α = 2.3 for k > k kink .We assume that our source test particles are coming from the magnetosheath represented by the upper region in Fig. 1 above x = 0.5 L plane.This plasma is carried by the solar wind magnetic field (according the frozen law condition) and when it meets the obstacle of the Earth's magnetic field, the plasma flows around the magnetosheath, that is the plane xy in Fig. 2 with a streaming velocity in the same plane.
The streaming velocity vector was taken in the present calculations as a fixed parameter, placed in the x = 0.5 L plane, having a value of u ≈ 100 km/s, and making a fixed angle θ (taken to be equal to θ = π/4 in the present calculations) with the magnetic field in this plane.The influence of the variation of the streaming velocity, and hence, the variation of the ratio Larmor radius over the correlation length in the normal direction x, l x , will be reported later.The relatively cold ion distribution in the magnetosheath can be described as a shifted Maxwellian: Here ν th is the thermal velocity, and the temperature was estimated as T = m i ν 2 th ≈ 100 eV (Le et al., 1996;Song et al., 1990Song et al., , 1993)).Particles are injected with random velocity components in such a way as to reproduce the distribution function given by Eq. ( 4).The details of the injection scheme in velocity space are given in Appendix A. The level of fluctuations is characterized by the parameter δB/B 0 , where δB is the root mean square of the fluctuation field given by Eq. ( 2), averaged over the whole simulation box.
Without fluctuations, δB/B 0 = 0, the constant magnetic field B = B 2 0 + B 2 y of the upper plane x = 0.5 L acts as an impenetrable boundary for the ions, since their thermal Larmor radius is significantly smaller than the cross layer scale L. The situation changes substantially in the presence of the fluctuations δB/B 0 = 0: some particles, helped by the complicated perturbed field, manage to traverse the MCS and cross the lower plane x = −0.5 L, thus penetrating into the magnetosphere.We call them penetrating particles.The sample of such a trajectory in xz (upper panel) and xy (lower panel) planes is presented in Fig. 3.A larger number of particles is not crossing the current sheet while being scattered by fluctuations, going back to the magnetosheath, and we call them reflected particles (Fig. 5).There are also particles that enter into the MCS, perform some stochastic motion and then exit from the sides of the simulation box in the y and z directions (flanks), because their orbits are still strongly controlled by the guiding, non-fluctuating sheared magnetic field.We call them escaping particles (Fig 4).
Figure 3 displays the trajectory of a penetrating particle performing a complicated random and bouncing motion in the MCS.It covers a rather large distance both in the z and y directions before crossing the MCS in the x direction and reaching the surface x = −0.5 L. Note that the trajectory becomes more stochastic close to the central plane x = 0, where the fluctuations are stronger while it demonstrates relatively "regular" behavior, typical for a magnetized particle, at the edges of the MCS, where the influence of the nonfluctuating magnetic field on the particle dynamics is prevailing.
The trajectory of an escaping particle shown in Fig. 4 is much more "regular" than the trajectory of the penetrating particle, while the distance it traverses in the z and y directions is considerably less.Note that the distance that the particle traversed in the x direction is also small, so that the motion occurs close to the outer edge of the MCS, and that is why it is more "regular".After some bouncing, mostly in the z direction, this particle crosses the flank-border of the simulation box, which, in the displayed case, corresponds to z = 0 plane.Finally, from the sample of a reflected particle displayed in Fig. 5, we notice that it enters the simulation box from the outer edge and its trajectory consists of Larmor gyration and B × ∇B gradient drift motion, with velocity  components both in the z and y directions.After few gyrations and short traversed distances in all directions, it leaves the simulation box from the same injection plane.

Simulation results
Figures 6-8 show the number of penetrating, escaping and reflected particles respectively, as a percentage of the total number of the launched particles, as a function of the fluctuation level δB/B 0 , for two different values of the shear parameter, b y = 0.5 and 1.0.
Without any fluctuations, δB/B 0 = 0, or when the fluctuations are relatively small, a constant, non-fluctuating magnetic field (or regular field as we call it), which is parallel to the yz plane, is an obstacle for the particles, by not allowing them to move in a cross layer (x) direction.This is clearly seen in Fig. 6, which shows almost no penetration of ions (only few tenths of percent) up to δB/B 0 = 0.3 .When the fluctuations become relatively large, the destruction of magnetic surfaces provides the particles with mobility across the average position of the distorted field line.This mobility is driven by: 1.The finite ion Larmor radius effect, allowing an ion to jump from one magnetic surface to another and eventually to be displaced in the cross layer (x) direction; 2. The perturbed component of the magnetic field in the x-direction, which bends the magnetic field lines, allowing particles, even with zero Larmor radii, to have mobility in the cross layer direction.
The larger the fluctuation level the more the distortion is of the magnetic field lines and the more the particles are able to penetrate and reach the opposite side of the simulation box (x = −0.5 L).However, this scenario only works up to a certain fluctuation level.The number of penetrating particles grows rapidly with the increase of the fluctuation level up to δB/B 0 ≤ 1, (see Fig. 6), but then, for strongly nonlinear fluctuations δB/B 0 > 1, the growth in the number of penetrating particles slows down and may even be reduced, as is demonstrated for the curve b y = 1 in the interval 2 < δB/B 0 < 2.5.This is because very high fluctuations strongly scatter particles in all directions and thus are an obstacle for the effective particle penetration in the cross layer direction as well.From this figure we see also that for a given fluctuation level the number of penetrating particles is slightly higher for smaller b y , due to the fact that with the increase in the regular (not fluctuating) magnetic field component, the effect of the finite Larmor radius, one of the drivers of the penetration process through the MCS, becomes weaker.The maximum penetration rate, obtained in our calculations, does not exceed 15% of the total number of particles.
In Fig. 7 the number of escaping ions is plotted: it decreases monotonically with the growth of the fluctuation level, being almost 95% for the smallest fluctuations and then rapidly going down for larger δB/B 0 .The reason for this is the following: it is a non-fluctuating sheared magnetic field that guides particles towards the flanks both in the y and z directions, and when the magnetic field becomes strongly distorted, its guiding role essentially weakens.Clearly, for the same reason, the escaped particle number is larger for larger b y .
The behavior of the number of reflected particles is shown in Fig. 8: it grows monotonically with the increase in the fluctuation level and for the same level of perturbations it is smaller for larger b y , because the growth of the regular magnetic field leads to the increase in the share of particles escaping from the flanks, as has been discussed above.But this difference vanishes for the smallest and largest fluctuations, which is due to the interplay of different influences of the scattering effect of the fluctuating field and the regular field on particle dynamics.For highest fluctuations the number of reflected particles reaches the level of 80%.It is as if a very dynamic magnetopause is "frying" and most particles are repelled off.Then, our modeling predicts an effective reflection of particles coming from the magnetosheath for high turbulence levels in connection with the high-latitude magnetopause and cusp region, where the high-level of turbulence is nearly permanently observed (S.Savin, private communications).
In Fig. 9 we plotted the percentage for all types of trajectories for the same value of the shear parameter.It is clear that the curves representing escaping and reflected particles look like mirror reflections of each other.The above mentioned results are confirmed by the calculations of the density profile of ions in the cross layer direction, plotted in Fig. 10 (here, density is normalized to its value at the border x = 0.5 L) for different fluctuation levels: density drops rapidly across the magnetopause, in agreement with the observations of lowshear and high-shear magnetopause (Paschmann et al., 1993;Phan et al., 1994;Song et al., 1990Song et al., , 1993)).Still, density at x = −0.5 L is larger for larger δB/B 0 , showing that turbulence allows an effective plasma flux across the magnetopause.
Fig. 10.The cross layer profile of normalized ion density for different δB/B 0 .The lowest solid curve corresponds to the smallest δB/B 0 = 0.1, the dotted curve to 0.2, the short-dashed one to 0.3, the long-dashed one to 0.5, the dot-dashed to 0.7 and the upper curve to the largest δB/B 0 = 0.9.

Conclusions
We performed a test particle simulation of ion motion in the turbulent magnetopause, modeling it as a sheared magnetic field configuration with additional stochastic static magnetic field fluctuations.By varying the level of the fluctuations and the magnetic field shear parameter we calculated the number of particles penetrating from the magneosheath to the magnetosphere side of the magnetopause.We found that the number of penetrating particles grows rapidly with the growth of the fluctuation level in the "linear" regime, δB/B 0 ≤ 1, while for highly nonlinear fluctuations, δB/B 0 > 1, the number of penetrating particles slows down or even sometimes is reduced.The penetrating particle number is slightly larger for smaller shear parameter b y values, but this dependence does not seem to be strong (though the latter statement has to be checked additionally, since so far we used only two values for the shear parameter).The maximum number of penetrated particles is ∼15% of the total number of particles injected from the magnetosheath.
We also counted the number of particles reflected back to the magnetosheath and escaping from the flanks in the direction of the sheared magnetic field component.These two kinds of ion motion exhibit opposite behavior: with the growth of the fluctuation level the number of reflected ions increases while the number of escaping ions decreases, so that they are almost mirror image of each other.The dependence on the b y parameter is not strong in this case either.
The analysis of ion density profile in the cross layer direction revealed its rapid drop when going in the direction from the magnetosheath towards the magnetosphere, which is in agreement with the observations.The higher the fluctuation level, the higher the density and the flatter the density gradient became in the cross layer direction.In this study we did not change the parameters other than the shear field value b y and the fluctuation level δB/B 0 .
Our results could slightly depend on the ratio between the correlation length in the x direction and the Larmor radius; in order to see how the turbulence effect can be changed for smaller or larger ratios, we have varied the ion injection velocity u, which, in turn, means varying the Larmor radius.The results, shown in Fig. 11 and Fig. 12, are almost the same, with a slight increase in the penetrated particles and a slight decrease in the reflected ions for increased energy.Still, the dimensions of the simulation box and the thickness of the current sheet could affect our results; in order to check this dependence, we have tried to inject the particles from a plane closer to the current sheet, and we have obtained percentage values slightly different from the case described above.We have also changed the profile of the turbulence, distributing it uniformly in the entire simulation box, and in this case the results remain quite unchanged, showing that the main new ingredient that allows for particles to penetrate the magnetic field in the cross layer direction is the presence of magnetic turbulence.
Finally, we would like to make some remarks on the applicability of our results to the actual situation in the magnetopause.Of course, static magnetic fluctuations only approximately describe real magnetic turbulence, which is always time dependent.The non-stationary character of the magnetic fluctuations could be neglected in the case when particle thermal velocity is larger than Alfvén velocity, and thus, the influence of induced electric fluctuations on the particle dynamics may be ignored (this was discussed in more detail by Veltri et al., 1998).But unfortunately, this is not always the case for cold magnetosheath plasma.However, we think that, when the phases are random, qualitatively, the physics does not change much, whether the turbulence is time dependent or not, because in any case, particles "see" different values of the magnetic field.
At least by neglecting the time dependence and inductive electric fields, we, of course, underestimate the rate of diffusion in our analysis.Indeed, by including the time dependent field, the rate of diffusion would only be enhanced by electric fluctuations (additional random force acting on particles).In the same way we do not take into account slow diffusive energization of ions in the course of their penetration through the magnetopause in the cross layer direction, though this could be a small effect, because a particle spends a relatively short amount of time in the turbulent layer.

Fig. 1 .
Fig. 1.The schematic of the magnetopause sheared magnetic field geometry (fluctuating field not shown).

Fig. 2 .
Fig. 2. The model spectrum of the magnetic fluctuations (axis have logarithmic scale).Magnetic field is normalized to B 0 , k to 1/L.

Fig. 3 .
Fig. 3.A sample of a penetrating ion trajectory in the simulation box.

Fig. 6 .
Fig. 6.The number of penetrating ions as a percentage of the total number of ions launched, as a function of the fluctuation level for two values of the shear parameter b y 0.5 (squares) and 1.0 (triangles).

Fig. 9 .
Fig. 9.The same as in Figs.6-8 for all types of trajectories and b y = 1.

Fig. 11 .
Fig. 11.The number of penetrating ions as a percentage of the total number of ions launched, as a function of the fluctuation level for one value of the shear parameter b y = 1.0 and for three values of the injection velocity u.

Fig. 12 .
Fig. 12.The number of reflected ions as a percentage of the total number of ions launched, as a function of the fluctuation level for one value of the shear parameter b y = 1.0 and for three values of the injection velocity u.