On the Generation of Hydrodynamic Shocks by Mixed Beams and Occurrence of Sunquakes in Flares

Observations of solar flares with sunquakes by space- and ground-based instruments reveal essentially different dynamics of seismic events in different flares. Some sunquakes are found to be closely associated with the locations of hard X-ray (HXR) and white-light (WL) emission, while others are located outside either of them. In this article we investigate possible sources causing a seismic response in a form of hydrodynamic shocks produced by the injection of mixed (electron plus proton) beams, discuss the velocities of these shocks, and the depths where they deposit the bulk of their energy and momentum. The simulation of hydrodynamic shocks in flaring atmospheres induced by electron-rich and proton-rich beams reveals that the linear depth of the shock termination is shifted beneath the level of the quiet solar photosphere on a distance from 200 to 5000 km. The parameters of these atmospheric hydrodynamic shocks are used as initial condition for another hydrodynamic model developed for acoustic-wave propagation in the solar interior (Zharkov, Mon. Not. Roy. Astron. Soc.431, 3414, 2013). The model reveals that the depth of energy and momentum deposition by the atmospheric shocks strongly affects the propagation velocity of the acoustic-wave packet in the interior. The locations of the first bounces from the photosphere of acoustic waves generated in the vicinity of a flare are seen as ripples on the solar surface, or sunquakes. Mixed proton-dominated beams are found to produce a strong supersonic shock at depths 200 – 300 km under the level of the quiet-Sun photosphere and in this way produce well-observable acoustic waves, while electron-dominated beams create a slightly supersonic shock propagating down to 5000 km under the photosphere. This shock can only generate acoustic waves at the top layers beneath the photosphere since the shock velocity very quickly drops below the local sound speed. The distance Δ\documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$\Delta$\end{document} of the first bounce of the generated acoustic waves is discussed in relation to the minimal phase velocities of wave packets defined by the acoustic cutoff frequency and the parameters of atmospheric shock termination beneath the photosphere.


Introduction
During the past decade, it became well established that flares can significantly affect the solar interior, as first predicted by Wolff (1972). The first observation of this impact was made with the Michelson Doppler Imager (MDI) onboard the Solar and Heliospheric Observatory (SOHO) and was seen as ripples travelling from a flare location called a 'sunquake'. This was reported by Kosovichev and Zharkova (1998) during the flare of 9 July 1996. The original explanation for sunquakes suggested that these ripples are reflections from the solar surface of acoustic waves induced during flares by some sharp deposition of momentum (and energy) (Kosovichev and Zharkova, 1995).
The most likely agents delivering this momentum are considered to be hydrodynamic shocks, which arise in the chromosphere as a result of pressure transients driven by a hydrodynamic response of the ambient plasma to the precipitation of energetic electrons into the chromosphere Zharkova, 1995, 1998;Donea and Lindsey, 2005;Kosovichev, 2006). Numerous hydrodynamic simulations have shown that the injection of energetic-particle beams leads to fast heating of the ambient plasma in the corona and chromosphere and pushes it towards the photospheric levels (Somov, Spektor, and Syrovatskii, 1981;Fisher, Canfield, and McClymont, 1985c). This process causes either mild (Allred et al., 2005) or explosive evaporation into the corona (Somov, Spektor, and Syrovatskii, 1981;Fisher, Canfield, and McClymont, 1985c;Zharkova and Zharkov, 2007) of the sweptup chromospheric plasma, depending on the intensity of the electron beams, at the same time producing low-temperature hydrodynamic shocks moving towards the photosphere and lasting up to several minutes (Zharkova and Zharkov, 2007).
However, the first sunquake (Kosovichev and Zharkova, 1998) indicated essential discrepancies in both timing and energy of the deposited disturbances (momentum) recorded from observation of the sunquake and those delivered by a shock derived from the standard model of a flaring atmosphere heated by electron beams. This provided an additional stimulus to investigate the physics of solar flares and the redistribution of energy from the primary energy release site to the lower atmospheric levels and the solar interior. In turn, it led to the introduction of protons as the agents causing hydrodynamic shocks with the required momentum and travel time (Zharkova and Zharkov, 2007). Eleven locations were detected that showed evidence of hydrodynamic shocks in Dopplergrams of the flare of 28 October 2003, which had termination downward velocities from 1.0 to 2.1 km s −1 , each lasting up to several minutes. However, only three of them with termination downward velocities above 1.5 km s −1 produced detectable sunquakes on the solar surface (Zharkova and Zharkov, 2007).
Since the first identification of a sunquake by time-distance methods, the development of the helioseismic holography method led to the identification of many more seismic sources during flares of class M and X observed with SOHO/MDI and the Solar Dynamic Observatory (SDO) instruments, with the total number growing towards a hundred (see for example Donea et al., 2006;Zharkova and Zharkov, 2007;Martínez-Oliveros, Moradi, and Donea, 2008;Zharkova, 2008;Matthews, Zharkov, and Zharkova, 2011;Donea, 2011;Alvarado-Gomez et al., 2012;Martínez-Oliveros et al., 2012). An increasing number of sunquakes combining observations in different ranges of the electromagnetic spectrum, from high-energy hard X-ray (HXR) and γ -ray emission to lower energy optical and infrared radiation, provides unique opportunities to study the excitation of acoustic oscillations induced by flares in detail. This will lead to a better understanding of various physical conditions in flares and of the role of different agents that are capable of causing sunquakes.
In many cases the occurrence of sunquakes is accompanied by strong high-energy emission in HXR, γ -rays, and in the extreme ultra-violet (EUV), which some authors suggest as the alternative source of pressure transients causing sunquakes. This aspect is related to backwarming of the photosphere by the enhanced coronal and chromospheric radiation, leading to white-light and seismic emission (Lindsey and Braun, 2000;Donea, Braun, and Lindsey, 1999;Donea and Lindsey, 2005). Some sunquakes are found in locations with little or no HXR and white-light emission , which indicates that at least for some flares, backwarming heating is not essential.
The shocks formed by a hydrodynamic response of the ambient plasma to precipitation of electron or proton beams appear to be capable of causing seismic events associated with solar flares because they carry sufficient mass and momentum to the solar photosphere (Kosovichev and Zharkova, 1998;Zharkova and Zharkov, 2007). However, it is still a question how exactly the energy deposition by these shocks into the solar interior takes place (depths, timescale) and how it compares with observations. These are closely related to the other implications of particle precipitation into a flaring loop, such as the formation of magnetic transients and non-thermal plasma ionization leading to prolonged white-light emission from solar flares (Zharkova, 2008).
In addition, there are essential changes to the photospheric magnetic field often recorded during a flare well above the noise levels, on both short timescales of minutes (i.e. the duration of the impulsive phase) and longer timescales of hours (before and after a flare) Zharkova, 1999, 2001;Sudol and Harvey, 2005;Petrie and Sudol, 2010;Petrie, 2012Petrie, , 2013. The magnetic field changes in flare locations show magnitudes of 50 -120 G that are well above the background noise of any instruments used, as was confirmed by the investigation of magnetic changes in general, and variations of a magnetic inversion line in particular .
The first type of changes, transient (or V-type) changes to the line-of-sight (LOS) magnetic field, occur in the small pixel-size locations of loop footpoints during the impulsive phase of flares , and thus appear to be closely related to the initial energy release (Zharkova et al., 2011b), often being close in time and location to the appearance of enhanced WL continuum and HXR emission. The second type, the step-type changes, of the LOS magnetic field changes occur over larger areas including a magnetic inversion line (MIL), or polarity inversion line (PIL) Petrie and Sudol, 2010). These changes appear just before the impulsive phase and continue for the whole duration of the flares (a timescale of one or two hours) .
The second type of changes has been shown to be related to the whole energy deposited during a flare Petrie and Sudol, 2010) and is most likely associated with magnetic reconnection producing the flare. The V-type magnetic field changes, or magnetic transients, can appear in small areas of flares associated with strong HXR emission as a result of the changes in the magnetic field due to precipitation of particle beams (e.g. Zharkova and Kosovichev, 2002;. These beams carry a strong electric field (Zharkova and Gordovskyy, 2006;Siversky and Zharkova, 2009b), which in turn induces an electromagnetic field in the ambient plasma (van den Oord, 1990). This modifies the resulting magnetic field of the loop where the particles precipitate, inducing the magnetic field seen as a transient magnetic field (Zharkova et al., 2011b).
Sunquake locations often coincide spatially with one of these magnetic field changes (Zharkova, 2008;Donea, 2011). This led some researchers to suggest that the magnetic changes can directly induce magneto-acoustic waves in the solar interior that are seen on the solar surface as sunquakes (Hudson, Fisher, and Welsch, 2008;Fisher et al., 2012). However, a close investigation of magnetic field changes in some recent flares did not provide confirmation of a direct association between sunquakes with step-type changes of the magnetic field (Martínez-Oliveros and Donea, 2009;Matthews, Zharkov, and Zharkova, 2011;. Similar to Zharkova (2008), Lindsey and Donea (2008) highlighted that in the case of magnetic changes, the transient component is the most relevant to acoustic emission. These changes are shown to occur on a timescale of τ ≈ 2H/c or less, where H is the density scale height and c is the sound speed. In the photosphere, this timescale is about 40 s. It is still debated whether or not these transient changes, which take place during the impulsive phase, represent genuine change in the photospheric magnetic field. These localized sign reversals, or 'magnetic anomalies' (Patterson and Zirin, 1981), are often regarded as being caused the sudden heating of the ambient plasma or the direct bombardment by high-energy particles, which lead to changes in the line profile.
However, in the case of the Ni 6768 Å line used by the instruments of the Global Oscillation Network Group (GONG) and MDI onboard the SOHO satellite to make magnetic measurements, the non-LTE simulations have shown that thermal heating occurring at the start of the impulsive phase is insufficient to turn this line into emission when the transient magnetic changes are recorded. Instead, it can be caused by a strong increase of non-thermal excitation and ionization caused by precipitating electrons, i.e. intense particle bombardment (Zharkova and Kosovichev, 2002).
Good correspondence between HXR sources and magnetic anomalies was found by Qiu and Gary (2003), providing additional support to the hypothesis that magnetic field changes are associated with energetic particles. Further evidence for that was found in the simulations of GONG and MDI observations by Edelman et al. (2004), who showed that magnetic measurements are less sensitive to the changes in the line profile than Doppler measurements.
Similar conclusions can be drawn for the iron line Fe I 6173 Å, which is used for measuring the magnetic field in the Helioseismic and Magnetic Imager (HMI) instrument onboard the Solar Dynamic Observatory (SDO) satellite following the comparative investigation of this line by the HMI instrument team and the Ni line used for SOHO/MDI (Norton et al., 2006; see also http://hmi.stanford.edu/doc/Tech_Notes/HMITN-2003-002Line_Choice/ line6173_6768_V1.4.pdf).
Further investigation is required into the processes of energy transfer from the magnetic field, which is reconnected during solar flares to energetic particles and by these particles to the flaring atmospheres. In this article we investigate the role of mixed beams (electrons and protons) in the heating of the ambient plasma and formation of hydrodynamic shocks, which can cause sunquakes. In Section 2 we discuss the hydrodynamics of the ambient plasma, in Section 3 we probe the shock parameters with the conditions required for detection of ripples on the solar surface by the theory of acoustic-wave propagation in the solar interior, and we conclude in Section 4.

Magnetic Reconnection and Particle Acceleration
The step-type changes occurring during flares in a LOS magnetic field are likely to reflect the process of 3D magnetic reconnection, the primary energy release mechanism initiating flares and leading to a change of the magnetic field directivity and/or shear seen in the abrupt changes of a LOS magnitude (Somov, 2000;Priest and Forbes, 2002;Petrie and Sudol, 2010;Petrie, 2012Petrie, , 2013Wang et al., 2012).
These changes occur over a large area of the whole active region including a magnetic inversion line (MIL) separating the LOS magnetic field with opposite polarity signs. The MILs were consistently derived from the magnetograms using the Gaussian gradient  in the vicinity of a few different flares, e.g. 14 July 2000 or 23 July 2002 Zharkova, 1999, 2001;.  have also shown that the step-type temporal variations of the magnetic field correlate rather closely with the variations of HXR emission measured during the flare of 23 July 2002. However, the step-type magnetic field changes only occurred at the locations of some of the HXR sources, while in the others only the transient changes were registered in the MDI magnetograms Martínez-Oliveros and Donea, 2009;Matthews, Zharkov, and Zharkova, 2011); these are also clearly seen in some of the recent HMI observations by Petrie and Sudol (2010) and Petrie (2012).
The difference in magnetic field magnitude between the step positions is found to account for the total energy deposited during a flare . In fact, the estimations show that the magnetic energy calculated from the difference of step-type changes can be two or three times higher than the total energy of HXR, soft X-ray, and EUV emissions added together (Holman et al., 2003;Petrie, 2012Petrie, , 2013. In this sense, all the energy deposited during a flare comes from a magnetic reconnection process, or Lorentz force (of step-type changes), as suggested by  and Hudson, Fisher, and Welsch (2008). These step-type changes were recently measured for different flares by Petrie and Sudol (2010), Petrie (2012Petrie ( , 2013 and Wang et al. (2012) and were quantified with an MHD model (Fisher et al., 2012). The latter simulations are in line with the previous 2D and 3D simulations of magnetic reconnection carried out by various authors (see for example books by Somov, 2000;Priest and Forbes, 2002, and references therein).
Part of the energy released by magnetic reconnection is converted into the energy of energetic particles. Particle acceleration occurs in current sheets at the intersection of interacting loops in a diffusive region where magnetic field lines reconnect. Particles gain their energy from a reconnection electric field directed along the guiding field (Litvinenko, 1996;Zharkova and Gordovskyy, 2004;Zharkova and Agapitov, 2009) while travelling through a 3D magnetic field of a current sheet. During their passage, particles are shown to gain subrelativistic energies with power-law distributions Zharkova and Agapitov, 2009). Hence, the particle energy gains are strongly affected by the magnetic field topology in the reconnecting current sheet, e.g. the longer the time particles spend inside a current sheet, the higher the energy they gain (Zharkova et al., 2011a).
Particles with opposite charges, protons and electrons, are simultaneously accelerated as they are dragged from the neutral ambient plasma into a diffusion region during the process of magnetic reconnection. The particles undergo a complex energisation process at the current sheet midplane until they gain sufficient energy to escape from a given 3D magnetic topology. Electrons are accelerated very quickly (within 10 −6 s), so that they can break quickly from a current sheet and precipitate into a loop leg. Proton acceleration is slower by three orders of magnitude (about 10 −3 s), causing the protons to gyrate much longer about the midplane than electrons. They gain energy at the midplane before they are eventually ejected. Very often protons are ejected into a loop leg with opposite polarity to the leg where electrons are ejected (Zharkova and Gordovskyy, 2004).
In the magnetic topologies with weaker guiding fields, accelerated electrons cannot fully break from a current sheet. Instead, they continue to return to its midplane and create electron clouds, which are probably observed as coronal sources (Zharkova and Siversky, 2011). This happens because the positive charge from protons, which undergo a much longer acceleration from gyrating about the midplane, creates an electrostatic electric field, which returns electrons to the diffusion region. Only after the protons gain sufficient energy to break from the current sheet are they ejected from this current sheet together with the electrons. In this case, the ejection would therefore contain mixed beams of protons and electrons in each leg. The ratio of proton-to-electron numbers has been shown to be dependent on the magnetic field topology in the current sheet (Zharkova and Gordovskyy, 2004;Zharkova and Agapitov, 2009).
Moreover, a full kinetic particle-in-cell (PIC) approach for current sheets in the corona (Siversky and Zharkova, 2009a) or in the heliosphere (the heliospheric current sheet, HCS) (Zharkova and Khabarova, 2012) indicates that full separation of the particles are the rare events applicable to the topologies with large magnitudes of the guiding field. For many other magnetic field topologies this separation is only partial. In the other words, there can be a larger number of electrons than protons in one loop leg and a larger number of protons than electrons in the other leg.
It is evidently rather difficult to distinguish signatures of accelerated particles from the acceleration site in the corona from those modified during their precipitation into loop footpoints (or loop legs). In the HCS, in contrast, it is possible to obtain direct in-situ measurements of accelerated particle parameters. These measurements show a close agreement of the observed parameters (ion and electron density, velocity profiles, pitch angle distributions) of solar wind particles crossing the sector boundary of the HCS with those derived from the PIC simulations carried out for low-density current sheets using a 3D magnetic topology similar to the HCS (Zharkova and Khabarova, 2012). These direct in-situ measurements confirm the PIC simulation results in the HCS magnetic field topology and strongly support the results derived from PIC simulations for current sheet conditions of the solar corona that are suitable for flares (Siversky and Zharkova, 2009a).
Particles of the ambient plasma pass through a current sheet during reconnection and gain subrelativistic energies and negative power-law energy distributions, with spectral indices being essentially different for electrons and protons or ions. Hence, there are four particle populations that appear during a magnetic reconnection process : (1) power-law subrelativistic electrons, (2) power-law subrelativistic protons (produced by direct acceleration in a current sheet) combined with (3) quasi-thermal (up to 10 -20 keV) electrons, and (4) quasi-thermal protons (produced by Alfvén and other waves in the current sheet exhausts). The third and fourth populations have energies much higher than the thermal energies gained from a magnetic reconnection process (Somov, 2000;Gordovskyy et al., 2005). All four types of particles gain subrelativistic velocities from magnetic field reconnection. Together with high velocities, these particles also gain high momenta Zharkova and Zharkov, 2007), which they deposit into the footpoints of a flaring loop.
It is therefore evident that magnetic reconnection energy is delivered to flaring atmospheres by mixed electron and proton beams, which most likely become only partially separated at the current sheet they pass and are injected into the loop legs. In addition, these mixed beams can be accompanied by quasi-thermal protons or electrons appearing because of Alfvén waves generated at the current sheet exhausts.
As an example, there can be 70 % electrons and 30 % protons (in particle numbers) ejected into one leg and 70 % protons and 30 % electrons into the other one. The protons in the first leg have much lower velocities than the electrons, and their energy is not high enough for them to be detected in γ -ray emission. In the second leg, the electrons have gained lower energies, and their number is rather small, therefore they produce weak HXR emission. The charge neutralisation process will include ambipolar diffusion in both cases and not only a formation of return currents from the ambient plasma and lower energy electrons in the beam, as in the case of pure electron beams (Zharkova and Gordovskyy, 2006).

Transient Changes of Magnetic Field
We now evaluate the appearance of transient magnetic fields induced by energetic electron beams precipitating into the footpoints of the flaring atmospheres. These particles deposit their energy by collisions and Ohmic heating of the ambient particles, causing a hydrodynamic response of the ambient plasma that we discuss in Section 2.2.
It is evident that electron beams should carry a strong electrostatic electric field that sometimes exceeds the local Dreicer field by two orders of magnitude (Zharkova and Gordovskyy, 2006). This electric field strongly affects the spectral indices of HXR bremsstrahlung emission in solar flares (Siversky and Zharkova, 2009a). At the same time, these beam electrons are turned back to the corona by their own electrostatic electric field.
This occurs at stopping depths of the electrons with a lower cutoff energy where they start travelling across the magnetic field of a loop (Siversky and Zharkova, 2009a). Then, following Faraday's law, these electrons can induce a strong magnetic field in the ambient plasma, which we have called a transient magnetic field in the Introduction (van den Oord, 1990). This field can exist as long as the beam electrons are injected into the atmosphere, as shown by the time-dependent Fokker-Planck simulations (Siversky and Zharkova, 2009a). This transient magnetic field is exactly what is observed in many flares (Sudol and Harvey, 2005;Petrie and Sudol, 2010), indicating the precipitation of electron beams.
From the difference in photon spectral indices at higher and lower energies, the induced electric field can be derived as (Zharkova and Gordovskyy, 2006). From Faraday's law it then follows that 1 where c is the speed of light. To produce an induced magnetic field, the electrons therefore need to move across the existing magnetic field. This is easily achievable in the case of returning electrons with a lower cutoff energy, which move across the loop at the point where they are turned back to the source by their own electrostatic electric field, as shown by Fokker-Planck numerical simulations (Siversky and Zharkova, 2009b).
This induced transient magnetic field is of opposite polarity to that originally present in a given footpoint. The magnitude of this transient magnetic field can be in the range of 30 -120 G as measured from the magnetic field variations in flare locations . The energy released in this process is linked to the energy of energetic particles, heating the ambient plasma, which we discuss in Section 2.2.

The Radiative Processes and Backwarming Heating
Since energetic particles produce high-energy HXR bremsstrahlung and γ -ray emission, these emissions can play some role in the additional heating of lower atmospheric levels (called backwarming heating), leading to white-light flares (Hudson, 1972;Machado, Emslie, and Avrett, 1989;Metcalf et al., 2003;Donea et al., 2006;Zharkova, 2008), which often accompany seismic emission in flares. A substantial fraction of the continuum emission of white-light flares comes from the overheated upper atmosphere, which can result in radiative backwarming of the photosphere by recombination radiation (Hudson, 1972;Machado, Emslie, and Avrett, 1989). This was simulated using radiative hydrodynamics by Allred et al. (2005).
The photosphere is assumed to absorb the part of coronal and chromospheric radiation that is emitted downward. Therefore, it becomes heated not only by collisions of beam electrons or protons with the ambient particles, but also by the HXR and ultra-violet (UV) or EUV radiation coming from the upper layers. The immediate effect of this absorption in the visible spectrum is mostly dissociation of neutral hydrogen, leading to white-light emission caused by long-lasting recombination of free electrons into neutral hydrogen atoms and creation of H − ions that define the photospheric opacity. This process lasts until the conditions reach the pre-heated level, after which the heating is stopped. The radiation energy estimated from white-light flares is found to be close to that measured using seismic emission (Donea et al., 2006). This allowed Donea et al. (2006) to suggest backwarming as an additional source of delivering energy to the photosphere and led to the assumption that the acoustic emission is closely associated with the continuum emission.
However, Allred et al. (2005) have pointed out that in their hydrodynamic models, which considered heating by electron beams of a very moderate power, the backwarming heating is not very effective during the impulsive phase when the hydrogen emission is formed mostly by non-thermal ionization by beam electrons (Zharkova and Kobylinskii, 1993). Moreover, it was recently shown that the location of white-light emission in a limb flare is very close to the location of low-energy (< 100 keV) HXR emission . This significantly reduces the likelihood that this WL emission will affect the seismic response of the flare, which occurs at much deeper atmospheric levels.
After beam heating is stopped (in about 10 -30 seconds), the thermo-conductivity continues to heat the lower atmosphere further, thereby increasing its emission in all lines and continua, which leads to further increases of the line and continuum opacity. This was reported by Zharkova and Kobylinskii (1993) and has been confirmed by Allred et al. (2005). Moreover, Allred et al. (2005) reported that any additional heating of the ambient plasma does not increase the backwarming heating of the photosphere, which remains modest since it is governed by the high opacity of the continuum radiation.
On the other hand, Zharkova (2008) pointed out that the hydrogen ionization degree is increased by up to six orders of magnitude by collisions with beam electrons, approaching nearly full ionization level, while the kinetic temperature of hydrogen atoms remains as low as it should be at the photosphere. The increased ionization caused by collisions with beam electrons produces strong emission in hydrogen lines and continua (Zharkova and Kobylinskii, 1993;Zharkova, 2008). This means that the opacity of continuum radiation is increased, which in turn keeps this radiation locked in the atmosphere for a timescale comparable with the difference (reaching sometimes tens of minutes) between the impact excitation rates of hydrogen atoms by beam electrons and the recombination rates of the ambient electrons with these hydrogen atoms, as discussed by Zharkova (2008).
Hence, the presence of white-light emission does not necessarily mean that there is backwarming heating of the plasma. It might be just a sign of over-ionization of the ambient plasma by subrelativistic particles (electrons and protons) (Zharkova, 2008).

Delivered Momenta and Timing of Seismic Signatures
The detailed evaluation of the momenta required to produce the sunquakes observed in the flare of 28 October 2003 (Zharkova and Zharkov, 2007) has shown that 11 shocks are detected in the Dopplergram datacubes taken around the flare location during the event. The shocks lasted 1.5 -2.5 minutes each, with downward velocities ranging between 1.0 and 2.1 km s −1 . However, ripples associated with sunquakes were not detected in all of the locations of these downward motions (Zharkova and Zharkov, 2007) but only in the three locations where the downward velocities of the shocks exceeded 1.5 kms −1 . This means that only a certain amount of momentum delivered to the photosphere within a very short timescale can initiate observable sunquakes.
Although the occurrence of seismic signatures in flares is often close to the locations of HXR and WL emission (Donea, 2011), they do not always coincide with them, as reported for the flare of 28 October 2003 (Zharkova and Zharkov, 2007). A similar picture was detected in GONG data for the flare of 14 December 2006 , where only two seismic sources were associated with HXR sources, while the other two did not show any or very little HXR emission. The first sunquake of Cycle 24, which occurred on 15 February 2011, revealed no significant HXR emission associated with one of the sunquake locations, while the other showed some HXR emission, which was much lower than that measured during the main flaring event.
From this point of view, even if backwarming does exist, which is often the case, its starting time after the impulsive phase and extended duration (for up 40 minutes) cannot deliver energy and momentum fast enough to the photosphere. For this reason, it cannot account for the timing and energy required to produce seismic ripples close to the start of HXR emission (Kosovichev and Zharkova, 1998;Zharkova andKosovichev, 2000, 2002) or within a minute or two after the HXR onset (Zharkova and Zharkov, 2007).
The points discussed above increase the requirements for delivering energy and momentum to the photosphere and lead to consideration of hydrodynamic shocks produced not only by electron (Kosovichev and Zharkova, 1998;Kosovichev, 2007), but also by proton beams (Zharkova and Zharkov, 2007) or magnetohydrodynamic shocks produced by the Lorentz force (Fisher et al., 2012). These cases emphasize the importance of considering mixed beams of electrons and protons, which most likely occur during the flare onset, as discussed in Section 2.1.1.

Hydrodynamic Response
To emulate the dynamics of a flaring atmosphere, we simulated the hydrodynamic response of the two-temperature ambient plasma to the injection of energetic particles (either electrons or protons) considering the ion viscosity (Somov, Spektor, and Syrovatskii, 1981). We solved the two energy (for electrons and protons), continuity and momentum conservation equations with the ion viscosity term (Somov, Spektor, and Syrovatskii, 1981;Zharkova and Zharkov, 2007), including the radiative cooling by optically thin coronal emission and optically thick hydrogen emission (Kobylinskii and Zharkova, 1996). The initial conditions at the time t = 0 were considered to be the photospheric ones: T (ξ 0 , 0) = 6700 K, n(ξ 0 , 0) = 10 13 cm −3 , v(ξ 0 , 0) = 0 cm s −1 , and the original atmosphere was considered to be the simple exponential one.
This approach is essentially different from the simplified one-temperature approach developed by Fisher (1989) that does not consider the ion viscosity. The hydrodynamic conditions in this approach were applied to a pre-heated semi-empirical atmosphere using quiet-Sun physical parameters (Vernazza, Avrett, and Loeser, 1981). The outcome of this hydrodynamic simulation was then used by Fisher et al. (2012) to explain some sunquake properties. The differences in the results of these two different models is discussed in the next sections. Zharkova and Zharkov (2007) evaluated the ambient heating caused by electron and proton beams or by a mixed beam. They also studied the implications for the ambient plasma that responds to this heating (variations of temperature, density, and macrovelocities); we use a similar method in this article. We briefly recall the main points.

Heating Functions by Energetic Particles
Similar to Somov, Spektor, and Syrovatskii (1981) and Zharkova and Zharkov (2007), for the heating function by beam electrons we used the formula derived from a continuity equation by Syrovatskii and Shmeleva (1972). This heating function is essentially different from the heating function used by Nagai and Emslie (1984), Allred et al. (2005) and Fisher (1989). Their heating function was derived by Emslie (1978Emslie ( , 1980 from the flux conservation equation, where the maximum energy deposition by electrons occurs much higher in the chromosphere compared to that derived by Syrovatskii and Shmeleva (1972). Furthermore, the flux conservation approach was shown by Mauas and Gómez (1997) to produce a strong singularity in the electron heating function (in fact, it becomes infinity) at the stopping depth of low-energy electrons. To avoid this singularity, Emslie (1978Emslie ( , 1980 selected the shape of the heating function with the maximum heating occurring above this stopping depth, which is not the best approximation, as shown by Mauas and Gómez (1997).
In this article, similarly to Zharkova and Zharkov (2007), the heating of the ambient plasma by electron beams occurs due to Coulomb collisions and Ohmic heating by precipitating particles (electrons, protons, and their mixture). The heating by proton and/or mixed beams is considered via Coulomb collisions and generation of kinetic Alfvén waves (KAWs) with their subsequent dissipation in Cherenkov resonance at the atmospheric depths where their velocities are higher than the local Alfvén ones .
As expected, the heating by electron or proton/mixed beams with power-law energy distributions is found (Zharkova and Zharkov, 2007) to be strongly dependent on the initial beam parameters: softer and weaker beams deposit their energy mainly in the corona and upper chromosphere, whilst harder and more powerful beams deposit more energy deeper in the lower chromosphere. The other heating mechanism for electron beams, Ohmic dissipation in a self-induced electric field, contributes significantly to the heating of the plasma at coronal levels, but it does not significantly affect the heating of the lower atmosphere (Zharkova and Gordovskyy, 2006;Siversky and Zharkova, 2009b).
Evidently, before returning to the source in the corona, some of lower energy electrons lose the bulk of their energy in Coulomb collisions contributing to heating of the upper layers in the flaring atmosphere and depositing momentum into lower atmospheric levels (Zharkova and Zharkov, 2007). At the same time, these electrons create their own strong electrostatic electric field. This electric field leads to a large number of lower energy electrons returning to the corona, or moving nearly across the magnetic field at some depth in the chromosphere (electrical stopping depth). This is caused by the significant reduction in their energy by the self-induced electric field and anisotropic scattering in the presence of this field (Siversky and Zharkova, 2009b). This process would definitely increase the energy deposited at higher atmospheric levels and reduce the amount of energy deposited at lower atmospheric levels via a hydrodynamic response.

Hydrodynamic Response of a Flaring Atmosphere
We are here mainly interested in hydrodynamic shocks as the sources of sunquakes, produced in the ambient plasma that is heated by energetic particles. We calculate hydrodynamic responses and resulting shocks produced by different mixtures of particles, electrons and protons, injected into the loop legs from the top (see Section 2.1.1).
The physical conditions in a flaring atmosphere are described by a plasma density n, electron T e and proton/ion T i temperatures, and a vertical velocity v as a function of a vertical coordinate z, or a column density ξ , and time t . The initial conditions at the time t = 0 are considered to be the photospheric ones: T (ξ 0 , 0) = 6700 K, n(ξ 0 , 0) = 10 13 cm −3 , v(ξ 0 , 0) = 0 cms −1 . Then we calculate these parameters for different times as a hydrodynamic response of the 1D solar atmosphere to the injection of electrons and/or protons followed by radiative cooling using the continuity, momentum, and energy equations for the ambient electrons and protons/ions following Zharkova and Zharkov (2007) (see also Somov, Spektor, and Syrovatskii, 1981;Fisher, Canfield, and McClymont, 1985a,b,c;Allred et al., 2005 andFisher, 1989).
Similar to Somov, Spektor, and Syrovatskii (1981), in our simulations we consider the momentum equation for a two-temperature plasma approach, in addition to the two energy (for electrons and for ions) equations and continuity equation for injected particles (Somov, Spektor, and Syrovatskii, 1981;Zharkova and Zharkov, 2007). This allows us to include the momentum deposited by these particles, which is shown to be very important for particle dynamics (Brown and Craig, 1984). Fisher (1989) also included the momentum equation in a very simplified form without considering the energy exchange between the ambient electrons and ions included in our analysis (Zharkova and Zharkov, 2007). For the radiative cooling term we consider the losses from optically thin emission of all elements with the coronal abundance, similar to Fisher (1989), and the emission in optically thick hydrogen lines and continuum calculated for the 5 plus continuum model of the hydrogen atom (Kobylinskii and Zharkova, 1996).
As a result of the hydrodynamic simulation for a flaring atmosphere, one can derive (within the first seconds after the beam injection onset) a sharp increase of the temperature and decrease of the density in the corona, formation of a transition region just above the chromosphere, and the formation of a low-temperature condensation in the lower chromosphere (Somov, Spektor, and Syrovatskii, 1981;Fisher, Canfield, and McClymont, 1985a,b,c;Zharkova and Zharkov, 2007). These processes are accompanied by the bulk plasma motion: the fast explosive upward motion of the photospheric plasma into the chromosphere and corona, called an 'explosive evaporation', and a supersonic motion of the low-temperature condensation towards the photosphere.
These two motions of the ambient plasma (upwards and downwards) are reported in all the hydrodynamic simulations, and the velocities of this bulk motion are called macrovelocities (Somov, Spektor, and Syrovatskii, 1981;Nagai and Emslie, 1984;Fisher, Canfield, and McClymont, 1985a,b,c). These macro-motions are widely investigated from the blueand red-shifted spectral measurements in UV and Hα emission, respectively (Zarro, Slater, and Freeland, 1988). For some events, or some beam parameters, the momenta of these two motions can be nearly equal . Some other events reveal only the upward Figure 1 Snapshots of the macrovelocity profiles for different times after beam injection onset (colour-marked description on the right) as a function of the column depth (X-axis) with the shocks formed at different times after the injection of (top plot) a pure electron beam with the parameters derived from HXR emission (initial energy flux 5 × 10 11 erg cm −2 s −1 , spectral index γ = 4.9), and (bottom plot) a mixture of 80 % of an electron beam and 20 % of a proton beam (particle number) (an initial energy flux 3 × 10 12 erg cm −2 s −1 , γ = 3.0). The macrovelocity profiles are given on the Y -axis in logarithmic scale in cm s −1 .
motions with blue-shifts indicating explosive evaporation without noticeable downward motions (to be seen as red-shifts). For still some other events, stronger red-shifts and weaker blue-shifts are observed, indicating a more modest chromospheric evaporation.
The hydrodynamic heating of flaring atmospheres caused by injected electron and proton beams has been discussed in great detail in our previous articles (Zharkova and Zharkov, 2007;Zharkova, 2008). Zharkova and Zharkov (2007) simulated and presented the full set of plots for electron temperature, ambient density, and macrovelocity in response to the injection of electron and proton beams. They showed that the low-temperature condensation formed in response to either of these beams has a temperature below 10000 K and the density exceeds 10 14 -10 15 cm 3 . This condensation moves downward to the photosphere with a speed of 10 7 cm s −1 in the chromosphere and slows down to a few units of 10 5 cm s −1 before its full termination at the photosphere. Furthermore, Zharkova and Zharkov (2007) showed in their Figure 9 that the temporal profile of the shock termination produced by protons was very close to the temporal profiles directly measured from Dopplergrams for the flare of 28 October 2003.
In this article we decided to explore the heating and hydrodynamic response in flaring atmospheres extended towards the photosphere and beneath in an attempt to detect the times when or if the lower-temperature shocks can reach the photosphere and how deep they can push their photospheric level into the solar interior compared to this level in the quiet atmosphere.
In Figure 1 the upward and downward macrovelocities of the bulk of the ambient plasma are plotted versus a column depth similar to the plot presented by Zharkova and Zharkov (2007), reflecting the macromotion of the ambient plasma in response to the injection of a powerful electron beam with an initial energy flux of F 0 = 5×10 11 erg cm −2 s −1 and spectral index γ = 4.9 (top plot) and mixed beam (80 % of electrons and 20 % of protons) with the initial energy flux of F 0 = 3 × 10 12 erg cm −2 s −1 and spectral index γ = 3 (bottom plot). Figure 1 shows a very significant difference in the generation of macrovelocities by different agents, which define different velocities for chromospheric evaporation into the Figure 2 Snapshots of the macrovelocity profiles (Y -axis, cm s −1 ) for different times after the beam injection onset (colour-marked description in the plot insets) as a function of the column depth ξ (the top plot) and the linear depth Z (the bottom plot) with respect to the quiet-Sun photosphere (zero level). The shocks are formed in response to the injection of a beam with 70 % electrons and 30 % protons with an initial energy flux of F = 5 × 10 12 erg cm −2 s −1 , γ = 3.0. We note that the pre-flare photospheric level is at a column depth of 10 23 cm −2 (or the linear depth is zero), while during a flare this level moves below the level of the quiet photosphere by 5000 km.
corona as well as different depths for the formation of low-temperature shocks and their termination velocities. The macrovelocity variations caused by a pure electron beam with the initial flux of 2 × 10 12 erg cm −2 s −1 , and spectral index ≈ 5 (top plot) reveal that very high velocities of chromospheric evaporation (up to 10 8 cm/s) occur from the upper chromosphere (column depth up to 5 × 10 19 cm −2 ). In addition, a lower temperature shock forms immediately below this depth, reaching velocities of up to 10 7 cm s −1 in the upper chromosphere while terminating at the column depth just before 10 21 cm −2 with a velocity of up to (1 -2) × 10 5 cm s −1 . This reflects the fact that pure electron beams have a maximum heating function in the upper chromosphere (Syrovatskii and Shmeleva, 1972;Somov, Spektor, and Syrovatskii, 1981), therefore the shock forms directly below the transition region.
For the beams with 20 % protons and 80 % electrons the initial energy flux is noticeably increased (see Figure 1, bottom plot), leading to a shift of the energy deposition maximum from the lower corona to the upper chromosphere. This reduces the macrovelocities of chromospheric evaporation into the corona from 10 8 cm s −1 in the first ten seconds to a factor of 10 5 cm s −1 after 130 seconds. At the same time, the lower temperature shock is now split into two shocks. One moves from the lower corona towards the chromospheric column depth of 10 19 cm −2 and terminates at depths of 10 20 cm −2 , which is a chromospheric level. The second shock appears at a column depth of 4 × 10 20 cm −2 and terminates at a depth of 1 × 10 21 cm −2 . The maximum velocity approaches 10 6 cm s −1 in the first shock and a few units of 10 5 cm s −1 in the second.
The macrovelocities of the bulk of the ambient plasma induced by a hydrodynamic response to the injection of mixed beams are plotted in Figure 2 for the electron-dominated beam with F = 5 × 10 12 erg cm −2 s −1 and in Figure 3 for the proton-dominated beam with the initial energy flux F = 5 × 10 13 erg cm −2 s −1 .
A further increase of the initial energy flux and momentum ( Figure 2) for a slightly more intense mixed beam with an initial flux of 5 × 10 12 erg cm −2 s −1 (containing 70 % electrons We note that the pre-flare photospheric level is located at a column depth of 10 23 cm −2 (or the linear depth is zero), while during a flare this level moves below the level of the quiet photosphere by 500 km. and 30 % protons) leads to formation of the double shocks, similar to those described in Figure 1 (bottom plot). For a larger contribution of protons, the upper shock spreads out much more from a column depth of 8 × 10 19 cm −2 to 10 22 cm −2 , with the second shock starting at a depth of 5 × 10 22 cm −2 and terminating at the depth 10 23 cm −2 with a velocity of a few units 10 5 cm s −1 .
The upper shocks can produce Moreton waves in the chromosphere, which are occasionally observed in association with sunquakes (Kosovichev and Zharkova, 1998;Warmuth and Mann, 2011), while the bottom shocks can be considered to be strong candidates for producing sunquake ripples (Zharkova and Zharkov, 2007).
The simulation plotted in Figure 3 in function of column depth (top plot) considering the beams of protons (70 % of a total abundance) and electrons (30 % of a total abundance) shows a reduction of the velocities of chromospheric evaporation (up to 10 6 cm s −1 ) while still keeping the wide first low-temperature shock spreading from the whole chromosphere with the second low-temperature shock starting at a depth of 5 × 10 22 cm −2 and terminating at the depth 10 23 cm −2 with the velocity of a few units of 10 5 cm s −1 .
To show the linear depth of the formation of this second shock, we plot in Figures 2 and 3 the macrovelocities of the ambient plasma versus a column depth (top plots) and a linear depth (bottom plots), with the zero level corresponding to the undisturbed photospheric level (defined at the start of simulations by the boundary condition for the hydrodynamic equations).
In the case of the mixed beam dominated by electrons (see Figure 2) the second shock can travel from the photosphere and beneath to a distance of about 5000 km (bottom plot). Since the two shocks are initiated in the upper and middle chromosphere of a flaring atmosphere, the latter is swept below the quiet-Sun photospheric level. This leads to large linear depths (up to 5000 km) in which these shocks travel below the quiet-Sun photosphere before they finally terminate. In fact, the shocks are both created and terminated smoothly at increasingly deeper levels that approach 5000 km.
With a stronger contribution of protons in hydrodynamic heating, the second shock generated by this mixed beam is terminated in the solar interior at depths of about 300 km below the quiet-Sun photosphere (see Figure 3, bottom plot). This happens because the proton-dominated beam creates the shock in a much denser plasma at the bottom of the chromosphere (Zharkova and Zharkov, 2007). As a result, this shock sweeps much denser plasma below the quiet-Sun photospheric level, and it travels a much shorter linear distance before terminating. In this particular simulation, the shock moves only 300 km below the quiet photosphere level before it terminates, as shown in Figure 3 (bottom plot).
The results of the hydrodynamic simulations in general confirm within the limits of the simplified model used (Fisher, 1989) the estimations by Fisher et al. (2012) of the momentum and energy delivered by hydrodynamic shocks. However, our hydrodynamic results have the great advantage that we know the initial condition for the quiet photosphere (before flare onset) and consider a two-temperature plasma with the energy exchange between the ambient electrons and protons (see Section 2.2). This allows us to reliably calculate a linear depth within a flaring atmosphere compared to the quiet-Sun atmosphere level.
In contrast, the model by Fisher (1989) and Fisher et al. (2012) ascribed the initial condition to the semi-empirical quiet-Sun model derived from fitting the observed solar emission in multiple lines (Vernazza, Avrett, and Loeser, 1981), which is dependent on settings of the non-LTE model and is known to underestimate the height of the quiet photosphere level (Somov, Spektor, and Syrovatskii, 1981). This immediately moves the zero level of a flaring atmosphere to a level below the photosphere in the models used by Somov, Spektor, and Syrovatskii (1981) and Zharkova and Zharkov (2007).
In addition, the mild heating functions by electron beams with initial fluxes below 10 11 erg cm −2 s −1 used in the model by Fisher (1989) and Fisher et al. (2012) do not produce the hydrodynamic responses with the shocks formed in the lower chromosphere and the photosphere like those reported in this article. These numerical limitations impose some restriction on the outcomes in the hydrodynamic models used by Fisher (1989) and Fisher et al. (2012), allowing only detection of weak shocks in the chromosphere above the photospheric level, which are unable to produce seismic signatures in flares similar to those reported in this study.
The occurrence of either single or double hydrodynamic shocks and their termination either in the solar atmosphere or the interior well below the undisturbed photospheric level can explain the observational features related to seismic responses, or sunquakes, that occur in the solar interior and their occasional links with Moreton waves, which occur in the chromosphere, and coronal waves measured by the SOHO Extreme ultraviolet Imaging Telescope (EIT) (Warmuth and Mann, 2011). Only mixed beams (with a substantial proportion of protons) with high energies and momenta can form two shocks: one in the upper atmospheric levels, leading to Moreton waves, and the other one at the photospheric level, which terminates deeply in the solar interior and leads to sunquakes (see Section 3 below). Table 1 Estimates of momenta P hd and kinetic energies E k of hydrodynamic shocks in the three seismic sources terminated with the areas A with downward (vertical) velocities V vert , with duration T of the sunquake ripples, and the sunquake average horizontal velocity, V horiz (for details see the article by Zharkova and Zharkov, 2007

The Momentum and Energy Delivered by a Hydrodynamic Shock
The momentum delivered by a hydrodynamic shock can be evaluated using the simple formula (Zharkova and Zharkov, 2007) where the summation is made over the time from 0 to τ , τ is the duration of the impact causing the seismic waves, m is the mass of the plasma delivering the momentum related to the flaring area A where the momentum is deposited, V vert is a starting velocity of the plasma at the impact time, and t is the duration of the impact. For the known plasma mass density ρ = m H · n, where n is the particle density per volume defined from hydrodynamic solutions, this equation can be rewritten as follows: where ρ is the average density of the plasma delivering the momentum, A is the flaring area where the momentum is deposited, v vert is the averaged vertical velocity of the plasma propagation at the impact time, and τ is the duration time of the impact. Then, the kinetic energy delivered by a hydrodynamic shock can be written as Of course, this energy is deposited over a depth (column depth ξ or linear depth z) of a flaring atmosphere, into which this shock moves during time t . Only a small part of the shock is deposited into acoustic waves; this fraction is left at the moment when the hydrodynamic shock terminates into dense plasma, which causes the acoustic waves in the solar interior.
To estimate this fraction of the energy, we need to estimate the plasma density, the flare area where the sunquake occurs, the duration of the shock propagation (or the existence of the sunquake), and a measured Doppler velocity of the downward motion (hydrodynamic shock) in the location of the sunquake. Then, the momentum and energy deposited by a shock in a flare with the given physical conditions can be estimated following the techniques described step-by-step by Zharkova and Zharkov (2007). An example of estimating the momentum for the sunquakes detected for the flare of 28 October 2003 (Zharkova and Zharkov, 2007) is given in Table 1. The energy deposited by the hydrodynamic shocks generated in a flaring atmosphere presented in Table 1 is much higher than that derived from the acoustic waves produced by flares, which normally ranges in 10 27 -10 28 erg (Donea, 2011;Zharkov et al., 2013). The first reason is that we only measure in the sunquakes the energy of the first bounce of the acoustic waves that are induced by these shocks, and they are derived as a ridge in the timedistance diagrams while many more bounces may occur beyond the datacube considered. The second reason is that some of the deposited energy can be spent on a drag force inside the solar interior, and, third, the properties of the photosphere (density, temperature, etc.) where the acoustic waves reflected from are likewise unknown. All these properties need yet to be investigated in more detail from observations of sunquakes.

Hydrodynamic Response of the Solar Interior and Formation of Acoustic Waves
Sunquakes are relatively rare events that occur in active regions during flares when a flare deposits sufficient momentum and energy into the lower atmosphere and, under certain conditions, can generate acoustic waves in the solar interior. These waves are seen in photospheric Dopplergrams as circular shaped surface ripples accelerating away from the source region where the impact occurs. However, random background oscillations often make such ripples difficult to distinguish. Therefore, helioseismic methods such as a time-distance diagram analysis (Kosovichev and Zharkova, 1998) and helioseismic holography (see for example Lindsey and Braun, 1997;Donea, Braun, and Lindsey, 1999, and references therein) are usually applied to detect sunquakes in suitably processed series of 2D line-of-sight velocity observations (Dopplergrams), called datacubes (varying from 180 Mm×180 Mm to 250 Mm ×250 Mm), obtained with cadences of 45 seconds to 1 minute for the duration of two or more hours.
The time-distance diagram technique is a direct observational method using the position of the quake source in the 3D datacube (two horizontal dimensions and time) to remap it to a 2D time-distance image while deriving the surface velocity variations as a function of time and radial distance from the assumed source. Sunquake ripples in this wave are seen as a ridge on the time-distance diagram (2D image).
On the other hand, the helioseismic egression measurements are based on the theoretical modelling of acoustic signals propagating from a point source and emanating outwards from the flare location. This allows producing 2D egression power maps of the acoustic sources and sinks, where the quake signatures are represented by compact bright kernels of the enhanced emission surrounded by the acoustically absorbing interior (Lindsey and Braun, 1997). This procedure was used for many flares (e.g. see Donea et al., 2006;Matthews, Zharkov, and Zharkova, 2011;.
Both helioseismic methods show the very localized nature of the origin point of seismic sources (where the ripples start from), validating the point source assumption used in the theoretical interpretation (Kosovichev and Zharkova, 1995;Zharkov, 2013) and in the observational analysis of Dopplergrams for detecting ripples (see for example Kosovichev and Zharkova, 1998;Matthews, Zharkov, and Zharkova, 2011;. We here consider that these point-like sources are the hydrodynamic shocks of the ambient plasma formed in flaring atmospheres by the mixed beams, which release their energy and momentum into the production of acoustic wavepackets in the solar interior. In fact, if the sources are the hydrodynamic shocks caused by mixed beams, as shown in Figures 2  and 3, then it becomes evident that these shocks move below the level of the quiet-Sun photosphere for a distance of a few hundred to a few thousand kilometers. With every second Figure 4 The individual acoustic rays (generated at the depth of 500 km and travelling to the bottom of the plot) by a moving supersonic source depositing a momentum below the photosphere (the origin) under some angle to the local vertical for the times shown above the plots. The photosphere is denoted by the top line, Y -axis shows a depth in Mm beneath the photosphere and the X-axis denotes a distance in Mm from the point of momentum deposition. The red arrow shows the direction of the momentum deposition. The rays are reflected from the photosphere level after 20 minutes (the bottom right plot), which can be observed as a ripple on the surface, or a sunquake, generated at about the central point of the momentum deposition. these shocks move downwards with a supersonic speed (up to 12 km s −1 for a mixed beam with a larger portion of protons or up to 6 -8 km s −1 for an electron-dominated mixed beam, with the adiabatic soundspeed, c being around 7 km s −1 near the surface and growing to around 20 km s −1 at 5000 km depth).
These hydrodynamic shocks from the atmosphere above can be used as the initial condition for another hydrodynamic model developed for acoustic wave propagation in the solar interior. As the initial hydrodynamic shock terminates within relatively shallow depths and with strongly supersonic velocities, the generated waves are formed at the point of deposition as a closed cone around the velocity vector in the solar interior, which, in accordance to Fermat's principle, propagates deeper into the interior, refracting as a result of the increasing temperature, and reflecting back to the photosphere (see Figure 4, bottom plots).
The model of the interior can be assumed to be either a polytrope (Zharkov, 2013) or the standard interior model (Christensen-Dalsgaard et al., 1996) used in this study. By solving the hydrodynamic equations for acoustic wave propagation in the interior (Zharkov, 2013), we can evaluate the parameters of the generated acoustic wave packets and the condition of their detection from Doppler observations.
We note that the wave packet comprises a large number of individual rays with different wavelengths or frequencies. Then, as shown by Zharkov (2013), the vertical shock perturbation moving with a supersonic velocity can generate the set of multiple acoustic waves, from which only the waves whose phase speed exceeds a certain threshold (Zharkov, 2013, see Equation (5.8)) can produce the observable acoustic waves. An individual ray characterized by a constant (along the ray) frequency, ω, and horizontal wavenumber, k h , initialized at a given depth, will have two, upper and lower, turning points (see Figure 1 in Zharkov, 2013). The first upper turning point along the ray defines its first surface appearance (as a first ripple), the lower turning point indicates where the wave changes the direction of its motion in the interior by being reflected back to the surface. Then the propagating ripples correspond to a sequence of the source-generated acoustic rays from the packet, successively reaching their upper turning points.
The source of the deposited impulse, depending on its properties, generates a family of the rays that provides the solution to the ray equations in the phase space and defines the generated wave front. As the source is located in the interior, the first ray (out of all generated by the source) to reach its upper turning point defines the minimal distance where the ripple is formed. This distance will depend on the source depositing the momentum, the depth where this momentum is deposited, and the interior model as described below.
For a near-surface source, i.e. for the ray initiated near its upper turning point, the first surface appearance, or the minimal distance, can be approximated by the ray skip distance, , the distance between its surface bounces. This distance depends on the ray horizontal phase speed, ω/k h (see for example Christensen-Dalsgaard et al., 1996;Zharkov, 2013). For the polytrope model of the solar interior used in Zharkov (2013), the minimal skip distance, , or the distance from the point of the initial impulse deposition to the first ripple occurrence, is derived from Equations (A1) in the Appendix A1 (Zharkov, 2013) as follows: where g is the gravitational constant, g = 2.67 × 10 −4 Mm s −2 , m is the polytrope index, and V ph = w/k h is the horizontal speed of wave propagation. For the supersonic (non-oscillatory) source (see Section 5.3 in Zharkov, 2013, and the corresponding plots), the waves of the packet are generated in a cone/sphere around the velocity vector, with the ray frequencies increasing with the angle measured from the velocity vector (see Equation (5.5) in Zharkov (2013)). In this case, the observations of highfrequency waves will be also limited by the Nyquist frequency and cadence of the series (8.33 mHz for MDI, 11.11 mHz for SDO/HMI). Thus, if the observations are made at a certain Nyquist frequency, ω N , the high-frequency waves may not be observed, therefore further restrictions for 'observable ripples' are considered, namely ω ≤ ω N holds.
This leads to the threshold condition for the minimal phase speed, v min ph , defining the condition for registering the first ripples on the surface by the following relation (see Equation (5.8) in Zharkov (2013)): where ω αc is the acoustic cutoff frequency and c is the speed of sound. The propagation of these waves is tracked by the rays, shown in Figure 4. The abscissa defines the horizontal distance in Mm about the locations of supersonic disturbances (shocks) and the ordinate shows the depth, z s , below the photosphere. The moving source shown in Figure 4 introduces an anisotropy in the produced acoustic wave packet. The waves are generated in a cone around the velocity vector, with frequencies increasing with the increase of the angle between the initial shock velocity vector and the ray take-off direction.
The phase speed of the ripple appearance increases with the distance from the source of the original impact, as seen in all the time-distance diagrams of sunquakes (see, for example Kosovichev and Zharkova, 1998;. Since the propagation of surface ripples from a near-surface source can be determined by the phase speed (Zharkov, 2013), the minimal distance, , where these ripples are observed can be estimated from Equation (4) after substituting the minimal phase speed given in Equation (5). For the horizontal velocities from 15 to 45 km s −1 reported by Zharkova (1995, 1998), the first bounce normally comes to the surface at a skip distance about = 40 Mm that is seen in about 20 -25 minutes after the initial impulse deposition in the source.
The other ripples, which are observed farther from the source, are the first bounces of the other rays in the wave packet. The linear phase speed of rays increases with a decrease of the ray frequency (or an increase of the ray velocity) (see Equation (5)). Therefore the ripple produced by the next ray is seen much farther from the first one. This process will continue until the location of the next ripple occurrence exceeds the size of the datacube. This defines the conditions under which the ripples can be detected depending on the initial deposition velocity.
For example, for the initial supersonic shock produced by an electron-rich beam and deposited with the velocity of v(z s ) = 10 km s −1 in the interior with the sound speed of c(z s ) = 7.5 kms −1 , ω αc (z s )/2π = 5.4 mHz and for ω N = 2π · 8.3 mHz (corresponding to the Dopplergram cadence), the condition that the cyclic frequency along the ray is no greater than Nyquist frequency ω N , gives us the minimum phase speed estimate, v ph , over of 60 km s −1 . After substituting this speed into Equation (4), the minimal skip distance where we can observe the first ripple is found to exceed 80 Mm.
It is evident that the ripples moving with a speed of 60 to 80 km s −1 (or an average speed of 70 km s −1 ) reach a distance beyond 180 Mm (70 · 45 · 60 = 189 Mm) within less than 45 minutes. Hence, this ripple propagation will be rather difficult to detect from the timedistance diagrams because it will be seen far away from the source. This would suggest that for a reliable observation of ripples on the solar surface associated with acoustic waves, the limit of their horizontal phase velocity should not exceed about 50 km s −1 .
These analytical estimations are in line with the numerical simulations showing the generation of acoustic waves by convective vortices in the quiet-Sun interior (Kitiashvili et al., 2011;Moll, Cameron, and Schüssler, 2011), reporting also supersonic movements at the vortices, similar to those derived from the observations of sunquakes (e.g. see Kosovichev and Zharkova, 1998;Kosovichev, 2007;Matthews, Zharkov, and Zharkova, 2011;Zharkov et al., 2013). These similarities strongly support the analytical theory of acoustic waves induced by the hydrodynamic processes in flares (Zharkov, 2013).
In addition, the geometry of wave propagation will largely depend on the angle under which the initial shock enters the photosphere, and, because the rays are generated within a cone, their frequencies can be growing away from the initial velocity vector. In other words, if at the sunquake locations a supersonic movement is observed with a strong anisotropy in the amplitude in one direction, this indicates that the initial atmospheric shock was deposited under some angle to the surface in this direction. Therefore, it is logical to suggest that this mechanism plays an important if not the key role in the sunquake wave front generation, which is linked to the parameters of hydrodynamic shocks induced in flaring atmospheres by energetic particles.
In the application to the models considered in this study, the atmospheric hydrodynamic shocks generated by a proton-dominated beam (Figure 3, bottom plot) are found deposited consistently at a shallow depth range between 200 and 300 km under the photospheric level. The high supersonic speed of this shock v > c (v = 12 km s −1 and c = 6 km s −1 ), the shock propagation creates a sharp set of waves shown in Figure 4, derived from the hydrodynamic theory about a movement of the supersonic source in the solar interior (see Section 5.3 in Zharkov, 2013).
The shocks created by electron beams (see Figure 2, bottom plot) are at each instant deposited over a wide range of depths from 0 to 5000 km below the quiet-Sun photosphere. The shocks generated by electron-dominated beams have a speed (8 km s −1 ) that only marginally exceeds the sound speed at the top of the interior (6 km s −1 ).
Similar to the previous case with a proton-dominated beam, the acoustic waves generated by electron-dominated shocks are located within a cone of the velocity vector with the angular cone width determined by the ratio v/c. This means that the shocks can still generate individual rays at the upper depths of the solar interior, shown as separate lines in Figure 4. However, because these shocks are terminated smoothly at the deeper interior layers, this would result in the larger upper skip distances (Equation (4)) for their ripples to be observed. This means that these shocks can be reflected from the surface at much greater distances than the standard datacube size of 180 × 180 Mm 2 used in helioseismology of sunquakes.
In addition, in deeper layers of the interior the speed of the shocks generated by electronrich beams drops very quickly below the local sound speed, which would eliminate their ability to generate acoustic waves shown in Figure 4. Therefore, the summary set of acoustic waves (rays) produced by electron-dominated shock would be rather scattered in the interior. Thus, in this particular occasion the shock cannot not generate significant visible effects on the photosphere seen as ripples within the given datacube sizes.
For the cases presented here, we find that the hydrodynamic theory in the solar interior with acoustic cutoff frequency (Zharkov, 2013) can logically explain the generation of the observable acoustic waves, or sunquakes, in solar flares observed in the datacubes of 180 × 180 Mm 2 . These detectable acoustic waves are shown to be linked closer to the shocks induced by the mixed proton-dominated beams. In contrast, the acoustic waves induced by the electron-dominated shocks can have a much wider spectrum and can therefore be reflected from the surface at much larger distances than the standard datacubes used in the helioseismology of sunquakes.
Of course, when observers register sunquakes within a given datacube, they can only derive the energy of the first bounce of these acoustic packets, which excludes the majority of waves, which travel much further in the solar interior. This limits the measurements of acoustic wave energy to approximately 1 -10 % of their total energy. In addition, to estimate the fraction of the energy transferred from the hydrodynamic shocks deposited from a flaring atmosphere, presented in Section 2.2.3, to the acoustic waves in the solar interior, discussed in this section, the drag force of the waves during their propagation in this interior needs to be taken into account; this force leads to its heating. The properties of the photosphere where the waves are reflected from, while producing enhancement of the plasma above the surface seen as ripples, also need to be further investigated from the observations of sunquakes. These are the important tasks for future studies before any modelling for a oneto-one reproduction of acoustic waves generated during flares can be carried out.

Conclusions
We have combined for the first time the two hydrodynamic models associated with solar flares to derive the conditions for a favourable observation of acoustic waves, or sunquakes. The first hydrodynamic model was developed for a two-temperature flaring atmosphere heated by mixed particle beams (protons plus electrons) (Somov, Spektor, and Syrovatskii, 1981;Zharkova and Zharkov, 2007). The heating produces an increase of large macromotions of the ambient plasma upwards (explosive evaporation) and downwards (lowtemperature hydrodynamic shocks), moving to the photosphere and below it.
The velocities and directivity of these shocks were considered to be the initial condition for the second hydrodynamic model (Zharkov, 2013) used for the investigation of acoustic wave generation and propagation in the solar interior below a flare. This model defines the acoustic cutoff frequency and the minimal phase velocity above which the acoustic waves can be measured from the surface observations. This approach is different from the one with a single kinetic temperature and simplified momentum equation used by Fisher et al. (2012). We also considered the initial condition used for the undisturbed photosphere with a density of 10 13 cm −3 and a temperature of 6700 K that sets the quiet-Sun photospheric level. In the other hydrodynamic models the authors considered a semi-empirical model of a flaring atmosphere as the initial condition.
The heating functions for beam electrons and protons were derived from the continuity equations (Syrovatskii and Shmeleva, 1972;Somov, Spektor, and Syrovatskii, 1981;Gordovskyy et al., 2005), allowing us to increase their initial energy fluxes of precipitating beams to very high magnitudes (up to 5 × 10 13 erg cm −2 s −1 , for a proton-rich beam). This is in contrast to the heating function derived from the flux conservation equation used by Fisher et al. (2012), which is known to have well-known limitations on the energy deposited by lower energy electrons at the chromosphere (Mauas and Gómez, 1997). This, in turn, imposes essential limitations on a maximum energy flux of injected beams that can be used in hydrodynamic models.
These hydrodynamic simulations allowed us to calculate the macrovelocity profiles to vary not only with column depth, but also with a linear depth measured from the photospheric level of the quiet Sun. We compared the linear depths of deposition of the hydrodynamic shocks generated by hydrodynamic responses of the ambient plasma to the injection of pure electron beams, electron-rich (with some fraction of protons) beams, and proton-rich (with some fraction of electrons) beams.
We found that the shocks induced by pure electrons are terminated either above the photosphere or just at its surface. The shocks induced by proton-rich beams were found to terminate at the narrow band of 200 -300 km below the quiet solar photosphere, while the shocks produced by electron-rich beams were terminated smoothly over depths from 0 down to 5000 km below the quiet-Sun photosphere.
The shocks induced by mixed beams gain a downward macrovelocity at the surface of about 8 km s −1 (electron-rich beam) and 12 km s −1 (proton-rich beam), which exceeds the sound speed of 6 km s −1 at the depths just below the photosphere. This increases to 20 km s −1 in the deep interior. The shocks were found to terminate at some deeper layers with a speed of a few km s −1 .
By using these hydrodynamic shocks from the above atmosphere as the initial condition for another hydrodynamic model developed for acoustic wave propagation in the solar interior (Zharkov, 2013) following a standard polytrope model, we evaluated the parameters of the generated acoustic wave packets and the conditions of detecting these waves from Doppler observations. A vertical shock perturbation moving with a supersonic velocity was shown to generate a cone-like multiple acoustic wave-packet in the direction of the shock motion. From this wave packet only the waves with the horizontal phase speed exceeding a certain threshold can produce observable acoustic waves, or sunquakes (Zharkov, 2013).
According to this second hydrodynamic model of acoustic wave propagation in the solar interior, the atmospheric shock caused by a proton-rich mixed beam propagates with a high supersonic speed and deposits its momentum at linear depths of 200 -300 km below the quiet-Sun photospheric level. This means that this shock generates a set of acoustic waves with very similar parameters, which propagate in the solar interior.
In contrast, the shock produced by the electron-rich beam considered in this article is deposited nearly continuously over the depth from below the surface down to 5000 km. In this case, acoustic waves are reflected from the surface at large distances of several hundred Mm, which are not included in the datacubes used in local helioseismology of flares. The continuous deposition over a depth of 5000 km also means that the acoustic waves produced by this shock can smoothly travel in the interior and become smoothly reflected from the photosphere at large distances from the source without producing any noticeable signatures (ripples) on the surface within the datacube used.
We noted in this study that the energy deposited by hydrodynamic shocks is much higher than that derived from the observations of sunquakes in flares. This reflects the facts that so far, only the energy of the first bounce of acoustic waves induced by these shocks is measured in sunquakes. This limits the measurements of the acoustic wave energy to approximately 1 -10 % of their total energy. In addition, some energy is spent on the drag force inside the solar interior, whose properties need yet to be investigated. Similarly, the properties of the photosphere need to be defined. Here the waves are reflected, which enhances the plasma above the surface and is seen as ripples.
Therefore, further theoretical investigation is required to determine the conditions for the generation, propagation, and directivity of acoustic waves in the solar interior that are induced by supersonic shocks of different nature (hydrodynamic or magneto-hydrodynamic), which are generated in flaring atmospheres and deposited below the solar surface of the quiet Sun with different velocities and at different angles.

Declaration of Conflict of Interests
The authors declare that they have no conflicts of interest.