Numerical Simulations of Magnetoacoustic-Gravity Waves in the Solar Atmosphere

We investigate the excitation of magnetoacoustic-gravity waves generated from localized pulses in the gas pressure as well as in vertical component of velocity. These pulses are initially launched at the top of the solar photosphere that is permeated by a weak magnetic field. We investigate three different configurations of the background magnetic field lines: horizontal, vertical and oblique to the gravitational force. We numerically model magnetoacoustic-gravity waves by implementing a realistic (VAL-C) model of solar temperature. We solve two-dimensional ideal magnetohydrodynamic equations numerically with the use of the FLASH code to simulate the dynamics of the lower solar atmosphere. The initial pulses result in shocks at higher altitudes. Our numerical simulations reveal that a small-amplitude initial pulse can produce magnetoacoustic-gravity waves, which are later reflected from the transition region due to the large temperature gradient. The atmospheric cavities in the lower solar atmosphere are found to be the ideal places that may act as a resonator for various oscillations, including their trapping and leakage into the higher atmosphere. Our numerical simulations successfully model the excitation of such wave modes, their reflection and trapping, as well as the associated plasma dynamics.


Introduction
The complicated magnetic field configuration of the Sun plays a key role in various types of dynamical plasma processes in its atmosphere, including all the significant plasma dynamics of the lower solar atmosphere. The resulting magnetic structures channel energy from the photosphere into the upper atmosphere, in form of magnetohydrodynamic (MHD) waves, and such waves experience mode conversion, resonances, trapping and reflection which results in the complicated dynamical processes in the lower solar atmosphere, the details of which depend on the plasma properties as well as strength of the magnetic field.
The complex magnetic field and plasma structuring in the lower solar atmosphere support the excitation of various kinds of MHD waves, and their propagation, reflection and trapping has been studied extensively both on theoretical and observational grounds (e.g. Murawski et al., 2011;Srivastava & Dwivedi, 2010;Srivastava, 2010;Fedun et al., 2009;Srivastava et al., 2008;Hasan et al., 2005;McAteer et al., 2003;Gruszecki et al., 2011, and references therein). The evolving magnetic fields of the lower solar atmosphere also lead to transient processes across a wide range of spatial-temporal scales in form of the eruption and associated phenomena. For example, various types of solar jets are formed at short spatial-temporal scales which play a significant role in mass and energy transport and also couple the various layers of the solar atmosphere (e.g. Shibata et al., 2007;Katsukawa et al., 2007;Srivastava & Murawski, 2011, and references therein). In addition, the magnetic activity and injections of helicity into the lower solar atmosphere result in large-scale eruptive phenomena, including solar flares and coronal mass ejections (CMEs) in the outer part of the magnetized solar atmosphere (e.g. Shibata & Magara, 2011;Zhang et al., 2012, and references therein). Therefore, the coupling of the complex magnetic field in various layers of the Sun due to waves and transients is one of the most significant areas of contemporary solar research.
In the quiet-Sun magnetic networks, cavities are important locations where the magnetic fields are sufficiently inclined due to their well evolved horizontal components. They are formed over the granular cells in form of the field-free regions due to the transport of plasma at their boundaries and overlaid by bipolar magnetic canopies. The vertical magnetic fields, however, reside mostly in the core of such magnetic networks (Schrijver & Title, 2003;Centeno et al., 2007). Such magnetic cavity-canopy systems are thought to be ideal resonators of the various MHD waves that can be trapped in the cavity, as well can also leak in form of the magnetoacoustic-gravity waves upward through the core of such magnetic networks. It is thought that the field-free cavity regions underlying the bipolar canopy can trap the high-frequency acoustic oscillations, and the low-frequency components may leak into the higher atmosphere in form of magnetoacoustic-gravity waves Srivastava et al., 2008;Srivastava, 2010). Therefore, such magnetic structures in the lower solar atmosphere are the regions that may play an important role in the wave filtering (e.g. McIntosh & Judge, 2001;Krijger et al., 2001;McAteer et al., 2002;Vecchio et al., 2007Vecchio et al., , 2009Srivastava, 2010, and references therein).
Alongside recent high-resolution observations of MHD waves in the lower solar atmosphere, extensive efforts have been made in the area of analytical and numerical modeling of such waves: Fedun et al. (2009) have investigated the 3D numerical modeling of the coupled slow and fast magnetoacoustic wave propagation in the lower solar atmosphere; and recently Fedun et al. (2011a) have reported the first numerical results of the frequency filtering of torsional Alfvén waves in the chromosphere. Apart from general numerical modeling of the waves in the lower solar atmosphere, models have also investigated the acoustic wave spectrum in the localized magnetic structures of the lower solar atmosphere, e.g., magnetic cavity-canopy system (e.g. Kuridze et al., 2008;Srivastava et al., 2008;Kuridze et al., 2009, and references cited therein).
It is also noteworthy that Bogdan et al. (2003), Fedun et al. (2009) andFedun et al. (2011a) discussed in detail the excitation, propagation and conversion of magnetoacoustic waves in a realistic 3D MHD simulation. However, in these references, the waves driven by a periodic driver were discussed, whereas we numerically simulate the excitation of magnetoacoustic-gravity waves generated by pulses in the gas pressure and vertical component of velocity, mimicking an isolated solar granule. We aim to investigate and understand this simpler (albeit complex enough) system before tackling the more realistic, multiple granule system. Our philosophy is to build up our models incrementally, with a clear focus on the underlying physical processes at each step.
In this paper, we investigate the excitation of magnetoacoustic-gravity waves generated from localized pulses in the gas pressure as well as in vertical component of velocity, modelling the effect of an isolated solar granule. These pulses are initially launched at the top of the solar photosphere that is permeated by a weak magnetic field. We investigate three different configurations of the background magnetic field lines: vertical, horizontal, and oblique to the gravitational force.
We aim to show that small amplitude perturbations that are associated with such a granule are able to trigger large amplitude, complicated oscillations in the solar corona, which exhibit periodicities within the detected range of 3 − 5 min.
The structure of the paper is as follows: In Sect. 2 we describe the numerical model. We report the numerical results in Sect. 3 and present the discussion and conclusions in Sect. 4 .

Numerical model
We consider a gravitationally-stratified solar atmosphere that is described by the ideal two-dimensional (2D) MHD equations: Here ̺ is mass density, V is the flow velocity, B is the magnetic field, p = k B ̺T /m is the gas pressure, T is the temperature, γ = 5/3 is the adiabatic index, g = (0, −g, 0) is the gravitational acceleration, where g = 274 m s −2 , m is the mean particle mass and k B is Boltzmann's constant. Throughout this paper, we use the Cartesian coordinate system with the vertical axis denoted by y and the horizontal axis x. Henceforth, we assume that the medium is invariant along the z-direction with ∂/∂z = 0 and set the z-components of plasma velocity and magnetic field equal to zero, i.e. V z = 0 and B z = 0. The latter assumption removes Alfvén waves from the system but still allows magnetoacoustic-gravity waves to propagate freely.

Equilibrium configuration
We assume that at its equilibrium the solar atmosphere is static (V e = 0) and threaded by a straight magnetic field, whereŝ is a unit vector that is either vertical, horizontal or oblique. We choose B 0 by requiring that at the reference point (0, 10) Mm, the Alfvén speed and sound speed satisfy the following constraint: c A (y = 10 Mm) = 10 c s (y = 10 Mm). This constraint reproduces typical conditions in the solar corona where, typically, c s = 0.1 Mm s −1 and c A = 1 Mm s −1 . As a result, the solar corona is magnetically dominated with B 0 ≃ 14.5 × 10 −4 T. The plasma β, in the solar corona attains a value of β(y = 10 Mm) = 0.012. It grows slowly with depth within the chromosphere and below y = 2 Mm reaches abruptly a value of β(y = 0 Mm) ≃ 4 × 10 5 at the bottom of the photosphere (Fig. 1, bottom).
This large value of β evidences that in these low regions of the solar atmosphere the effect of magnetic field is negligibly small. The straight magnetic field of Eq. (5) is a simplified model of curved magnetic field lines which are located in low regions of the solar atmosphere. However, as the plasma β grows fast with depth below y = 2.5 Mm reaching large values there (see Fig. 1, bottom) we actually model the very weakly-magnetized low layers of the atmosphere, which do not correspond to flux-tubes and sunspots.
Therefore, in this first approximation, the curved magnetic field lines can be replaced by a straight magnetic field, and our models are justified.
In Eqs. (6) and (7), ̺ e (y) and p e (y) denote equilibrium mass density and gas pressure, respectively. They are specified by the hydrostatic constraint which in the context of Eq. (5) states that the pressure gradient is balanced by the gravity force, − ∇p e + ̺ e g = 0 .
With the ideal gas law and the y-component of Eq. (9), we arrive at where is the pressure scale-height, and p 0 denotes the gas pressure at the reference level that we choose in the solar corona at y r = 10 Mm.
We adopt an equilibrium temperature profile T e (z) for the solar atmosphere that is close to the VAL-C atmospheric model of Vernazza et al. (1981): see Higher up, T e (y) grows gradually with height up to the transition region which is located at y ≃ 2.7 Mm. Here T e (y) experiences a sudden growth up to the coronal value of 1.5 MK at y = 10 Mm. Then with Eq. (10) we obtain the corresponding gas pressure and mass density profiles. Both p e (y) (not shown) and ̺ e (y) (Fig. 1, right-top panel) experience a sudden fall-off from photosphere to coronal values at the transition region.

Initial Conditions
At t = 0 s, we initially perturb the equilibrium impulsively using localized Gaussian pulses simultaneously in the gas pressure and vertical component of velocity, viz., Here A p and A v are the amplitudes of the perturbations, (x 0 , y 0 ) is their initial position and (w x , w y ) denotes their widths along the xand y-directions. We set and hold fixed x 0 = 0 Mm, y 0 = 0.5 Mm, w x = 0.5 Mm, w y = 0.5 Mm, A p = 0.02 p e (y 0 ) and A v = 0.2 × 10 −3 Mm s −1 . These magnitudes of A v and A p represent the recently detected flow and temperature in solar granulation (Baran, 2011). We separately consider three orientations for our unidirectional equilibrium magnetic field lines: (a) horizontal, (b) vertical, (c) oblique to the gravitational force.

Numerical results
Equations (1)-(4) are solved numerically using the FLASH code (Lee & Deane, 2009).This code implements a second-order unsplit Godunov solver with various slope limiters and Riemann solvers, as well as adaptive mesh refinement (AMR). The main advantage of using AMR technique is to refine a numerical grid at steep spatial profiles while keeping a grid coarse at the places where fine spatial resolution is not essential. Such AMR technique usually introduces interpolation errors at different-size numerical cells. These errors can result in vertical flow that, albeit initially small, can grow with height to an unacceptable magnitude.
A remedy of this inherent phenomenon is to refine the whole region below the transition region, which was adopted in our simulations. We used the Roe solver and minmod flux limiter, and set the simulation box for the horizontal and vertical equilibrium magnetic field cases in the x-direction as −5 ≤ x ≤ 5 Mm. In order to trace plasma structures which extend more horizontally in the case of the oblique magnetic field we use −2 ≤ x ≤ 8 Mm in this case. Along the y-direction the   shock (and subsequent shocks which arrive respectively after about 160 s and 240 s) results from the nonlinear wake, which lags behind the leading shock (Sterling & Hollweg, 1989;. These times can be compared with the acoustic cut-off period (Roberts, 2006), which for y = 0.5 Mm attains a value of P ac (y) ≃ 180 s. This value differs by 70 s from the time-span between arrivals of neighboring shocks, which is 250 s. However, the 2D model we discuss here is more complex than the 1D scenario which is described by the Klein-Gordon equation (Roberts, 2006).The wave-period we detected is altered by the interaction between up-going waves from the launching place of the initial pulse and the reflected one from the transition region signals. As wave reflections result at the large temperature gradients, therefore, the chromosphere sustains a cavity for these waves which are represented by the oblique stripes located at altitudes 0.7 Mm < y < 1.4 Mm (Fig. 5). In the neighborhood of the point x = 0.4 Mm, y = 0.4 Mm, we observe the formation of vortices (Fig. 5, top panel) which experience energy cascade into smaller scales. These vortices are present till the end of our simulation runs (Fig. 5, bottom panel). The first vortex results from the initial pulse in V y , which is a characteristic feature of velocity perturbations. This first vortex seeds convection in the convectively unstable plasma layers. According to the Schwarzschild's instability condition, a medium is convectively unstable if the squared buoyancy (or Brunt-Väisälä) frequency, is negative (e.g., Roberts, 2006). As this criterion is satisfied for y < 0.6 Mm, therefore, the convection sets in there. Such vortices were also theorized by Konkol, Murawski, & Zaqarashvili (2011)in a similar context.

Horizontal equilibrium magnetic field:ŝ =x
In this section, we investigate the horizontal equilibrium magnetic field, which corresponds toŝ =x in Eq. (5). In this case, the spatial profiles of T drawn The corresponding movie can be found in the file fig6.avi online.
at t = 250 s and t = 10 3 s (Fig. 6) reveal small amplitude oscillations of the transition region without the presence of a jet which was observed for the case of the vertical magnetic field. On comparing with Fig. 3, the waves resulting from the initial pulses already arrived to the solar corona at t = 250 s (Fig. 6, left panel). As the magnetic field lines are horizontal, therefore, these waves are essentially fast magnetoacoustic-gravity waves in the region of low plasma β that occurs for y > 2.6 Mm (Fig. 1, bottom panel). It is interesting that the orientation of magnetic field plays very crucial role in determining the amplitude of oscillations of the transition region.
Similar to the time-signatures of Fig. 4, V y collected at the detection point (x = 0, y = 5) Mm, reveals the shocks. However, for the horizontal magnetic field, we now observe only two shocks (as seen in Fig. 7). The first shock arrives to the detection point at t = 250 s and the second shock reaches this point at  3.3. Oblique equilibrium magnetic field:ŝ = (x +ŷ) / √ 2 In this section, we discuss the case of the oblique magnetic field, for whicĥ s in Eq. (5) is at an angle of π/4 to the horizontal axis.  (Figure 9, bottom-right), we observe surface waves propagating along the transition region. However, for the oblique equilibrium magnetic field, we find that time signatures of V y reveal shocks, and the waveperiods are about 250 − 300 s (Fig. 10). The jet, which is associated with the shock, lasts for several hundred seconds and its velocity reaches up to maximum limit of ∼ 30 × 10 −3 Mm s −1 . This is a result of complex interaction of waves with the background plasma in upper regions of the solar atmosphere.
This complexity is seen at t = 10 3 s in the velocity profiles of Fig. 11, which clearly illustrates the trapped and reflected waves in the solar chromosphere as well as well-developed vortices.

Discussion and conclusion
We have investigated the impulsive excitation of magnetoacoustic-gravity oscillations and have compared and contrasted their resultant propagation under different orientations of equilibrium magnetic fields. We performed 2D numerical simulation of the velocity and gas pressure pulses, which mimic a solar granule.
These pulses were initially launched at the top of solar photosphere, in a stratified solar atmosphere utilizing the VAL-C temperature profile. We find that in the cases the background magnetic field possesses a non-zero vertical component the amplitude of the upwardly-propagating perturbations rapidly grows with height due to the rapid decrease in the equilibrium mass density. Therefore, the perturbation quickly steepens into shocks in the upper regions of the solar chromosphere that launches the cool material behind it.
We find that the excitation of magnetoacoustic-gravity waves in a non-horizontal equilibrium magnetic field configuration is channeled upwards and is able to cause large amplitude transition region oscillations. Meanwhile, the horizontal equilibrium magnetic field configuration results in comparatively smaller amplitude transition region oscillations, due to fast magnetoacoustic-gravity waves which spread their energy across magnetic field lines in contrast to slow magnetoacoustic-gravity waves which in a strongly magnetized plasma region   are guided along magnetic field lines. This indicates that the wave energy is transfered into the upper atmosphere in larger amounts where the magnetic field is more vertical. However, the complexity and therefore the evolution of a horizontal equilibrium magnetic field allows the reflection and trapping of the waves in the lower solar atmosphere and can influence the localized dynamics and heating of the atmosphere. The reflection of the waves in the lower solar atmosphere and its trapping may also result in the presence of chromospheric cavities.
By analyzing the time-signature of V y collected at a fixed spatial point, we also find that the period of oscillation is longer when comparing the vertical magnetic field case (here, the characteristic period of oscillation was ≈ 250 − 300 s) with the horizontal equilibrium magnetic system (in which the characteristic period of oscillation was ≈ 150 − 200 s). For the vertical magnetic field system in low plasma β regions, the magnetoacoustic-gravity waves are well-described as slow magnetoacoustic-gravity waves (propagating along the magnetic field lines at approximately the sound speed, c s (y), see Eq. 7). For the horizontal magnetic field system in the strongly magnetized plasma, the oscillations transverse to the magnetic field are well-described as fast magnetoacousticgravity waves (propagating across the magnetic field lines at the fast speed, c 2 f (y) = c 2 A (y) + c 2 s (y), see Eqs. 6 and 7). As for y = 0.5 Mm we have c A ≪ c s , therefore the cut-off frequencies of the fast and slow magnetoacoustic-gravity waves are very close to each other, which results in very close values of the cutoff waveperiods. Therefore, we could expect that the detected waveperiods for the vertical and horizontal magnetic fields exhibit similar values. However, the cut-off waveperiods are derived on the basis of the linear theory which is valid for small amplitude oscillations only. As such small oscillations are present in the case of the horizontal magnetic field, therefore, the obtained numerical data lies close to the analytical prediction. In the case of the vertical magnetic field, the upwardly propagating waves interact with the waves which become reflected from the transition region. As a result of that larger amplitude oscillations originate, which significantly alter the background plasma. The waves, which propagate through this strongly modified medium exhibit modified velocities and they get reflected from the transition region that is locally largely curved.
As a consequence of that waveperiods within the range 250 − 300 s result, which differ from P ac .
It is known that radiation is an effective mechanism of wave damping in the low photosphere ( (Mihalas & Toomre, 1982).It might not radically alter the system's behavior, but radiative damping is at least likely to reduce the amplitude of the waves reaching the transition region, leading to shorter jets than what we see in our models, as well as lower amplitudes of the coronal shocks. In this paper, we do not invoke the radiative cooling nor thermal conduction in our model atmosphere as we aim to model the small-scale atmospheric regions above the solar photosphere where such effects are not believed to be dominant. However, we intend to include these effects in our future studies.
It is noteworthy that the solar atmosphere is structured by convective overshoot which is absent in the model we devised. Instead, we isolated a single granule-like perturbation by aiming to mimic a magnitude of flow and plasma temperature that are associated with the solar granulation (Baran 2011), and considered the complex scenario which results by this simple model to understand the physics of wave phenomena above such magnetic structures in the solar atmosphere.
We intend to develop more advanced models in our future studies.
In conclusion, our numerical simulations clearly demonstrate that small amplitude initial pulses in vertical velocity and gas pressure are able to trigger a plethora of dynamic phenomena in the upper regions of the solar atmosphere with waveperiods within the range of 150 − 300 s, a value which depends on orientation of the background magnetic field. However, it should be noted that the performed 2D simulations are idealized in the sense that they do not include radiative transfer and thermal conduction along field lines. The magnetic field configuration and the equilibrium stratification are simple and we modeled a single granule only. These limitations require additional studies which we intend to carry on in near future.