Breathing of heliospheric structures triggered by the solar-cycle activity

Abstract. Solar wind ram pressure variations occuring within the solar activity cycle are communicated to the outer heliosphere as complicated time-variabilities, but repeating its typical form with the activity period of about 11 years. At outer heliospheric regions, the main surviving solar cycle feature is a periodic variation of the solar wind dynamical pressure or momentum flow, as clearly recognized by observations of the VOYAGER-1/2 space probes. This long-periodic variation of the solar wind dynamical pressure is modeled here through application of appropriately time-dependent inner boundary conditions within our multifluid code to describe the solar wind – interstellar medium interaction. As we can show, it takes several solar cycles until the heliospheric structures adapt to an average location about which they carry out a periodic breathing, however, lagged in phase with respect to the solar cycle. The dynamically active heliosphere behaves differently from a static heliosphere and especially shows a historic hysteresis in the sense that the shock structures move out to larger distances than explained by the average ram pressure. Obviously, additional energies are pumped into the heliosheath by means of density and pressure waves which are excited. These waves travel outwards through the interface from the termination shock towards the bow shock. Depending on longitude, the heliospheric sheath region memorizes 2–3 (upwind) and up to 6–7 (downwind) preceding solar activity cycles, i.e. the cycle-induced waves need corresponding travel times for the passage over the heliosheath. Within our multifluid code we also adequately describe the solar cycle variations in the energy distributions of anomalous and galactic cosmic rays, respectively. According to these results the distribution of these high energetic species cannot be correctly described on the basis of the actually prevailing solar wind conditions. Key words. Interplanetary physics (heliopause and solar wind termination; general or miscellaneous) – Space plasma physics (experimental and mathematical techniques)


Introduction to solar cycle -induced variations
Since Parker (1963) published the first description of the heliospheric interface configuration, it was already obvious that the locations of the solar wind termination shock and the heliopause would vary under different inner heliospheric pressure conditions. However, these early models resulted from stationary equilibrium considerations, requiring that the boundary conditions do not change in time. Model calculations with these restrictions are completed, as soon as an asymptotic state in the locations of the heliospheric structures is reached. Truly dynamical variations of the heliosphere could only be studied in recent years, after refined numerical codes had been developed to model the solar wind -interstellar medium interaction. Following the course of the publications in recent years one can identify a gradual increase in numerical and physical perfections of the models (for a review, see Zank, 1999a). Based on the virtues of these models, studies of time-dependences have meanwhile become feasible.
First studies in this respect showed that the response reaction of the heliospheric shock properties very much depends on the dimensionality in which the problem is described. Barnes (1993Barnes ( , 1994Barnes ( , 1995, Grzedzielski and Lazarus (1993), and Naidu and Barnes (1994) and Belcher (1998, 1999) have considered the motion of the termination shock in response to a jump in the solar wind density or an interplanetary shock passage, using a one-dimensional planar approximation. Ratkiewicz et al. (1996) then showed that the resulting responses viewed in compression ratios and propagation speeds of the shock, are much less pronounced in a spherically symmetric approach. Similar studies, however, already based on measured perturbations of solar wind properties, were carried out by Whang and Burlaga (1993) and Whang, Lu and Burlaga (1999). These latter authors conclude, that the location of the upstream termination shock is anti-correlated with the sunspot number, i.e. moving inward or outward at the rising or the declining phase of the solar cycle. Since then it is recognized that the heliospheric response to inner heliospheric solar wind variations can only be ade-quately described if the whole interaction system is consistently estimated, including all involved time-dependent processes, especially, for example, the back reactions from the heliopause and the bow shock.
In this sense Steinolfson (1994) and Karmesin et al. (1995) have run the first two-dimensional, time-dependent interaction codes to study solar cycle effects and could show that the resulting migrations of the shock location in their timedependent studies turn out to be much more moderate compared to migrations achieved if the shock could instantly react to the actual inner solar wind conditions. This broadened view has then still been substantially improved by Zank (1999b), emphasizing the fact that the interstellar neutral Hatoms, a dynamically important component of the interaction process, was not yet taken into account. By including Hatoms in a consistent form into a two-dimensional, two-fluid interaction code he showed that the earlier arythmic motion of the termination shock is now reduced as a consequence of a reduction of the upwind to downwind shock layer asymmetry and of the shock distance from the Sun.
Appreciating the relevance played by the H-atom fluid in this context, we now feel strongly encouraged to revisit the problem of heliospheric interface responses to solar cycle variations, however, including now, in addition to the approach presented by Zank (1999b), three more fluid components of dynamical relevance within this context. These latter components are pick-up ions (PUIs), anomalous cosmic rays (ACRs) and galactic cosmic rays (GCRs) (see .

The solar-cyclic response of the 5-fluid heliosheath
The theoretical basics of the 5-fluid model has already been published by Kausch (1998 and Fahr (2000). Thus, we do not want to describe all the physical and numerical details here again. We only briefly summarize that besides protons of solar and interstellar origin, we consistently include interstellar H-atoms, heliospheric Hpick-up ions, H-ACRs and GCRs as separate fluids which are dynamically and thermodynamically coupled to build the multifluid interaction system. PUIs are co-moving with the local solar wind, but behave thermodynamically as an independent fluid with separate values of density, temperature, sound speed and pressure. ACRs are injected at the shock with a rate dependent on the local strength of the shock into the Fermi-1 acceleration process. This rate corresponds to a specific fraction of the PUI flux which locally passes over the shock. This injection constitutes a local ACR source at the shock which is taken into account in an energy-averaged ACR transport equation describing convection and spatial diffusion of this high-energy 10 MeV/nuc fluid. Finally, we include the GCRs by coupling to the other fluids an energyaveraged GCR transport equation. Hereby, the coupling occurs due to convective GCR changes connected with the local plasma bulk flow, and due to a modulation of the plasma motion by ACR and GCR pressure gradients. The parameters for a stationary heliosphere are shown in Table 1. These model calculations, after reaching its asymptotical state, are used as input for modeling the solar cycle variations. All different solar cycle models discussed below use this stationary state as a starting point for the new calculation.
The high-energy components, i.e. the ACRs and the GCRs, also represent two additional fluid species diffusing in space due to being scattered at MHD-wave turbulences. These species, in a mathematical sense, enforce the upgrading of the system of coupled partial differential equations to a system of second order equations. The inclusion of PUI-, ACR-and GCR-fluids in a physical sense opens up new reaction modes for the whole interaction system, if timedependent inner boundary conditions modeling typical solar cycle variations are installed.
As is well known, such variations have been reported on the basis of deep space probe data obtained with PIONEER and VOYAGER spacecraft by Kayser, Barnes and Mihalov (1984) and Lazarus and McNutt (1990). They essentially show cycle-periodic variations of the solar wind mass flow and dynamical pressure roughly by a factor of 2, with maxima delayed by roughly 4 years with respect to the activity cycle maximum. As can be extracted from data presented by Lazarus and Belcher (1987), Gazis (1994) Table 1 is shown, while in the lower panel the dynamic heliosphere is presented, in which the solar wind speed and density at the inner boundary varies according to a sinusoidal law between 300 km/s and 500 km/s. For the varying distances of the termination shock and the bow shock see Fig. 4. wind mass flux, since solar wind velocity and density appear to be appropriately anti-correlated. In view of these findings we prefer here models in which the mass flow is constant during the solar cycle, as described in the next section.

Solar cycle models
At a distance of 5 AU, taken as the inner boundary of our simulation model, the radially symmetric inner solar wind velocity becomes time-dependent, with velocity magnitudes varying periodically with the solar cycle. The solar wind mass flow, e.g. ρ sw v sw , thereby is kept constant, meaning that the density ρ sw is anti-correlated with the velocity v sw . This also means that the velocity varies by the same factor, as does the ram pressure ρ sw v 2 sw . We modeled two types of solar cycle activity: (1) a pure sinusoidal variation and (2) a more complicated but perhaps a more realistic one, with an exponentially modulated amplitude, given by Fahr et al. (1987) and Fahr and Scherer (1990): where ω s is the solar cycle frequency. v sw0 ∼ = 359 km/s and v sw ∼ = 162 km/s are chosen in such a way that the solar wind velocity during maximum conditions reaches v max = 800 km/s while during minimum conditions it is v min = 300 km/s. The phase φ is needed to adapt to a stationary heliosphere with a solar wind speed of v sw stationary = 400 km/s.
We also used a sinusoidal law to model the solar cycle variations: where the parameters v sw 0 , v sw , ω s , and φ have the same meaning as in Eq. (1).
In the preceding models in the literature, the authors varied either the solar wind speed or the solar wind density. In our models we kept the solar mass flux constant, requiring the variability of both the solar wind speed and the solar wind density simultaneously. Therefore, the associated solar wind density at the inner boundary is adopted such that the resulting mass flux at 5 AU is constant during the course of the solar cycle and is equal to m0 (5 AU) = v sw · n sw (5 AU)m p = 8 · 10 6 m p cm 2 s −1 , which is approximately the value measured at the Ulysses spacecraft m0,U = 1.2 · 10 7 m p cm 2 /s −1 (McComas et al., 2000).

Snapshots of the large-scale heliospheric structure
In Figs. 1 to 3 the proton and hydrogen density of the heliosphere are shown in the form of color-coded plots. In the following displays we keep the set of the static heliosphere parameters as standard and only vary in time the solar wind speed v sw = v sw (t) as discussed above, according to Eqs. ( 1) and (2), where the solar mass flux is constant and hence, the solar wind number density n sw = n sw (t) is anti-correlated to the solar wind speed. The period of the solar cycle is kept constant and is adopted with 11 years.
For a pure sinusoidal solar-cyclic variation of the solar wind speed between V sw min = 300 km/s and v sw max = 500 km/s, density contours for protons at a specific phase of the solar cycle (approximately in the middle of the cycle) are shown in Figs. 1b, 2 and 3b (i.e. a phase snapshot of the dynamic heliosphere). As one can easily identify by eye, the main structure features of the dynamic heliosphere appear to be a little blown up compared to those of the static heliosphere, both when viewed in proton and hydrogen densities. Furthermore, in the dynamic heliosphere one can clearly see (Fig. 1b) two proton density waves propagating outwards from the heliopause towards the bow shock. The waves traveling outwards, as can be controlled in a cinematographic sequence of successive heliospheric states, clearly are pointed out as subsonic waves running into the bow shock. This most likely should also have implications for a post-terminationshock acceleration of ACRs and for a modulation of GCRs and ACRs downstream of the termination shock.
Since the proton density in this contour plot varies over several orders of magnitude, we present in addition to the linear scale a logarithmic scale (see Fig. 2), to enable the visibility of all the wavy features, especially in the heliotail, which can be observed in the inner and the outer heliosphere, as found in our models. An animation of these features can also be found under http://www.aef-ev.de/Kom03/heliopuls. Similar features can be observed for the hydrogen density, as presented in Fig. 3 within a logarithmic scale, to better pronounce the details in the heliotail. Moreover, to enhance the effects we allow for a maximal solar wind speed of 800 km/s. The waves in the hydrogen wall cannot be resolved in this scale, but the "chaotic" behavior of the hydrogen distribution inside the heliosphere is clearly seen.
Alternatively, we also studied the heliospheric solar cycle variations resulting from the time-dependence given by Eq. (1), which seems to be more realistic and better confirmed by data showing that, after an extended quiet period, the solar wind speed then rapidly increases towards its maximum and then again rapidly decreases from this maximum towards a new extended minimum. These exponential features, rather than a sinusoidal triggering, resemble more the   Fig. 1 for the hydrogen distribution. To enhance the effects the solar wind speed variations reach from 300 km/s to 800 km/s for a pure sinusoidal law. Moreover, a logarithmic scale is chosen. case of a large Global Merging Interaction Region (GMIR) running into the termination shock. We have modeled solar wind speed variations between v sw min = 300 km/s and v sw max = 500 km/s or v sw max = 800 km/s for both solar cycle models.

Distances to termination and bow shock in different solar cycle models
In Fig. 4 five different representations of solar cycle variations are used by us to model the solar-cyclic heliosphere variation. The differences in the results can then be com-pared. The minimum solar wind speed v sw min = 300 km/s is kept fixed in all of the five models. The black line represents a sinusoidal modeling with a maximal solar speed of v sw max = 500 km/s. In this model no PUIs, ACRs, and GCRs are taken into account, i.e. their dynamical role and relevance is suppressed (model 1). In contrast, in the other four models all five energetic particle species (i.e. protons, hydrogen atoms, PUIs, ACRs, GCRs) are consistently taken into account with their mutual dynamical interactions. While the red curve (model 2) and green curve ( sults for an analogous amplitude of the solar wind speed, but now for a sinusoidal variation. In the case of model 1 (in which only hydrogen and protons are taken into account), we used a similar parameter set, as discussed in Table 1, but without any interaction of the PUIs, ACRs, and GCRs as input for the sinusoidal solar cycle model (v sw max = 500 km/s). It turns out that the distance of the termination shock is larger than compared to the standard stationary model where all 5-fluid components are used. Furthermore, the TS in model 1 migrates only a small distance outwards and then more or less purely oscillates around its average position. After almost four solar cycles the minimum termination shock distance (Fig. 4c) is constant for the following cycles, which also holds true for the bow shock. In contrast, the bow shock distance for models 2-5 is still slightly growing after 10 solar cycles.
The phase lag between the sinusoidal and exponential models (Fig. 4a) is caused by the condition that they have to fit to the stationary model at t = 0. While the amplitudes are largest for the models with v sw max = 800 km/s, the form of the curves is much different: while the motion of the TS in the exponential model shows slow growth, with a flat saddle in the slope, the sinusoidal models increase fast and mono-tonically towards their maximal position. The decrease towards the minimum condition of the TS is slow and (almost) monotonic for the sinusoidal models, while the exponential models decrease rapidly and monotonically towards the minimum position. In all four models an increasing in the minimal and maximal distance of the TS can be observed. This is shown in detail in Fig. 4c for the relative minimal distance of the termination shock and in Fig. 4d for its relative maximal distance. The relative distances are calculated as the respective deviations from the static model.
Depending on the amplitude of the solar wind speed variation in five of the solar cycle models, the distance of the bow shock migrates without any oscillation. While no equilibrium location seems to be achieved in model 1, for models 2-5 an equilibrium position is reached almost after 6-8 solar cycles. It should be mentioned here that the bow shock has a nearly stable location over one solar cycle. The jumps are occuring due to numerical uncertainties and the handling of large data sets. Only the tendency of a systematic increase in the bow shock distance over several consecutive cycles expresses a physically relevant fact.
From the above mentioned, it becomes clear that the heliosphere has a memory and recognizes at least a couple of solar cycles. Thus, after an extended period of low solar activity, as was the case in the Maunder Minimum (1645-1715, where almost no solar cycle activity was observed, e.g. Beer 2000), some time is needed for the termination shock and bow shock to reach well adapted positions. The behavior of and the importance for the cosmic ray fluxes in the outer heliosphere will be discussed elsewhere by Scherer et al. (2003).
In the following, we will concentrate on model 5, i.e. the exponential law with a solar wind speed variation from 300 km/s to 800 km/s, because it shows the most pronounced features. All what will be discussed below holds true for the other models 2-4.
Instead of showing a snapshot in time, as in Figs. 1-2, we now show in Figs. 5-8 an isochronic cut through space, i.e. along one direction (in this paper we chose the stagnation line), as a function of time. There, it is nicely seen that the maximum distances to the termination shock, the heliopause and the bow shock systematically increase for about 6 solar cycles. Only after about that time do these structure features purely oscillate in distances with the solar cycle.

Hydrogen waves
In Fig. 5 the hydrogen density variation is shown, while the corresponding features for the proton density are discussed in Scherer and Fahr (2002). As expected. the region ahead of the heliopause appears as a region of enhanced Hatom densities. However, in the dynamic heliosphere the socalled "hydrogen-wall"-structure turns out to be a substantial time variable feature. The changes in the formation of the wall configuration are associated with H-atom density waves propagating inwards and outwards from this wall, with about equal velocities of the order of 25 km/s. In addition, it can be recognized that the thickness (extent) of the hydrogen wall is periodically varying according to the frequency of the solar cycle. This very interesting phenomenon could be directly controlled by observation of hydrogen-wall modulated absorption features in stellar hydrogen spectra near H-Lyman-α (see Linsky and Wood, 1996;Linsky et al., 1993), because the spectral imprints to these spectra from the hydrogen wall should periodically vary over the solar cycle.
But not only can heliospheric hydrogen be detected by observing Lyman-α line features. In addition, inconsistencies in the Lyman-α spectra intepretations from other stars are also explained as resulting from their astrospheres, i.e. as their "stellar heliospheres" (Wood et al., 2000). Thus, stellar cycle induced astrophere changes are also reflected in the variations of their hydrogen wall, and the observation of the resulting Lyman-α emission will allow for an indirect observation of time-varying stellar winds around Sun-like stars.

Shock waves and compressible waves
To study the characteristic features of the waves traveling towards the bow shock, we estimated the proton Mach-number. In Fig. 6 it is plotted along the stagnation line for almost 11 solar cycles. The heliopause is recognized here as the black (dark blue) area between 100 AU and 139 AU. It can also be seen that subsonic regions extend towards the termination shock, as well as towards the bow shock. But, nevertheless, large areas in the heliosheath are filled with traveling subsonic waves.
The radial velocities of the protons between 100 AU and 180 AU (Fig. 7) also show a systematic change from a small inflow to a small outflow of material between the heliopause and the bow shock. This should provide an increase in the MHD-turbulence. Its estimation must be done elsewhere with a different model as the 5-fluid approach. Nevertheless, from our model we can presently conclude, that the periodic increase and decrease in the Mach-number, as well as in the   7. The variation of the radial proton velocity in the sheath between termination shock and bow shock proton velocities, causes similar increases and decreases in the turbulence levels. The latter affects the transport of energetic particles, as well as the effectiveness of the acceleration of ACRs. Based on our model data presented here, Scherer et al. (2003) discuss some of the consequences for GCRs.

Merging waves in the tail region
Figure 8 also shows a time versus distance plot for the proton density in the heliotail. Here, the densities are very small, but, nevertheless, the solar activity cycles blow even more "holes" into the density patterns. These depression areas start to expand and then propagate further out where they slowly decay, i.e. they induce proton density increases at larger distances. Furthermore, density waves rising from different solar cycles merge at even larger distances, which also should increase the MHD-turbulence levels. These merging waves may accelerate newly ionized particles, which are created by charge exchange processes between the protons and hydrogen atoms in the heliotail, to some keV and which can be observed on SOHO (Czechowski et al., 2001).
The merging wave regions are moving toward the termination shock, but did not reach it after 10 solar cycles. Therefore, we will have to perform longer runs of our models to better study the long-term behavior in the tail.

More dimensions needed?
The heliosphere is an essentially three-dimensional structure, not only due to the frozen-in magnetic fields carried by the solar wind, but also due to the three-dimensional asymmetric distribution of the solar wind momentum flow. The effects of the three-dimensional distribution of the solar wind momentum flow during solar minimum may be handled by separate static axisymmetric models, under the assumption that as observed by Ulysses (McComas et al., 2000), there exists a sharp transition region between the slow and fast solar wind. But the principle physical effects of the fast solar wind and the slow ecliptical solar wind can only be studied by dynamically fully developed 3-D models, because the transport of the cosmic rays depends on the different turbulence levels caused by the three-dimensional structure of the frozen-in magnetic field as well as that of the 3-D solar wind momentum distribution. Hence, the solar cycle dependence studied here is a first attempt to understand the difference between static and dynamical models in a 2-D approximation for different ecliptical boundary conditions simulating an ecliptic cut through the three-dimensional structure of the heliosphere. Thus, the models 2-4, discussed above, give a first glimpse of the nature of the three-dimensional heliosphere.
On the other hand, the three-dimensional magnetic field induces different polarities in both hemispheres, which can be observed in the modulation of cosmic rays (e.g. Potgieter, 1998;. These effects cannot be modeled in our present five-fluid code. Nevertheless, axisymmetric models, with a kinematic treatment of the magnetic field to calculate the diffusion coefficients, describe the CR-transport in the ecliptic quite well. Therefore, the modulation of CRs in our model is a first attempt to include these important species in heliospheric models. In addition to the above discussed three-dimensional structure of the heliosphere from the inside, an outside LISM magnetic field may have to be expected, giving an imprint on the heliospheric configuration (Ratkiewicz et al., 1996). At this time, no three-dimensional models exist which include the frozen-in solar wind magnetic field as well as the magnetic field of the LISM, nor do they handle sophisticated cosmic ray transport models. Thus, our present five-fluid code can be understood as a first step in the development of more sophisticated heliospheric models.

Summary
We presented time-dependent models for the large-scale heliosphere, in which the solar cycle dependence was described by a periodic variation of the solar wind speed and density, respectively. Keeping the solar wind mass flow constant leads to moderate variations, which causes an overall migration of the termination shock, as well as the bow shock to larger distances.
In the heliosheath region solar cycle induced plasma waves and waves of neutral atoms (i.e. hydrogen) are merging into the bow shock. These waves affect the propagation of cosmic rays, as well as the acceleration of its anomalous component.
It turns out that periodic variations of the solar wind change the structure of the heliosphere remarkably. For the first time we could model the "heliospheric memory" (see, also, Scherer and Fahr, 2002), i.e. the fact that not only the last solar cycle influences the structure of the outer heliosphere but up to 6-7 previous cycles (and even more in the tail region) must be taken into account to correctly model the dynamic heliosphere and its hysteresis effects. Because the large-scale structure and physics of the heliosphere is used in related fields, such as the propagation of cosmic rays, production of energetic neutral atoms, the interpretation of Lyman-α line from other stars, etc., the dynamical heliosphere should be taken into account in future applications.