© ESO 2021

1. Introduction

The solar atmosphere can be divided into three layers with different plasma properties: the photosphere (with the conventional surface of a star) extending up to only 500 km above the surface, the chromosphere ranging from the top of the photosphere to about 2500 km in height, and the solar corona from 2500 km and expanding into the solar wind. In the photosphere, the temperature decreases with height from about 5600 K at the bottom until it reaches its minimum of 4300 K at the top of the photosphere (or the bottom of the chromosphere), which is about 500 km above the solar surface (Athay 1976). Higher up, the temperature starts rising again with height, first slowly in the lower chromosphere and then faster in the upper chromosphere until the transition region, which separates the chromosphere from the solar corona. There, the temperature experiences its sudden increase to 1 million K and from there upwards it increases steadily with height in the low corona. A clear explanation of this temperature increase with height above the solar surface remains to be found. As a result of the lower temperatures, the plasma in the low atmospheric layers is only partially ionized (Avrett 2003). However, the corona, where the average temperature is about 1 − 3 MK, is fully ionized (Aschwanden 2005).

Multiple solar missions, such as the Solar Dynamics Observatory (SDO) and the Interface Region Imaging Spectrograph (IRIS), have shown that a diversity of waves occur in the solar atmosphere (Jess et al. 2009; McIntosh et al. 2011; Okamoto & De Pontieu 2011). The various wave types that occur include Alfvén waves (Alfvén 1942). These waves are transverse magnetohydrodynamic (MHD) waves that travel along the magnetic field lines. Alfvén waves were reported to be present in both the photosphere and chromosphere (Srivastava et al. 2017; Baker et al. 2021). As they pass by, they modify the transverse magnetic field and velocity components but, do not alter the gas pressure or the mass density (Nakariakov & Verwichte 2005), at least in the linear limit and in a homogeneous background plasma. A thorough understanding of Alfvén waves is essential because they could be a part of the solution to the major problems of heliophysics, such as the solar coronal heating and wind acceleration (Uchida & Kaburaki 1974; Ofman 2010). Recent theoretical research revealed that Alfvén waves can carry enough energy to heat the solar corona (Yang & Xiang 2016). However, the details of the mechanism(s) of the thermal energy release related to their dissipation remain unknown. One potential candidate for that may be associated with ion–neutral collisions (Soler et al. 2017). Piddington (1956), Osterbrock (1961), and Haerendel (1992) were the first to study ion–neutral collisions, but they did not find that this interaction affects the chromospheric temperature. Ballester et al. (2018) showed that ambipolar diffusion leads to substantial chromospheric heating, and Zaqarashvili et al. (2013) derived a dispersion relation for two-fluid Alfvén waves and confirmed that the damping of Alfvén waves resulting from the ion–neutral collisions is quite significant. Khomenko (2017), based on a two-fluid model, stated that the presence of neutrals affects the solar atmosphere. The effect of ion–neutral interactions is expected to influence the energy balance of the chromosphere. Zaqarashvili et al. (2013) also confirmed that low- and high-frequency photospheric Alfvén waves might not reach the solar corona because ion–neutral collisions damp them very efficiently in the upper chromosphere. According to Song & Vasyliūnas (2011), the rate of Alfvén wave damping varies with magnetic field strength and wave frequency. For a strong magnetic field, wave damping is low. Low-frequency waves are also weakly damped, and so there is a chance to detect low-frequency Alfvén waves in the solar corona under the condition of a strong magnetic field.

According to Hollweg (1978, 1981), Hollweg et al. (1982), Kudoh & Shibata (1999), and Matsumoto & Suzuki (2012), Alfvén waves can potentially also induce plasma outflows. More recent studies by Yang & Xiang (2016) and Shestov et al. (2017) indeed show that Alfvén waves manage to produce fast plasma outflows. As Tu et al. (2005) stated, we can detect outflows in the solar atmosphere between 5 and 20 Mm altitude, and at the height of 20 Mm, the outflow speed reaches about 10 km⋅ s⁻¹. A notable problem is to examine such outflows in the lower solar atmosphere, specifically in the chromosphere and transition region (Marsch et al. 2008; Tian et al. 2009; McIntosh et al. 2011; McIntosh 2012; Kayshap et al. 2015).

The goal of the present paper is to study impulsively generated two-fluid Alfvén waves in the partially ionized lower layers of the solar atmosphere. These studies are performed numerically in the context of the observed atmospheric heating and plasma outflows, which in the higher coronal regions may result in nascent solar wind.

The organisation of the remainder of this paper is as follows. Section 2 presents the governing two-fluid equations, an initial equilibrium model of the solar atmosphere, and the impulsive perturbations we applied in the simulations. In Sect. 3, the results of the numerical simulations are presented and Sect. 4 contains the conclusions of this paper with a discussion and summary of the results of the numerical experiments performed.

2. Two-fluid numerical model of the partially ionized solar atmosphere

In the dens(er) lower layers of the solar atmosphere, the interactions between the different plasma components are commonly considered to be collective, which means that each particle simultaneously affects several other particles. Therefore, it is justified to consider the plasma in these atmospheric layers as a continuum, that is, as a fluid (Roy & Pandey 2002). The temperature depends on the altitude of the considered layer in the solar atmosphere. There, where the temperature is high enough, as in the transition layer and the corona mentioned above, the plasma is fully ionized. However, as mentioned in the introduction, in the photosphere and chromosphere the temperatures are lower and the plasma there is not fully ionized. In the photosphere, the ionization degree is even as low as 0.01%. As discussed above, it is therefore necessary to take into account the dynamics of neutral particles in these atmospheric layers. As a simple conceivable approach, we use the two-fluid plasma model to describe the partially ionized atmosphere, which consists of ions+electrons and neutrals, treated as two separate fluids (Zaqarashvili et al. 2011). Each of these two fluids has its own mass density, flow velocity, and gas pressure and the fluids interact with each other through ion–neutral collisions.

2.1. Two-fluid equations

The two-fluid approach combines the fluid mechanics theory (governed by the Navier-Stokes equations) and electromagnetism (governed by Maxwell’ s equations). The evolution equations in the framework of a two-fluid model can therefore be written in the following form (Ballester et al. 2018; Wójcik et al. 2019):

(1)

(2)

(3)

(4)

(5)

(6)

with

(7)

(8)

where

(9)

(10)

(11)

(12)

Here, 𝜚_i, n is the mass density of the ions and neutrals, respectively, and similarly V_i, n are the velocity fields, p_i, n the gas pressures, E_i, n the total energy densities, while g = [0, −g, 0] denotes the gravitational acceleration with its magnitude g = 274.78 m⋅s⁻². Moreover, k_B is the Boltzmann constant, I denotes the identity matrix, μ the magnetic permeability, B the magnetic field, γ = 1.4 the adiabatic coefficient, and T_i, n corresponds to the temperatures of the two fluids. Also, m_H is the mass of a hydrogen atom, and μ_i = 0.58,  μ_n = 1.21 denote the mean mass of each species. Hence, the indices _i,e,n correspond respectively to ions, electrons, and neutrals, while ν_i, n is the ion–neutral collision frequency which is given as (Braginskii 1965; Ballester et al. 2018)

(13)

with σ_in being the collision cross-section, which varies with the energy of the colliding particles, that is, with temperature. Here, we adopt its classical value, σ_in = 1.4 × 10⁻¹⁹ m², (Vranjes & Krstic 2013).

2.2. Magnetohydrostatic equilibrium

The solar corona continuously expands into interplanetary space. Yet, for computational economy during numerical simulations, we assume that the solar atmosphere is in hydrostatic equilibrium. This is the state in which the radially inward force of gravity and the radially outward gas pressure gradient force are in balance. This state is described by the hydrostatic condition, expressing this force balance, that is,

(14)

We initialise the velocity, temperature, pressure, and mass density of both plasma components from the hydrostatic equilibrium state given by the solution of Eq. (14). This means that the velocity of both ions V_i and neutrals V_n is set equal to zero. Additionally, we set the initial temperature to be equal for both plasma components, T_i = T_n = T, according to the semi-empirical model of Avrett & Loeser (2008). The following equations describe the initial equilibrium values of gas pressure and mass density (e.g., Kuźma et al. 2017):

(15)

(16)

(17)

Here, Λ_in denote the ion and neutral pressure scale heights, while p_0in are the plasma and gas pressures at the reference height, y_r is the reference height taken as 50 Mm, with p_0i = 10⁻² Pa and p_0n = 3 ⋅ 10⁻⁴ Pa, and m_H corresponds to the mass of a hydrogen atom.

We also assume that a magnetic field penetrates this hydrostatic equilibrium state, but we consider a magnetic field that is force-free (i.e., (∇ × B) × B/μ = 0) and even current-free (i.e., (∇ × B)/μ = 0), so that the hydrostatic force balance is not disturbed. A uniform vertical magnetic field , with a magnitude B_y that is chosen to be equal to 30 G, satisfies these conditions or equations. Here, denotes a unit vector directed along the y-axis. The magnetohydrostatic equilibrium model assumes a uniform and vertical magnetic field, while in the real solar atmosphere the magnetic field expands with height. Recent works that include such magnetic field line expansion include those by Soler et al. (2017, 2019). This expansion is more important in the low chromosphere and our model requires a revision there.

Figure 1 (top panel) illustrates the vertical variation of the initial equilibrium temperature T(y). This plot reveals that the minimum temperature is located 100 km above the photosphere (occupying the region given by 0.0 Mm ≤ y ≤ 0.5 Mm), mainly at y = 0.6 Mm, and T(y) becomes slightly higher in the middle and upper chromosphere (0.5 Mm ≤ y ≤ 2.1 Mm). A sudden increase in temperature occurs at the height corresponding to the transition region (y ≈ 2.1 Mm). Above the transition region, in the solar corona, the temperature keeps rising with height until it reaches a magnitude of 1 MK at y = 20 Mm.

Fig. 1. Variation with height y of the equilibrium temperature T, given by Avrett & Loeser (2008) (top), and Alfvén speed (bottom).

The bottom panel of Fig. 1 displays the vertical variation of the Alfvén speed c_a, which is defined as

(18)

At the bottom of the photosphere, located at y = 0 Mm, c_a ≈ 200 m⋅s⁻¹. We note that c_a grows with height and the most drastic growth takes place at the transition region. However, in the solar corona, c_a(y = 20 Mm)≈10⁶ m⋅s⁻¹, and it grows steadily with increasing y (height).

2.3. Impulsive perturbations

In order to generate Alfvén waves, we perturb the magnetohydrostatic equilibrium by initially (i.e. at t = 0 s) launching a localized pulse in the transverse component of ion velocity, V_iz. The following equation describes this pulse:

(19)

Here, A denotes the amplitude of the pulse, w its width, and y₀ its vertical location. By varying A and y₀, but holding the width w fixed at 0.2 Mm, a similar pulse of different amplitude A is launched in the photosphere or chromosphere.

3. Numerical results

In order to study Alfvén waves in the inhomogeneous and only partially ionized lower layers of the solar atmosphere, we have to rely on numerical simulations. In our model, initially (at t = 0 s), an ion velocity pulse described by Eq. (19) is injected into the simulation domain specified as ( − 0.08, 0.08)×(−0.5, 60) Mm². We use the JOANNA code (Wójcik et al. 2020) to perform the numerical simulations, that is, to solve the two-fluid equations presented and discussed in the previous section. Indeed, this code solves the initial-boundary value problem for the specified two-fluid model. All plasma variables are initialised with their magnetohydrostatic values at the lower and upper boundaries of the simulation domain and held fixed in time throughout the simulation. Setting non-reflecting or absorbing boundary conditions along the gravity action is a formidable and still not fully solved task. We have implied the simplest conceivable boundary conditions with all plasma quantities set to their equilibrium values. Such conditions being supplemented by a well stretched grid along the y-direction is found to significantly reduce reflections of the incoming signal. The part of the simulation domain specified by −0.5 Mm ≤ y ≤ 4.62 Mm is divided into 512 grid cells, leading to a numerical grid cell size Δy = 10 km in that area. The region higher up is covered by 128 grid cells and in this part of the computational domain the grid is stretched, that is, the size of the grid cells steadily grows with height. At the side (x-) boundaries, so-called “open” boundary conditions are implemented, meaning that the x-derivatives of all the plasma quantities are set to zero.

3.1. Numerical verification test

Before starting any simulation, we first verify the numerical accuracy by trying to quantify the numerical errors. One test consists of running the code without any initial pulse and relaxing the initial (analytical) equilibrium state numerically, that is, allowing the physical system to evolve in time without any external influence. Although a zero–amplitude pulse cannot lead to the generation of Alfvén waves, as a result of the numerical approximation which is due to the discretization of the equations on a mesh with finite resolution, there is some signal present in the vertical (y-) component of the ion velocity and temperature in the system. This signal corresponds to magnetoacoustic waves, which are affected by gravity. We can trace these waves by plotting the V_iy component in the (y, t)plane (Fig. 2, top) and also plotting the relative perturbed ion temperature,

(20)

Fig. 2. Time–distance plots for Viy (top) and δTi/T (bottom) for the initial pulse-free system, i.e., A = 0 km⋅s−1. In other words, this plot shows the numerical relaxation of the initial equilibrium state.

see the bottom plot of Fig. 2. We note that the obtained max (V_iy)≈0.34 km⋅s⁻¹, which takes place in the corona at y ≈ 20 Mm, and it is very small compared to the local sound speed c_s ≈ 100 km⋅s⁻¹ there (not shown). The signal propagates upward at an average velocity of 125 km⋅s⁻¹, which corresponds to the average sound speed (the average Alfvén velocity in that area is an order of magnitude higher). Similarly, max (δT_i/T) = 0.0045, which is also a very small perturbation and it propagates upward at the same velocity as the V_iy perturbation. Hence, we infer that for the chosen numerical grid and discretization scheme, the generated numerical noise in the relaxation phase is rather small and does not affect the numerical results presented below.

3.2. Small-amplitude Alfvén waves

In this section, we consider small-amplitude Alfvén waves generated by an initial pulse. Figure 3 displays the evolution of the Alfvén waves that are excited by the initial pulse with an amplitude as small as A = 1 km⋅s⁻¹ and launched from the bottom of the photosphere, that is, y₀ = 0 Mm (top plots), and at the top of the photosphere, at y₀ = 0.5 Mm (bottom plots). The initial pulse splits into two counter-propagating waves, which are damped. This damping results from collisions between ions and neutrals (Leake et al. 2005; Arber et al. 2016; Shelyag et al. 2016; Erdélyi & James 2004; Khodachenko et al. 2004; Forteza et al. 2007; Soler et al. 2019). For the initial pulse at the bottom of the photosphere, the upwardly propagating waves experience partial reflection at y ≈ 0.5 Mm (upper plots), corresponding to the top of the photosphere. However, when the initial pulse is located at the top of the photosphere, the reflection takes place at about y ≈ 1.0 Mm (bottom plots).

Fig. 3. Time–distance plots for Viz (left) and Bz (right) for A = 1 km⋅s−1, y0 = 0 Mm (top) and y0 = 0.5 Mm (bottom).

In the case of y₀ = 0 Mm (Fig. 3, left-top), the reflection takes place after t ≈ 1500 s, so the upward propagating waves travel at about 0.33 km⋅s⁻¹, which is the average Alfvén velocity in the photosphere (see Fig. 1, bottom). However, for y₀ = 0.5 Mm (Fig. 3, left-bottom), as a result of higher values of c_a(y) in the low chromosphere (Fig. 1, bottom), the waves propagate upward much faster and the reflection therefore occurs much earlier, that is, at about t ≈ 200 s, which means that the waves indeed propagate with an average velocity of 2.5 km⋅s⁻¹. Some waves propagate higher up, reaching the transition region and the solar corona (left). In the photosphere and the low chromosphere, the amplitudes of the upward- and downward-propagating waves decay with the distance of their propagation. This is a result of ion–neutral collisions taking place in these layers. The left panels of Fig. 3 show that the maximum value of the transversal velocity component is smaller for y₀ = 0 Mm with max(V_iz)≈2.5 km⋅s⁻¹ than for y₀ = 0.5 Mm with max(V_iz)≈3.4 km⋅s⁻¹. These maxima are reached in the corona as the plasma density is much lower there.

The right panels of Fig. 3 present plots of a transversal magnetic field B_z, which is in phase with V_iz for the downward-propagating Alfvén waves and in phase opposition for the upward propagating waves. Indeed, here we also observe upward- and downward-propagating waves from the pulse height. The upward-propagating waves experience reflection in both cases (i.e. for y₀ = 0 Mm and y₀ = 0.5 Mm) after passing a distance of about 0.5 Mm, essentially simultaneously with the reflections seen in the V_iz component (Fig. 3, left panels). However, we note that the B_z signals evolve differently from the V_z perturbations. The B_z perturbation amplitude strongly falls off with height and even becomes evanescent above the reflection height, while the V_iz component reaches its maximum values there. This difference is well known in the MHD framework, where the evolution of the B_z component is described by the following wave equation:

(21)

while the transversal component of velocity, V_z, is governed by the classic wave equation (e.g., Wójcik et al. 2017):

(22)

Hence, the increase in the Alfvén velocity c_a with height y yields a relatively large damping term (∼∂B_z/∂y) in the wave equation for B_z which does not appear in the wave equation for V_z.

3.3. Heating and plasma outflows by large-amplitude Alfvén waves

In this section, we explore the dynamics of large-amplitude Alfvén waves generated by an initial pulse in the lower atmospheric layers of the solar atmosphere. We limit our discussion to the case of an initial pulse of A = 10 km⋅s⁻¹. The time–distance profiles of V_iz and B_z are displayed in Fig. 4 (top), and they globally look similar to those illustrated in Fig. 3 which corresponds to A = 1 km⋅s⁻¹. We note that the magnetic field perturbations now travel higher up in the solar atmosphere before being dissipated. The bottom panels of Fig. 4 show the corresponding velocity drifts: V_iz − V_nz (bottom-left) and V_iy − V_ny (bottom right). In both bottom panels, it can be seen that the velocity drifts grow with height towards maximum values of V_iz − V_nz ≈ 1 km⋅s⁻¹ and max(V_iy − V_ny)≈10 km⋅s⁻¹, respectively. Due to the fact that the max(V_iy − V_ny) value is higher than the max(V_iz − V_nz) value, we infer that magnetoacoustic waves, which are driven by the Alfvén waves, are responsible for the observed heating, which is described by the first term of the right-hand side of Eqs. (11) and (12). The nonlinear coupling of Alfvén waves with magnetoacoustic waves is due to the ponderomotive force. A study of the ponderomotive force in the multi-fluid framework was recently presented by Martínez-Gómez et al. (2018), who investigated impulsively generated Alfvén waves. On the other hand, the fact that the nonlinearly driven magnetoacoustic waves heat the plasma more efficiently has already been discussed; namely by Arber et al. (2016), who used the single-fluid approximation.

Fig. 4. Time–distance plots for: Viz (top-left), Bz (top-right), Viz − Vnz (bottom-left), and Viy − Vny (bottom-right) in the case of A = 10 km⋅s−1 and y0 = 0.3 Mm.

Figure 5 (top panels) presents time–distance plots for the perturbed ion temperature δT_i/T (left panel) and for the vertical component of the ion velocity V_iy (right panel) in the case of A = 10  km⋅ s⁻¹ and y₀ = 0.3 Mm (i.e. more or less in the middle of the photosphere). Some correlation between the velocity and temperature signals is discernible from this analysis. Mainly, in the range of 0 Mm  ≤ y ≤ 2 Mm the signal in δT_i/T experiences a similar trend to the signal in Viy with max(V_iy)≈65 km⋅s⁻¹ and max(δT_i/T)≈1 taking place for the launching level y₀ = 0.3 Mm. In the case of the same launching point y₀, but for a much smaller pulse amplitude, mainly A = 1 km⋅s⁻¹, the maximum value of the vertical component of ion velocity is almost 20 times lower than for this amplitude A = 10 km⋅s⁻¹ pulse. Similarly, the maximum value of the perturbed relative ion temperature is about 25 times lower for the simulation with amplitude A = 1 km⋅s⁻¹. As a matter of fact, for A = 1 km⋅s⁻¹ we obtained max(V_iy)≈3.6 km⋅s⁻¹ and max(δT_i/T)≈0.05 (not shown). The bottom panels of Fig. 5 illustrate the perturbed relative ion temperature δT_i/T averaged over time, and the vertical component of the ion velocity V_iy averaged over time. These two quantities are defined as

(23)

Fig. 5. Time-distance plots for δTi/T (top-left) and Viy (top-right), perturbed temperature of ions δTi/T averaged over time (bottom-left), and the vertical component of ion velocity Viy averaged over time (bottom-right) for y0 = 0.3 Mm and A = 10 km⋅s−1.

and

(24)

where t₁ = 0 s and t₂ = 5000 s. We note that the max(H(y)) takes place for y = 2 Mm (bottom-left). However, a small increase in temperature is discernible at y ≈ 0.3 Mm, the height at which the initial pulse is launched in this case. From the bottom-right panel of Fig. 5, it is inferred that a slow downflow takes place in the lower atmospheric layers up to y ≈ 5 Mm. Higher up, on the other hand, there is an upflow or outflow, with a magnitude that is growing with height y; it reaches about 0.35 km⋅s⁻¹ at y = 20 Mm. Below y ≈ 5 Mm, the downflow reaches a minimum value of V_iy ≈ −0.7 km⋅s⁻¹. This means that when going from the lower atmosphere to the corona, the downflow turns into an upflow. A similar scenario was reported by Kayshap et al. (2015) who found even faster downflows of 0 − 13 km⋅s⁻¹ and also faster upflows of 10 − 12 km⋅s⁻¹. Hence, our simple 2.5D numerical findings reveal only the general trend of the vertical flows that is also present in the observational data. This vertical flow can be generated by ponderomotive force originating from Alfvén waves (e.g., Murawski 1992).

Figure 6 displays the relative perturbed ion temperature averaged over height and time. This quantity is defined as (with H(y) given by Eq. (23)):

(25)

Fig. 6. Relative perturbed ion temperature averaged over time and height H vs. pulse launching y0 for A = 10 km⋅s−1 (top) and vs. pulse amplitude A for y0 = 0.25 Mm (bottom).

where y₁ = 20 Mm. The top panel illustrates H versus the launching height y₀. We note that the minimal heating occurs for y₀ = 0 Mm and the maximal heating H is obtained launching the pulse at y₀ = 0.35 Mm, which is close to the middle of the photosphere. We also note that H strongly falls off with y₀ > 0.35 Mm. Such behavior of H versus y₀ results from the magnetoacoustic waves that release thermal energy as a result of the ion–neutral collisions and which are driven by Alfvén waves. The ion–neutral collisions are less frequent in the upper layers as the ionization degree is close to one (100%) there.

Figure 6 (bottom panel) illustrates the dependence of the average heating H on the pulse amplitude A for a pulse located at y₀ = 0.25 Mm. Here, according to our expectations, it is clearly seen that the higher the amplitude of the pulse, the more heat is deposited in the atmosphere.

4. Summary and conclusion

In this paper, we present the results of 2.5D simulations of impulsively generated two-fluid Alfvén waves. We show that as a consequence of ion–neutral collisions (e.g., Ballester et al. 2018), large amplitude Alfvén waves generated in the photosphere and chromosphere can contribute to the heating of the chromosphere (e.g., Murawski et al. 2020). These waves can also drive plasma outflows, which higher up can become more substantial and constitute the origin of the solar wind. However, the magnitude of the flows is found to be rather moderate and much lower than the observed in- and outflows by Tian et al. (2009) and Kayshap et al. (2015). Hence, we conclude that impulsively generated two-fluid Alfvén waves with an initial amplitude A = 10 km⋅s⁻¹ are not able to explain the observational data, even though the obtained flow amplitudes are in the observed ranges. Nevertheless, the magnitude of the obtained plasma heating is proportional to the pulse amplitude. We find that the maximum heating occurs for the pulse launched from the middle of the photosphere. The maximum of the averaged relative temperature H in the chromosphere increases by about 50%. The pulses launched from the low photosphere and in the chromosphere do not substantially affect the chromosphere temperature (Fig. 6, top). Also, the magnitude of the obtained outflows grows with the amplitude of the pulse. In particular, in the case of the pulse launched from the middle of the photosphere and with an amplitude A ≤ 2.5 km⋅s⁻¹, the averaged heating H is negligibly small. For larger values of the amplitude A, the trend in H(A) is essentially linear with H ≈ 0.42 for A = 10 km⋅s⁻¹ (Fig. 6, bottom).

Our model is based on the two-fluid equations in which the only source of plasma heating are ion–neutral collisions. Hence, more realistic models that include other nonideal and nonadiabatic effects are required to more accurately describe the solar atmosphere and make the simulations more realistic. Such models are left for future studies and the obtained results will be reported elsewhere.

Acknowledgments

The JOANNA code has been developed by Darek Wójcik. This work was done within the framework of the projects from the Polish National Foundation (NCN) grants Nos. 2017/25/B/ST9/00506 and 2020/37/B/ST9/00184. Numerical simulations were performed on the LUNAR cluster at Institute of Mathematics of University of M. Curie-Skłodowska, Lublin, Poland. S. P. acknowledges support from the project (EUHFORIA 2.0) that has received funding from the European Union’s Horizon 2020 research and innovation programme under grant agreement No 870405. SP also received support from the projects C14/19/089 (C1 project Internal Funds KU Leuven), G.0D07.19N (FWO-Vlaanderen), SIDC Data Exploitation (ESA Prodex-12), and the Belspo projects BR/165/A2/CCSOM and B2/191/P1/SWiM.

References

All Figures

	Fig. 1. Variation with height y of the equilibrium temperature T, given by Avrett & Loeser (2008) (top), and Alfvén speed (bottom).
In the text

	Fig. 2. Time–distance plots for Viy (top) and δTi/T (bottom) for the initial pulse-free system, i.e., A = 0 km⋅s−1. In other words, this plot shows the numerical relaxation of the initial equilibrium state.
In the text

	Fig. 3. Time–distance plots for Viz (left) and Bz (right) for A = 1 km⋅s−1, y0 = 0 Mm (top) and y0 = 0.5 Mm (bottom).
In the text

	Fig. 4. Time–distance plots for: Viz (top-left), Bz (top-right), Viz − Vnz (bottom-left), and Viy − Vny (bottom-right) in the case of A = 10 km⋅s−1 and y0 = 0.3 Mm.
In the text

	Fig. 5. Time-distance plots for δTi/T (top-left) and Viy (top-right), perturbed temperature of ions δTi/T averaged over time (bottom-left), and the vertical component of ion velocity Viy averaged over time (bottom-right) for y0 = 0.3 Mm and A = 10 km⋅s−1.
In the text

	Fig. 6. Relative perturbed ion temperature averaged over time and height H vs. pulse launching y0 for A = 10 km⋅s−1 (top) and vs. pulse amplitude A for y0 = 0.25 Mm (bottom).
In the text