Plasma depletion layer: Magnetosheath flow structure and forces

The plasma depletion layer (PDL) is a layer on the sunward side of the magnetopause with lower plasma density and higher magnetic field compared to the corresponding upstream magnetosheath values. In a previous study, we have validated the UCLA global magnetohydrodynamic (MHD) model in studying the formation of the PDL by comparing model results, using spacecraft solar wind observations as the driver, with in situ PDL observations. In this study, we extend our previous work and examine the detailed MHD forces responsible for the PDL formation. We argue that MHD models, instead of gasdynamic models, should be used to study the PDL, because gasdynamic models cannot produce the PDL on the sunward side of the magnetopause. For northward interplanetary magnetic field (IMF), flux tube depletion occurs in almost all the subsolar magnetosheath. However, the streamlines closest to the magnetopause and the stagnation line show the greatest depletion. The relative strength of the various MHD forces changes along these streamlines. Forces along a flux tube at different stages of its depletion in the magnetosheath are analyzed. We find that a strong plasma pressure gradient force along the magnetic field at the bow shock and a pressure gradient force along the flux tube within the magnetosheath usually exist pushing plasma away from the equatorial plane to deplete the flux tube. More complex force structures along the flux tube are found close to the magnetopause. This new, more detailed description of flux tube depletion is compared with the results of Zwan and Wolf (1976) and differences are found. Near the magnetopause, the pressure gradient force along the flux tube either drives plasma away from the equatorial plane or pushes plasma toward the equatorial plane. As a result, a slow mode structure is seen along the flux tube which might be responsible for the observed two-layered slow mode structures.


Introduction
The plasma depletion layer (PDL) is a distinct structure on the sunward side of the magnetopause that differentiates itself from the usual magnetosheath structure by its increased magnetic field and decreased plasma density. The PDL usually exists when the interplanetary magnetic field (IMF) is northward with low magnetic shear across the magnetopause, in which case the magnetic reconnection is weak so that magnetic field lines can pile up on the sunward side of the magnetopause. However, the PDL has also been observed under southward IMF when the reconnection process is insufficient to relieve the pile-up of magnetic field lines (e.g. Anderson et al., 1997).
Observations of the plasma depletion near the magnetopause were first reported by Cummings and Coleman (1968), with many follow-on studies of the stationary PDL (e.g. Crooker et al., 1979;Fuselier et al., 1991;Anderson and Fuselier, 1993;Paschmann et al., 1993;Song et al., 1993;Anderson et al., 1994;Phan et al., 1994Phan et al., , 1997. The PDL is essentially a slow-mode structure in magnetohydrodynamic (MHD) terms because the magnetic field magnitude is anti-correlated with the plasma density. By analyzing magnetosheath wave structures near the magnetopause, Song et al. (1990) also found slow mode structures with decreased magnetic field and increased plasma density upstream of the PDL. Such slow mode structures, together with the PDL, are usually called two-layered slow mode structures. Lees (1964) presented one of the earliest theoretical studies of the magnetic field effects in the magnetosheath. By studying a very simplified axis-symmetric flow configuration, he predicted a density minimum at the stagnation point. The first comprehensive theoretical plasma depletion layer study was attempted by Zwan and Wolf (1976), who investigated the evolution of a single thin flux tube moving from the solar wind into the magnetosheath. They proposed two processes for depleting the flux tube, as shown in Fig. 1. First, by deflecting plasma around the magnetosphere, the bow shock pushes plasma out along the field lines away from the nose of the magnetosphere. Second, the compressional  Zwan and Wolf (1976). Idealized sketches of the squeezing process at successive times T 1 , T 2 , and T 3 . As plasma moves around the magnetosphere, flux tubes become draped around the nose of the magnetosphere. stress exerted on magnetosheath flux tubes near the nose tends to squeeze plasma out along the field lines, further depleting them. However, there are intrinsic discrepancies in Zwan and Wolf's model. Although they pointed out that their solutions correspond to the slow mode waves expanding along the flux tube, the function of the wave in moving plasma is contrary to the slow mode wave properties, i.e. moving plasma into magnetopause nose instead of moving it out, as pointed out by Southwood and Kivelson (1992).
Observations generally have insufficient spatial and temporal coverage of the PDL for the study of its dynamics and global geometry. Observations from single spacecraft also suffer from temporal-spatial ambiguity which makes it difficult to separate a local PDL structure from temporal solar wind variations. Theoretical PDL studies also have the difficulty of treating the PDL which usually involves nonlinear processes, has complex geometry, and dynamically evolves under fluctuating solar wind conditions. Numerical simulations usually do not suffer from these difficulties and they have proven to be a powerful tool to investigate the PDL. The first numerical MHD simulations of the plasma depletion layer were conducted by Wu (1992), who simulated the large-scale three-dimensional MHD flow in the idealized magnetosheath with a perfectly conducting sphere as the magnetopause. Wu found that currents extend far into the magnetosheath and that a smooth PDL is formed on the magnetopause. Lyon (1994) used MHD simulations with much less numerical dissipation than Wu's simulations, but with similar boundary conditions, to study the magnetosheath pattern and the PDL. He found that the slow mode waves can only exist for low solar wind Mach number cases and no evidence of the slow mode waves was obtained for normal solar wind cases. The latter is inconsistent with some magnetosheath slow mode wave observations during typical solar wind conditions (Song et al., 1990(Song et al., , 1992). Lyon was unable to determine whether the observed magnetosheath density enhancement is an ideal MHD phenomenon or not. Denton and Lyon (2000) studied the effects of pressure anisotropy on the magnetosheath structure using a twodimensional MHD simulation with anisotropic pressure and by assuming a flux surface magnetopause. They found that the exact form of the parallel pressure gradient force may not be crucial for the global dynamics of the PDL. The anisotropy leads to a larger bow shock standoff distance compared to the isotropic case, due to the difference in perpendicular pressure. Their results also imply that the effects of pressure anisotropy may be even less for a three-dimensional system than for a two-dimensional system. Siscoe et al. (2002) summarized four important MHD effects in the magnetopause boundary layer and the magnetosheath that cannot be produced by gasdynamic models, one of which is the PDL. They found a clear dependence of the PDL thickness on the IMF clock angle based on their ISM global model results.
Most recently, Wang et al. (2003) compared the UCLA global MHD model results, using IMP 8 and ACE solar wind plasma and IMF observations as the driver, with Wind in situ observations for two PDL events on 12 January 1996 and 1 January 1999. They found good visual consistency between PDL observations and their model results. Meanwhile, the average model departures were found to be generally smaller than the standard deviations of observations. Wang et al. concluded that an MHD description can produce the plasma depletion layer. Furthermore, they also showed that the PDL is stable for stable solar wind conditions and that small magnetosheath variations observed by single spacecraft were temporal rather than spatial variations during these two events.
Until now global simulations have not yet been used to address the detailed formation and underlying physics of the plasma depletion layer. Although the Southwood and Kivelson (1995) model gives a phenomenological description of the structure and the underlying physics for the formation of the PDL, the validity of their model depends on the propagation of the slow mode waves in the complex flow and field geometry of the magnetosheath. The question whether their proposed slow mode front exists in the magnetosheath and, if so, whether it can significantly change the property of plasma and field across it, remains open. Other processes, for example, a slow mode expansion fan, may exist and compete with a slow mode front or shock.
In this paper we use global MHD simulations to address the fundamental physics that is relevant to the formation of the plasma depletion layer. First, we will briefly introduce MHD theory with emphasis on the MHD forces that are responsible for the PDL formation. Then, we will give a brief introduction to the UCLA global MHD model used in this study, as well as the important parameters and some specific settings for our model runs. After that we will show and discuss the model results with emphasis on the forces that govern the magnetosheath flow and the depletion of flux tubes. Finally, we will summarize our results.

Theory
Spreiter's model is one of the earliest models to describe the flow pattern in the magnetosheath (Spreiter et al., 1966). In this model, the Earth's magnetosphere is treated as a blunt body and the gasdynamic equations are numerically solved for the bow shock and the magnetosheath plasma conditions given a set of solar wind conditions. However, it is a gasdynamic model and the only force controlling the plasma motion in the magnetosheath is the plasma pressure gradient force. Although the magnetic field can be obtained from Spreiter's model (Alksne, 1967), it is done in a kinematic way and no magnetic forces are considered in the calculation. In order to stop the plasma motion toward the subsolar magnetopause, there must be a pressure gradient force pointing away from the sunward magnetopause. This usually corresponds to an enhancement of the plasma density toward the magnetopause. Unless the plasma β (the ratio between the plasma thermal pressure and the magnetic pressure) is much larger than unity in the magnetosheath, a gasdynamic model will not give the correct pressure and density along the stagnation streamline. However, when β value is very large and the magnetic force is relatively weak, the behavior of the magnetosheath plasma should be similar to the gasdynamic results, maybe with the exception of the region very close to the subsolar point on the magnetopause.
In ideal MHD, the plasma motion is described by the momentum equation: where du/dt is the rate of change of a plasma parcel's velocity along its flow path (Lagrangian derivative). There are two forces controlling the motion of the plasma parcel: the pressure gradient force, −∇p, and the magnetic force, j ×B. The pressure gradient force is directly related to the plasma density and temperature, and there is no constraint for its orientation. However, the magnetic force is always perpendicu-lar to the local magnetic field. This implies that the acceleration/deceleration of a plasma parcel along a magnetic field line can only be caused by the pressure gradient force. For a given particular field line in the magnetosheath, since different plasma parcels along this field line take different paths to their current locations from the solar wind, the plasma pressures at those points will, in general, be different. Thus, pressure gradients along this magnetosheath field line will form, which drive the plasma flow along the field line. We will study the detailed forces along magnetosheath magnetic field lines in Sect. 4.3. Due to the pressure gradients along magnetic field lines no one-dimensional model can account for the plasma depletion layer formation. For two-dimensional models, let's first assume that the flow is in the x direction, and the magnetopause extends infinitely in the z direction. Then twodimensional models can principally either have a magnetic field perpendicular to the xz plane or along the z direction. In the former case, field lines cannot convect around the obstacle which leads to an unconstrained flux pile-up that is not realistic. Assuming a finite normal velocity component through the magnetopause might lead to a steady solution. However, the dependence on such a parameter still would lead to an unrealistic model. In the later case, where the field is along the z direction, no PDL can develop because it is similar to the gasdynamic case by replacing p in the gasdynamic case with p+B 2 /2µ 0 in the MHD case. The case of a magnetic field having an arbitrary angle also leads to flux pile-up. Thus, we conclude that the PDL formation is fundamentally a three-dimensional process that not only requires a three-dimensional model, but also the analysis of gradients and forces in all three coordinate directions.

Model
The UCLA/NOAA global geospace model is used in this study, which has also been used by Wang et al. (2003) to validate the model in studying the plasma depletion layer. This model solves resistive MHD equations with isotropic pressure in a large three-dimensional volume surrounding the Earth, such that the entire interaction region between the solar wind and the magnetosphere, as well as the ionosphere, is included (Raeder, 1999). Usually, the physical domain of the model ranges from ∼−300 to ∼20 R E in the GSE x direction, and about ∼−40 to ∼40 R E in the GSE y and z directions, which comprises the bow shock, the magnetopause, and the long magnetotail. A stretched Cartesian grid is used in this model. There are usually ∼300 grid cells in the GSE x direction and ∼80 grid cells in the GSE y and z directions. In this study, we have ∼ 22 grid cells between the bow shock and the magnetopause along the Sun-Earth line, corresponding to a spatial resolution of ∼0.15 R E in the subsolar magnetosheath region. This is sufficient to resolve the PDL with a typical observed thickness of ∼0.3−1.0 R E . Another control simulation with a spatial resolution of ∼0.1 R E in the subsolar magnetosheath region has also been performed and has shown very little difference. Because of the computational requirements of the large-scale global magnetosphere simulations, the global model is parallelized and runs on supercomputers (e.g. IBM/SP2, CRAY/T3E) or local Beowulf clusters. The Message Passing Interface (MPI) library is used to parallelize the code. The scalability of the model on parallel computers has been tested with excellent results . In particular, the model can be run in real time with 60 nodes and a million grid cells.
The model starts with a cold (∼5000 K) plasma with a uniform density of 0.1 cm −3 . The magnetic field starts from a mirror dipole configuration. It takes ∼0.5−1.0 h for the whole model to reconfigure itself from this initial condition into a realistic magnetosphere configuration. After that period the model results can be used to address the physics of the dayside magnetosheath and the PDL. The model results shown in this paper are obtained three hours after the start of the model run, when stable magnetosheath structures have developed.
Spacecraft observations or idealized solar wind conditions can be used on the sunward model boundary. Open (free-flow) conditions are used on all other model outer boundaries. In this study, we use idealized constant solar wind plasma and IMF conditions to avoid the complexities that arise from temporal solar wind variations. The solar wind plasma and IMF parameters used in this study are (in the GSE coordinate): V =(−450, 0, 0) km/s, B=(0, 0, 7) nT, N=6 cm −3 , and T i =T e =10 eV. These are typical values for the solar wind, except that the IMF is due northward. The inner boundary of the global model couples to the NOAA Coupled Thermosphere Ionosphere Model (CTIM) that handles the coupling between the magnetosphere and the ionosphere (Fuller-Rowell et al., 1996;Raeder et al., 2001). In this study, however, we only use a uniform ionosphere conductivity of 5 Siemens because the ionosphere is not likely to play a significant role for the PDL.
The MHD module of the model integrates the MHD conservative equations (mass, momentum, energy) in time using flux-limited spatial difference and a predictor-corrector time stepping scheme (Raeder et al., 2001). For the time integration of magnetic field the constrained transport method is used to ensure ∇·B=0 (Evans and Hawley, 1988). There is no dipole tilt in the model run presented in this paper. In the model, anomalous resistivity is set to be dependent on the current density when the current density is larger than a threshold value. This anomalous resistivity is used to account for the magnetic reconnections in the global magnetosphere system. However, it is not important for this study because there is no reconnection near the subsolar magnetopause for northward IMF. Reconnection occurs in the polar region of the magnetosphere for this situation. However, it is not likely to have major effects on the magnetosheath near the equatorial plane where the main physics of the PDL occurs.

Basic magnetosheath pattern
The plasma density and the magnetic field magnitude in the z=0 plane are shown in Fig. 2. The open-closed magnetic field boundary is shown as a red zigzag curve and streamlines are shown as white smooth curves in each of the panels. Plasma flow coming from the solar wind is first deflected at the bow shock, then it moves around the magnetopause to the downstream magnetosheath. The plasma density and the magnetic field magnitude in the y=0 plane are shown in Fig. 3 in the same manner as in Fig. 2, except that magnetic field lines are shown as white smooth curves. The magnetic field magnitude increases from the bow shock toward the magnetopause along the streamline closest to the stagnation line in the right panel of Fig. 2. This implies that magnetic field lines are piling up on the magnetopause. Afterward, these field lines are draped around the magnetopause with plasma flow. In both Fig. 2 and Fig. 3, a clear PDL feature with decreased plasma density and enhanced magnetic field is seen near the subsolar magnetopause.
Several radially outward straight lines are drawn from the center of the Earth in Fig. 2, which are marked 0-4. Some plasma and field parameters along these lines are shown in Fig. 4. From top to bottom in Fig. 4 are the flow velocity, the magnetic field magnitude, the plasma density, and the ratio between the plasma density and the magnetic field magnitude (N/B). The magnetosphere is on the left side of the figure, the solar wind is on the right side of the figure, and the magnetosheath is between the magnetosphere and the solar wind. The magnetopause is defined as the boundary where the magnetosheath flow velocity is close to zero. The magnetopause locations on lines 0-4 are shown as black dots. In panel (c) of Fig. 4, plasma depletion is seen clearly near the magnetopause in the magnetosheath on lines 0-3, which implies that the PDL extends longitudinally along the magnetopause. However, the patterns are different on each of these lines, which implies the longitude dependence of the PDL. Specifically, the closer a radial straight line is to the Sun-Earth line, the thinner the PDL is. Especially clear on line 0, there are two regions with different density decreasing trends in the magnetosheath separated at ∼10.25 R E : a region of strong density decrease on the magnetosphere side and a region of weak density decrease on the sunward side. Both of these regions show the PDL feature with plasma density decrease and magnetic field increase, but we only call the first region the PDL because of its strong plasma depletion.
The N/B ratio shows the degree of plasma depletion in a flux tube. The smaller the N/B ratio is, the more depleted the flux tube is. The N/B ratio in the bottom panel of Fig. 4 shows that flux tube depletion occurs in the whole subsolar magnetosheath. However, stronger flux tube depletion usually exists closer to the magnetopause along lines 0-3. Similar results are obtained for the straight lines in the y=0 plane of Fig. 3, which are shown in Fig. 5. From Fig. 5 we also see that the PDL extends with latitude along the magnetopause and its thickness also depends on the latitude on the magnetopause. These results are consistent with the results of our previous PDL event study (Wang et al., 2003).
There are two forces controlling the plasma motion in the magnetosheath: the pressure gradient force and the magnetic force. The latter depends on the current density in the magnetosheath. Figure 6 shows the current density in the z=0 plane (the left panel) and the y=0 plane (the right panel). The openclosed magnetic field boundary is shown as a white zigzag curve in each of the planes. Although the open-closed magnetic field boundary is close to the peak of the current density on the magnetopause, they do not coincide. The current extends well into the magnetosheath, which is consistent with the results of Wu (1992). Together with the magnetic field in the magnetosheath, the magnetosheath current produces a magnetic force which plays an important role in shaping the flow pattern in the magnetosheath, as we will show next.

Forces
In the MHD regime, plasma motion is controlled by the pressure gradient force and the magnetic force. Figure 7 shows the pressure gradient force, the magnetic force, and their combined force in the z=0 plane (the left panels) and the y=0 plane (the right panels). The open-closed magnetic field boundary is shown as a red zigzag curve and the plasma density is shown color coded as the background in each panel. The streamlines are shown as white and black smooth curves in the left and the right panels, respectively. Field lines are shown as white smooth curves in the right panels. The arrows in each panel show the force vectors. The force arrows are properly scaled within each panel and among different panels. Note that the streamlines are only drawn on a plane, so they may not represent the real three-dimensional streamlines. Contrary to Spreiter's model magnetosheath results (Spreiter and Alksne, 1968), in which the pressure gradient force causes plasma deceleration toward the subsolar magnetopause, the pressure gradient force in the MHD simulation is directed toward the magnetopause in most of the subsolar magnetosheath (see the top panels of Fig. 7). The  Fig. 2. From top to bottom are the flow velocity, the magnetic field magnitude, the plasma density, and the ratio between the plasma density and the magnetic field magnitude (N/B). The magnetosphere is on the left side of the figure, the solar wind is on the right side of the figure, and the magnetosheath is between the magnetosphere and the solar wind. The magnetopause is defined as the boundary where the magnetosheath flow velocity is close to zero. The magnetopause locations on lines 0-4 are shown as black dots. magnetic force is the dominant force that decelerates the magnetosheath plasma as it approaches the magnetopause (see the middle panels of Fig. 7). The combined pressure gradient force and magnetic force first decelerate plasma near the subsolar magnetopause, and subsequently bend and accelerate the plasma flow around the magnetopause toward the downstream magnetosheath (see the bottom panels of Fig. 7). Note here that the bending of the streamlines does not occur in a layer close to the magnetopause in both the z=0 and y=0 planes, but occurs closer to the magnetopause near the Sun-Earth line and moves away from the magnetopause farther away from the Sun-Earth line.
To better understand plasma acceleration and deceleration in the magnetosheath, we show in Fig. 8 the combined pressure gradient force and magnetic force along the flow direction, (−∇p+j ×B) ||v , as the color-coded background in the z=0 (left) and y=0 (right) planes. The open-closed mag-netic field boundary is shown as a black zigzag curve and the streamlines are drawn with black or white thick smooth curves marked 0-3 in each of the panels. A positive background value means that the combined force is along the flow direction, thus flow is being accelerated. A negative background value means that the combined force is opposite to the flow direction, i.e. flow is being decelerated. In the z=0 plane a large deceleration region exists along the bow shock which causes the velocity drop across the shock. Inside the magnetosheath, deceleration is constrained to a region just downstream of the bow shock and close to the Sun-Earth line. Outside of this region in the magnetosheath, acceleration dominates throughout the magnetosheath but usually with a larger acceleration force closer to the open-closed magnetic field boundary. In the y=0 plane a similar large deceleration region exists along the bow shock. Inside the magnetosheath deceleration is also constrained in a region just downstream  The pressure gradient force, the magnetic force, and their combined force in the equatorial or z=0 plane (the left panels) and the noon-midnight or y=0 plane (the right panels). The open-closed magnetic field boundary is shown as a red zigzag curve, and the plasma density is shown color coded as the background in each panel. The streamlines are shown with white and black smooth curves in the left and the right panels, respectively. Field lines are shown with white smooth curves in the right panels. The arrows in each panel show the force vectors. The arrows are properly scaled within each panel and among different panels. Note that the streamlines are only drawn on a plane, so they may not represent the real three-dimensional streamlines. of the bow shock and close to the Sun-Earth line. Outside of this region in the magnetosheath the flow is accelerated. In contrast to the case in the z=0 plane, the peak of the acceleration force in the magnetosheath is smaller and the acceleration force is distributed more or less evenly for different distances from the magnetopause.
For streamlines close to the Sun-Earth line, like streamline 0 in the left panel of Fig. 8, plasma first experiences deceleration after passing through the bow shock, then acceleration. For streamlines farther away from the Sun-Earth line, like streamlines 1-3 in the left panel of Fig. 8, plasma is accelerated immediately after passing through the bow shock. This acceleration/deacceleration pattern differentiates streamline 0 from streamlines 1-3. Figure 9 shows the flow veloc-ity, the magnetic field magnitude, the plasma density, and the N/B ratio along the four marked streamlines in the left panel of Fig. 8. The horizontal axis is the distance along a streamline from its starting point in the solar wind. The PDL, i.e., decreased plasma density and increased magnetic field toward the magnetopause, is only seen on streamline 0 from ∼6 R E to ∼9 R E . On the other streamlines adiabatic fast mode plasma expansion dominates and no distinct PDL structure is seen. Next we will only concentrate on streamline 0 to investigate the forces for the formation of the PDL. field magnitude, the N/B ratio, the pressure gradient force parallel to magnetic field, the forces perpendicular to magnetic field, the forces perpendicular to the flow velocity, and the forces parallel to the flow velocity. The unit of the forces in panels (e)-(h) is 10 −17 Pascal/m. The horizontal axis is the distance along the streamline from its starting point in the solar wind. Panels (f)-(g) are drawn with the same scale, but panels (e) and (h) are drawn with a smaller scale and thus the values are more noisy. Panel (e) shows that the pressure gradient force along B is negligible. This is because streamline 0 is on the equatorial plane where the plasma pressure forms a local maximum along the magnetic field. However, the pressure gradient force along B is usually not zero out of the equatorial plane, which we will show later in the paper. The regions bounded by dashed lines, I, II, III, IV, and V, are the regions on streamline 0 with different force features. Regions I, II, III, and IV are marked in the left panel of Fig. 8, while region V is further downstream along the streamline and is not shown in that figure.
In region I, the flow velocity decreases and the magnetic field magnitude increases. This is a typical feature for a field line pile-up process due to compression. However, the plasma density does not increase in unison with the magnetic field, as one would expect for adiabatic compression. On the contrary, the plasma density first basically remains constant, then decreases. This implies that there must be forces moving plasma along the magnetic field line. In this region, the pressure gradient force decelerates plasma and the magnetic force accelerates plasma along the streamline, as shown in panel (h). The net result of these two forces is the plasma deceleration, which decreases in this region until it is close to zero at the end of this region.
In region II, although the magnetic force is still trying to decelerate plasma, the pressure gradient force begins to dominate. Plasma acceleration begins to produce a net positive force along the flow direction, as can be seen in panel (h). Much stronger forces perpendicular to the flow velocity direction and the magnetic field direction exist in this region. These forces are responsible for bending the flow around the magnetopause and draping magnetic field lines around the magnetopause into the downstream magnetosheath. At the end of region II, the bending forces reach their peak values.
In region III, continuous plasma acceleration exists because of the dominant pressure gradient force acceleration in the first half of this region and the dominant magnetic force acceleration in the second half of this region. Within this region, the magnetic field stops piling up and begins to decrease.  Fig. 9. The flow velocity, the magnetic field magnitude, the plasma density, and the N/B ratio along streamlines 0-3 in the left panel of Fig. 8. The x axis is the distance along a streamline from its starting point in the solar wind. The PDL feature, decreased plasma density and enhanced magnetic field is only seen clearly on streamline 0, which is closest to the magnetopause, from ∼6 R E to ∼9 R E .

Bow Shock
In region IV, forces are usually much smaller compared to the previous regions. Plasma and field go through a fast mode expansion, i.e. the flow velocity increases, and the plasma density and the magnetic field magnitude both decrease, with the accelerating force coming from both the magnetic force and the pressure gradient force.
Finally in region V, the plasma and field slowly expands further and gradually returns to their solar wind values.
A similar discussion is appropriate for the streamlines in the right panel of Fig. 8, as shown in Fig. 11. The PDL signature, with a decreased plasma density and enhanced magnetic field, is only seen on streamline 0, which is closest to the magnetopause, from ∼6 R E to ∼9 R E . There are some structures along streamline 1; however, they only occur far downstream in the magnetosheath and far away from the subsolar point. Figure 12 shows forces in a streamline coordinate system along streamline 0 in the right panel of Fig. 8, together with plasma and field parameters. The display is the same as in Fig. 10. Regions I, II, and III are very close to the y=0 plane, but regions IV and V start to deviate from this plane. The force features are similar to those shown in Fig. 10. However, some significant differences exist. The positive pressure gradient force along magnetic field exists in regions I, II, and the first half of region III, as shown in panel (e) of Fig. 12. This force is the only force driving plasma flow along magnetic field line to deplete a flux tube. In contrast to the region III in Fig. 10, there is both plasma acceleration and deceleration in this region, as shown in panel (h) of Fig. 12. The difference exists in region IV, too, with a net plasma deceleration. Streamline 0 in the y=0 plane ends up inside the magnetosphere through the polar reconnection; thus, there is no infinite expansion as streamline 0 in the equatorial plane. There is little spatial difference for the starting points of these two streamlines, (17, 1, 0) R E and (17, 0, 1) R E in the left and right panels of Fig. 8, respectively. However, this small difference is responsible for substantial differences along these streamlines.  Fig. 11. Values along the streamlines in the right panel of Fig. 8. The other settings are the same as those in Fig. 9. The PDL feature, decreased plasma density and enhanced magnetic field is only seen clearly on streamline 0, which is closest to the magnetopause, from ∼6 R E to ∼9 R E .

Flux tube depletion in the magnetosheath
As magnetic field lines are piling up toward the magnetopause, the plasma density and the plasma pressure are also trying to increase. In a one-dimensional model this would be the only possible solution. However, in three-dimensional models, as the one that we use in this study, the resulting pressure gradient can drive plasma away from the subsolar region along magnetic field lines. This depletion overcomes the plasma pile-up process and eventually produces a plasma density decrease toward the magnetopause. Thus, it is necessary to study the detailed force and flow patterns along magnetic field lines in the magnetosheath, to understand the plasma depletion process. In order to do so, we choose three field lines in the y=0 plane, which are shown in Fig. 13. The open-closed magnetic field boundary is shown as a red zigzag curve and the plasma density is shown color coded as the background. From right to left, the field lines are marked 1-3.
The plasma and field parameters along field lines 1-3 in Fig. 13 are shown in Fig. 14. The horizontal axis is the distance along each magnetic field line from its starting point in the equatorial plane. In each panel, the shaded region is the bow shock, the region to the left of the bow shock is the magnetosheath, and the region to the right of the bow shock is the solar wind. From top to bottom we show: the plasma density, the magnetic field magnitude, the plasma pressure, the flow velocity along and perpendicular to magnetic field, the pressure gradient force along magnetic field, the pressure gradient force perpendicular to magnetic field, and the magnetic force (perpendicular to magnetic field). The unit of the forces in panels (f)-(h) is 10 −19 Pascal/m. Figure 14 shows that there is a strong pressure gradient force perpendicular to the magnetic field on the bow shock to decelerate solar wind plasma. This deceleration force decreases as the magnetic field line moves further into the magnetosheath. The magnetic force (perpendicular to magnetic field) on the bow shock has a similar pattern to the pressure The display is the same as in Fig. 10.
gradient force perpendicular to magnetic field, except that it decreases much faster and becomes more and more insignificant compared to the pressure gradient force as the field line moves further into the magnetosheath. On the bow shock, the pressure gradient force along the magnetic field line gives plasma a strong kick away from the equatorial plane. The bow shock pushes the newly shocked plasma away from the equatorial plane and this occurs during every passage of the flux tube in the magnetosheath. Also, it decreases when the flux tube moves further downstream into the magnetosheath. Inside the magnetosheath, forces perpendicular to magnetic field, i.e. the pressure gradient force perpendicular to magnetic field and the j ×B force, are responsible for field line deceleration, acceleration, and draping around the magnetopause. These two forces are relatively small in most of the magnetosheath (lines 1-2) until close to the magnetopause (line 3), where they have distinct peak magnitudes close to the equatorial plane. The pressure gradient force along the flux tube exists inside the entire magnetosheath and further drives plasma away from the equatorial plane. This result is different from the description of Zwan and Wolf (1976), who believed that this depletion can only occur at the bow shock and near the nose of the magnetopause. Different from the other two forces in panels (g) and (h), the pressure gradient force along the flux tube only increases slightly when the field line moves toward the magnetopause. Very close to the magnetopause, along line 3, a complex pressure gradient force pattern in both directions along the flux tube exists. This complex force pattern is responsible for the disturbed structures in the plasma density and the magnetic field magnitude profiles in panel 3. Especially between the two dashed lines in panel 3, there is an enhancement of the plasma density and a decrease in the magnetic field magnitude, which is the signature of a slow mode wave. However, such a feature does not exist along the radial lines on the y=0 plane, as shown in Fig. 3 and Fig. 5. Near the equatorial plane, along lines 0-3 in Fig. 14, the flow velocity parallel to magnetic field is close to zero. It increases away from the equatorial plane, which looks like a fast mode expansion wave because the plasma density and the magnetic field magnitude decrease concurrently. However, the expansion is also very closely aligned with the magnetic field direction which would also allow for a slow mode expansion fan. We speculate at this point that the features seen in Fig. 14, line 3 are those of a slow mode expansion fan. We will address the existence of such a wave in a forthcoming paper.

Discussion
The conventional wisdom that the plasma density and pressure decrease as a result of the magnetic pressure increase, to keep the total pressure balance in the PDL, is only an approximation. For example, as shown in panel (h) of Fig. 10, the combined force along the flow direction can be comparable to both the pressure gradient force and the magnetic force. There is usually some inertial force so that the magnetic and plasma forces cannot precisely balance. As shown in Fig. 7, flow acceleration and deceleration exist almost everywhere in the magnetosheath. Zwan and Wolf (1976) used Spreiter's model results in their MHD PDL model to supply a critical pressure boundary condition for the thin flux tube evolution in the magnetosheath. Thus, although both Zwan and Wolf's model and our model use the MHD approach, specific differences exist between our results. We present a new, more detailed model to describe the formation of the PDL. In the new description, the pressure gradient force is responsible for the plasma depletion in the flux tube. By accelerating newly shocked solar wind plasma along magnetic field line away from the equatorial plane, the bow shock always plays an active role in depleting a flux tube during its passage in the magnetosheath. This pressure gradient force at the bow shock acts on the newly shocked plasma by providing significant acceleration along the flux tube. Plasma depletion along a flux tube inside the magnetosheath is also provided by the pressure gradient force along the flux tube. Flux tube depletion occurs during all of the passage of a flux tube in the magnetosheath before it reaches far downstream of the magnetosheath. A little further upstream of the magnetopause in the magnetosheath, the pressure gradient force along the magnetic field line pushes plasma away from the equatorial plane. However, closer to the magnetopause, a more complex pressure gradient force pattern exists, where in some part of a flux tube the force is pointing toward the equatorial plane. This feature may be indicative of a slow mode front, as proposed by Southwood and Kivelson (1992). Depending on which way a virtual spacecraft goes along a field line close to the magnetopause, a slow mode feature with enhanced magnetic field and decreased plasma density, or a slow mode feature with decreased magnetic field and enhanced plasma density

Fig. 14.
Parameters along field lines 1-3 in the noon midnight meridian in Fig. 13. The horizontal axis is the distance along each magnetic field line northward from its crossing on the equatorial plane. In each panel, the shaded region is the bow shock, the region to the left of it is the magnetosheath, and the region to the right of it is the solar wind. From top to bottom are: the plasma density, the magnetic field magnitude, the plasma pressure, the flow velocity along and perpendicular to magnetic field, the pressure gradient force along and perpendicular to magnetic field, and the magnetic force (perpendicular to magnetic field). The unit of the forces in panels (f)-(h) is 10 −19 Pascal/m. may be observed. This explains PDL observations with both types of slow mode features inside the magnetosheath (Song et al., 1990(Song et al., , 1992. However, as we have shown in Fig. 4 and Fig. 5, only the slow mode feature with enhanced magnetic field and decreased plasma density exists toward the magnetopause along the radial lines for the stable solar wind conditions in this study. Southwood and Kivelson (1992) proposed that, if the flow is to be diverted away from the magnetopause, the pressure must rise across the slow mode front. The field must drop if the pressure rises. This leads to density enhancement and field decrease after the slow mode front. This slow mode structure would correspond to the outer layer of the observed two-layered slow mode structure, i.e. a density increase/field decrease, followed by the PDL. Such features have been reported by Song et al. (1990Song et al. ( , 1992 and were deemed characteristic of the magnetosheath. However, in order to understand the density depletion, one has to consider the velocity divergence in the transverse direction as well. Around the subsolar point for the northward IMF case in this study, in the y direction the flow divergence leads to a decrease in both the plasma density and the magnetic field magnitude. In this direction the flow pattern is a fast mode expansion. In the z direction, however, magnetic field is not immediately affected by the flow divergence. Thus, the density can decrease without affecting the field. In the x direction, the compression tends to increase both the magnetic field and the plasma density. Combining the effects in all three directions, the flow divergence conspire to decrease the density and to increase the magnetic field magnitude toward the magnetopause. As shown in the previous section, a rather complicated force pattern causes flow divergence in the magnetosheath. However, as shown in Fig. 7, there is very little flow bending close to the magnetopause, as proposed by Southwood and Kivelson (1992). Thus, we do not see any plasma density enhancement with magnetic field decrease in front of the magnetopause. Actually, the bending occurs far ahead of the subsolar point in the magnetosheath for the streamlines that we have chosen. Along the Sun-Earth line, the flow toward the magnetopause decreases, and the plasma density also decreases, as shown in Fig. 4. This is because plasma flow is diverted to the other directions. Thus, the flow deceleration toward the magnetopause does not ensure the enhancement of the plasma density.

Summary and Conclusions
In this paper we investigate the underlying physics of the PDL using global MHD simulations. Our detailed force analysis shows: 1. Specific MHD forces play different roles for the PDL formation. The pressure gradient force along a field line is primarily responsible for plasma depletion. Both the pressure gradient force and the magnetic force are responsible for the flow pattern in the magnetosheath, and divert the plasma and magnetic field around the magnetosphere.
2. Only streamlines that originate close to the Sun-Earth line exhibit clear plasma depletion and contribute to the formation of the PDL. This is consistent with the consideration of Kivelson (1992, 1995). On such streamlines, distinct regions exist with different force features. The PDL is the result of the plasma motion controlled by those forces with their complex roles along streamline.
3. The analysis of the forces in the MHD simulation results leads us to propose a more detailed description of flux tube depletion in the magnetosheath. In this new description the bow shock always plays an important role to drive the newly shocked plasma along the flux tube away from the equatorial plane. The pressure gradient force exists in the flux tube's entire magnetosheath passage before it is close to the magnetopause, which further depletes the flux tube. Near the magnetopause, a more complex pressure gradient force pattern exists along the flux tube. These results are in contrast to the flux tube depletion description by Zwan and Wolf (1976) which only qualitatively considers the depletion effects on the bow shock and close to the magnetopause.