Smoothing of MHD Shocks in Mode Conversion

TL;DR

This paper investigates how shock waves behave when passing through the Alfvén-acoustic layer in a stratified magneto-atmosphere, revealing that only the slow shock is smoothed while the fast shock persists, contrary to some predictions.

Contribution

It provides new insights into MHD shock smoothing and mode conversion, clarifying the behavior of fast and slow shocks in stratified magnetic atmospheres.

Findings

Slow shock is smoothed during mode conversion.

Fast shock persists without smoothing.

Partial mode conversion occurs as predicted by linear theory.

Abstract

Shock waves are simulated passing through the Alfv\'en-acoustic equipartition layer in a stratified isothermal magneto-atmosphere. The recent ray-theoretic calculations of N\'u\~nez (2019) predicted smoothing of the shock through this layer, causing both the fast and slow components to emerge as continuous waves. However, it is found that the partial mode conversion expected from linear theory for oblique incidence of the shock on the magnetic field is accompanied by a smoothing of the slow-shock only, whilst the fast-shock persists. Explanations are presented based on MHD mode conversion and shock theory.

Figures

Equations

$$\tilde{F}(f)=-f_0\frac{\cos (\pi f/2f_0)}{(f^2-f_0^2)\pi }$$

$$h_{\parallel }=\left|v_{\parallel ,\text{max}}/\left(\frac{d{v}_{\parallel }}{dz}\right)\right|,$$

$$\mathcal{T}\approx \exp {\left[-\pi h_s|\mathbf{k}|{\sin }^2\alpha \right]}_{a=c},$$

$$\begin{array}{rl}A(v_1,{\theta }_1,\gamma ,c_1,a_1)=& (v_1^2-X{\cos }^2{\theta }_1{a}_1^2{)}^2\left\{\left[(\gamma +1)-(\gamma -1)X\right]v_1^2-2X{c}_1^2\right\}\\ & -X{\sin }^2{\theta }_1{a}_1^2\left\{\left[\gamma +(2-\gamma )X\right]-\left[(\gamma +1)-(\gamma -1)X\right]X{\cos }^2{\theta }_1{a}_1^2\right\}=0.\end{array}$$

Full text

Sun: oscillations – Sun: chromosphere

1 Introduction

Magnetohydrodynamic (MHD) waves may change nature – between fast, slow and Alfvén – several times in moving through stellar atmospheres. Both fast-slow (Schunker & Cally, 2006) and fast-Alfvén (Cally & Goossens, 2008; Cally & Hansen, 2011) conversions have been explored in detail for linear waves (reviewed by Cally et al., 2016). Linear fast-slow conversion occurs where the Alfvén speed $a$ matches the sound speed $c$, and fast-to-Alfvén conversion is affected near where the fast wave reflects in a stratified atmosphere (provided the wave vector is not in the same vertical plane as the magnetic field lines). However, mode conversion by these processes is less well understood for nonlinear waves, in particular shocks.

Recently, Núñez (2019) presented an analysis based on the generalized ray theory of Cally (2006) and Schunker & Cally (2006) (see also Tracy et al. 2003) to argue that magnetoacoustic shocks passing through a layer where $a$ and $c$ coincide (the Alfvén-acoustic equipartition) not only split into fast and slow components as suggested by linear theory, but also that both the resulting fast and slow shocks are smoothed in the process.

Núñez’s analysis proceeds by examining the convergence of the Fourier components of the incoming wave. A shock discontinuity is associated with an $\mathcal{O}(1/n)$ asymptotic behaviour in wavenumber $n$ space. Faster convergence results in continuous (smooth) solutions. However, the argument is predicated on an assumption which is merely stated and not justified: ‘Assuming that the passage through the conversion zone acts linearly in the waves, which is reasonable since they spend very little time there, we obtain the surprising result that the shocks which may be present in the incident wave are smoothed out in the outgoing ones, except for the single case of a pure transfer from a fast to a slow wave or vice versa.’ (Núñez, 2019, Conclusions).

The purpose of this letter is to present the results of 1.5-dimensional (1.5D) simulations testing this assumption and the resulting finding that shocks are smoothed. We conclude that the assumption is in fact not valid – the mode conversion region is advected with the slow ‘shock’ for several seconds thereby extending the time over which conversion is operating – and that the smoothing is only partially realized; the slow shock is smoothed but the fast shock is not. Elementary explanations for this are presented in the discussion. These suggest that the smoothing of the slow shock but persistence of the fast shock may happen even without advection of the $a=c$ equipartition layer.

2 Model

We adopt the simplest suitable model to explore the phenomenon: an isothermal gravitationally stratified 1.5D atmosphere with uniform magnetic field ${\mathbf{B}}=B_{0}(\sin\theta,0,\cos\theta)$ inclined at angle $\theta$ to the vertical. The magnetic field is chosen such that the Alfvén speed $a$ equals the sound speed $c$ at a sufficient altitude that a fast (i.e. acoustic) wave injected at the bottom of the computational domain has enough space to shock before reaching the equipartition level $a=c$.

The computational box stretches $2.0$ Mm both above and below the $z=0$ height. The equilibrium density scale height is $h=175.731$ km, the uniform sound speed is $c=8.958$ $\text{km}\,\text{s}^{-1}$, the adiabatic index is $\gamma=5/3$, and the magnetic field strength is $B_{0}=0.9$ kG oriented at angle $\theta$ to the vertical. The $a=c$ equipartition level is situated at height $z_{0}=125$ km. The acoustic cutoff frequency is $\omega_{c}=c/2h$, or 4.0 mHz. There is nothing particularly special about these values, though they are broadly characteristic of the active solar chromosphere. Changing $B$ simply moves the equipartition level up or down. Similar results may be expected for other stellar atmospheric models.

The numerical code used to solve the non-linear MHD equations is Lare2d (Arber et al., 2001). The grid is comprised of 8192 cells in the vertical direction, giving a resolution of $0.488$ km.

A shock viscosity term is employed to avoid the Gibbs phenomenon whilst not artificially smoothing the shock front. It can be seen that the shock fronts remain sharp where predicted throughout the simulations. An artificial cooling term is applied to avoid thermal runaway from shock heating, allowing the atmosphere to remain close to its initial profile. This term has the form of an exponentially weighted moving average, where the degree of weighting is set at $\alpha_{w}=0.05$. When the cooling term is removed, the results show a $<1\%$ difference within a scale height either side of the $a=c$ layer, our main area of interest.

Exact solutions for mode conversion and transmission are known for the linear model in terms of ${}_{2}F_{3}$ generalized hypergeometric functions (Hansen et al., 2016). Transmission decreases with increasing attack angle between the wave vector and the magnetic field (Schunker & Cally, 2006). According to Núñez (2019), an acoustic shock with zero attack angle should pass through $a=c$ unchanged. On the other hand, if we orient the magnetic field away from vertical, substantial mode splitting will occur and the resulting fast and slow ‘shocks’ should be smoothed throughout this process if Núñez’s analysis is correct.

The fast wave above $a=c$ is predominantly magnetic in nature, and will not reflect due to the ever-increasing Alfvén speed gradient as the wave is exactly vertical. In any case, that is beyond the phenomenon of interest here. The acoustic slow wave also propagates upward unhindered and becomes progressively more aligned with the magnetic field as the plasma beta reduces towards zero, whilst the fast wave becomes ever more transverse with height. This provides a convenient way to distinguish them in $a\gg c$, as well as via their propagation speeds.

Specifically, we plot the transverse and parallel-to-the-field plasma velocities $v_{\perp}$ and $v_{\parallel}$, as well as the density $\rho$ scaled by the equilibrium density $\rho_{0}(z)$. It is important to note that these ‘$\perp$’ and ‘$\parallel$’ components are with respect to the actual (perturbed) magnetic field, not the original field. Since the shock front compresses the plasma vertically by up to a factor of 4, and the field lines are tied to the plasma, the perturbed field inclination is greatly increased.

3 Results

A wave driver is placed at the base $z=-2.0$ Mm, and delivers an angular-frequency-$\omega_{0}$ half-period sinusoidal vertical velocity burst of amplitude $V_{0}$. The driver creates a strictly vertical wave vector. The plasma velocity is also initially vertical, until the inclined magnetic field takes effect as the plasma beta decreases with increasing altitude. For concreteness, we adopt $\omega_{0}=2\omega_{c}=8.0$ mHz. This corresponds to a broad frequency spectrum

[TABLE]

centred at zero, where $f=\omega/2\pi$. However, the acoustic cutoff only lets frequencies above 4 mHz propagate upward. This is confirmed numerically. Using drivers lasting several half-periods allows us to localize the driving frequency more precisely, but introduces other numerical complications, particularly first shocks reflect off the top of the computational domain and complicate behaviour in the region of interest. For cleanest results, the half-period driver has been found to be best. In any case, it makes very clean shocks from a broad spectrum of frequencies above $f_{c}$.

Kalkofen et al. (1994) have discussed the evolution of linear and nonlinear waves in non-magnetic 1D isothermal atmospheres. The wakes produced by the initial pulse follow each other by typically one cut-off period, supporting the theory that the atmosphere rings at the cut-off frequency. We restrict our results to tracking the initial pulse produced by the 8.0 mHz driver.

The velocity perturbation increases exponentially with height until a shock is formed well below the $a=c$ layer. When the magnetic field is vertical, the shock front remains sharp throughout the entire simulation and travels unhindered through to the top of the computational box, as expected (Fig. 1). Increasing $\theta$, the angle of the magnetic field from vertical, produces a splitting of the wave into its slow and fast components around the $a=c$ layer. When $a>c$, the fast wave accelerates out in front of the slow wave and the two wave types can be easily distinguished.

We choose to look separately at velocities both parallel and transverse to the (actual) magnetic field. This separation allows an approximation of the two wave types, however it is not precise. This can be seen by the fast wave appearing in the parallel velocity for $\theta>0$, and vice versa for the slow wave.

We also look at the in-place compression ratio $\rho/\rho_{0}$, where $\rho$ is the density value at the given time and $\rho_{0}$ is the initial background density value at the same position. This is not the same as the classical shock-jump compression ratio $X=\rho_{2}/\rho_{1}$ from the pre-shock to post-shock regions (‘1’ and ‘2’ respectively), though they do correspond immediately after a shock front passes. For an oblique compressive shock, where the magnetic field is not completely perpendicular to the shock front, the limiting compression ratio for $\gamma=5/3$ is 4 (Priest, 1982). For these simulations, this remains true across the shock front as expected; however the values behind the front can exceed this ratio due to the transport of fluid via the propagating shock. Deeper, higher density fluid is carried upwards and is then related to the original density values of greater heights giving the appearance of a much higher compression ratio. To eliminate this, the movement of each fluid particle would need to be traced throughout the simulation and the densities compared to those at their original equilibrium height, which is an unnecessary complication.

A non-zero attack angle gives rise to the splitting of the fast and slow components around $a=c$. However both the fast and slow shock remain sharp in the simulations when $\theta$ is less than about $15^{\circ}$. (In reality, we believe that the slow shock will smooth for any non-zero attack angle, but this is computationally imperceptible at these very narrow incidences.) After this point, we begin to see smoothing of the slow component as it separates from the fast shock and traverses the layers around $a=c$. As the slow wave exits the $a=c$ layer, it begins to steepen again to ‘re-shock’, best seen in Fig. 2. Conversely, the fast shock remains sharp throughout and continues to propagate through the atmosphere unaltered in its form. It does not reflect because it is vertical.

Figures 2, 3, and 4 display the evolution of the $\theta=15^{\circ}$, $30^{\circ}$, and $60^{\circ}$ cases respectively. The greater the value of $\theta$, the more smoothing of the slow component is seen. We quantify this by measuring the steepness of the slow shock at its inflection point. The blue dots shown on the ‘shock’ front in each case is where $v_{\parallel}$ has its inflection point, i.e., where the ‘shock’ is at its steepest. The value of $h_{\parallel}$ given above each panel is the ‘velocity scale height’ there and is defined as

[TABLE]

which is a measure of shock steepness. A small $h_{\parallel}$ corresponds to a steep shock front (numerically limited).

Núñez (2019) assumes the passage through the conversion zone acts linearly on the waves, which is based on the supposition that they spend very little time in that area. However, the incoming shock front acts to drag the $a=c$ layer along with it (Fig. 5), which can increase in height by up to 400 km for low $\theta$. The shock front itself can spend multiple seconds within close proximity of the $a=c$ layer.

4 Discussion and Conclusions

1.5D simulations were conducted using a bottom-driven 8.0 mHz half-period pulse injected into an isothermal, gravitationally stratified atmosphere with the angle the magnetic field makes from vertical being varied from $\theta=0^{\circ}$ to $60^{\circ}$. Only power above the acoustic cutoff of 4 mHz propagates upward. The perturbation was allowed to freely evolve, forming a shock well below the Alfvén-acoustic equipartition ($a=c$) layer. If the attack angle $\alpha$ between the (vertical) wave vector $\mathbf{k}$ and magnetic field is non-zero, the shock splits into its slow and fast components around $a=c$. For $\theta\gtrsim 15^{\circ}$, the slow-shock in our simulations shows smoothing through this area, then begins to re-shock upon exiting. The fast-shock, however, continues to propagate freely throughout the atmosphere regardless of the value of $\theta$. This is in only partial agreement with Núñez (2019), who stated that both the fast and slow components of the shock should be smoothed. A possible reason for the discrepancy is that the $a=c$ equipartition region is advected with the shock for several seconds, making the interaction intrinsically nonlinear.

However, it is not clear that Núñez’s result is correct even without this effect. Figure 6 shows the exact linear transmission coefficient $\mathcal{T}$ for a vertical wave in an isothermal atmosphere (Cally, 2009) as a function of frequency with the four field inclinations discussed above. The complement, $\mathcal{C}=1-\mathcal{T}$ is the conversion coefficient, i.e., the fraction of energy going into the fast wave above $a=c$. Of course, it starts at the cutoff frequency $\omega_{c}$. Any shock wave discontinuity is necessarily constructed from high-order Fourier modes with amplitudes asymptotically proportional to the inverse of the vertical wavenumber $k_{z}$. In turn, $k_{z}$ is proportional to $\omega$ through a dispersion relation. Figure 6 shows that, especially at greater magnetic field inclinations, high-frequency high-wavenumber modes are not transmitted to any appreciable degree, thereby smoothing the incident acoustic shock. On the other hand, these small-scale Fourier modes are essentially perfectly converted to fast waves, which are therefore sharp. This is in accord with our simulations. The advection of the $a=c$ level with the shock could only amplify this effect.

The same conclusion follows from the approximate WKB-based formula for transmission (Schunker & Cally, 2006),

[TABLE]

where $h_{s}=[d(a^{2}/c^{2})/ds]_{a=c}^{-1}$ is the thickness of the conversion layer along the direction of the wave vector $\mathbf{k}$, and $\alpha$ is the attack angle between $\mathbf{k}$ and ${\mathbf{B}}$. In our case $\alpha=\theta$, $s=z$, and $\mathcal{T}\to 0$ exponentially with increasing $|\mathbf{k}|$, faster at greater $\theta$. Again, this is in accord with our findings. The advection of $a=c$ with the shock may effectively increase $h_{s}$ and thereby decrease $\mathcal{T}$ and increase $\mathcal{C}$, as suggested above.

A further perspective is afforded by the jump relations for MHD shocks in uniform media, and in particular the so-called shock adiabatic, a bi-cubic polynomial for the shock speed $v_{1}$ dependent on adiabatic index $\gamma$, magnetic field inclination $\theta_{1}$ in the pre-shock plasma, shock compression ratio $X=\rho_{2}/\rho_{1}$ (where $1<X<(\gamma+1)/(\gamma-1)$, as for hydrodynamic shocks), and pre-shock sound and Alfvén speeds $c_{1}$ and $a_{1}$ respectively (Fitzpatrick, 2014):

[TABLE]

The post-shock region is denoted by label ‘2’. For $\gamma=5/3$, compression is limited to $1<X<4$. This should give a good local description of shocks in a stratified atmosphere.

Consider the case $\gamma=5/3$, $\theta=30^{\circ}$, corresponding to Fig. 3. Adiabats for our canonical sound speed $c_{1}=8.958$ km $\rm s^{-1}$ and four Alfvén speeds are plotted in Fig. 7. For each Alfvén speed case, the upper curve corresponds to the fast shock. The loops beneath correspond to the slow shock (lower branch) and the intermediate or Alfvén shock (upper branch). The intermediate shock is not excited in our simulations. Clearly, the fast shock exists for each $1<X<4$, but the slow shock is absent if $X$ is to the right of its loop. As $a_{1}$ increases with height in the stratified atmosphere, the intermediate/fast loop extends to greater $X$ and the shock redevelops. For example, at the height where $a_{1}=12$, a slow shock may only exist for $X\lesssim 1.53$, but for $a_{1}=35$ there is a slow shock for all $X\lesssim 2.99$.

This is consistent with our simulations (Fig. 3), in that the compression is large across the slow wave front in the bottom two panels, where the somewhat-steep slow wave has passed beyond $a=c$. This means there is no available slow shock solution, providing another explanation for its absence in the simulations.

On the other hand, the fast shock always has a solution. Because $X$ is close to 1 (and gets closer with height) across the fast shock in the bottom panels, the fast shock travels at a speed only just above the Alfvén speed $a_{1}$.

The reason for the slow wave above $a=c$ inheriting most of the density contrast is difficult to see from the exact mathematics, but it is plausible conceptually. The slow wave above the equipartition height is the primarily-acoustic successor to the incident fast shock below $a=c$, which is also mainly acoustic. It is therefore unsurprising that density is largely tied to these waves. On the other hand, the fast wave in $a>c$ is more weakly coupled to density, and in fact becomes incompressive in the $a\gg c$ limit. Irrespective of these broad observations, the numerical solutions do suggest that the slow wave in $a>c$ accounts for most of the density contrast, and therefore cannot shock until $a$ increases sufficiently for a shock solution to become available (see the 35 km $\rm s^{-1}$ lobe in Fig. 7).

Implications for solar and stellar chromospheres might include reduced acoustic wave heating near and above the equipartition level, both because the acoustic shock is smoothed and because energy is removed from the acoustic wave. Any enhancement of energy in the fast shock due to equipartition-level advection might benefit Alfvén wave production by fast/Alfvén conversion higher in the atmosphere.