Philosophical Transactions of the Royal Society A: Mathematical, Physical and Engineering Sciences
Open AccessResearch articles

Characterization of shock wave signatures at millimetre wavelengths from Bifrost simulations



    Observations at millimetre wavelengths provide a valuable tool to study the small-scale dynamics in the solar chromosphere. We evaluate the physical conditions of the atmosphere in the presence of a propagating shock wave and link that to the observable signatures in mm-wavelength radiation, providing valuable insights into the underlying physics of mm-wavelength observations. A realistic numerical simulation from the three-dimensional radiative magnetohydrodynamic code Bifrost is used to interpret changes in the atmosphere caused by shock wave propagation. High-cadence (1 s) time series of brightness temperature (Tb) maps are calculated with the Advanced Radiative Transfer code at the wavelengths 1.309 mm and 1.204 mm, which represents opposite sides of spectral band 6 of the Atacama Large Millimeter/submillimeter Array (ALMA). An example of shock wave propagation is presented. The brightness temperatures show a strong shock wave signature with large variation in formation height between approximately 0.7 and 1.4 Mm. The results demonstrate that millimetre brightness temperatures efficiently track upwardly propagating shock waves in the middle chromosphere. In addition, we show that the gradient of the brightness temperature between wavelengths within ALMA band 6 can potentially be used as a diagnostics tool in understanding the small-scale dynamics at the sampled layers.

    This article is part of the Theo Murphy meeting issue ‘High-resolution wave dynamics in the lower solar atmosphere’.

    1. Introduction

    The solar atmosphere is highly dynamic at small scales at chromospheric heights, also under quiet Sun conditions with low magnetic-field strength [1]. A major contribution to the small-scale dynamics comes from the propagation of shock waves. Acoustic waves propagating upwards from the solar surface steepen into shock waves as a result of the decrease in gas density with height. The formation of shock waves and their propagation through the atmosphere have been studied by means of detailed one-dimensional (1D) numerical simulations (e.g. [215]).

    Three-dimensional (3D) simulations, for instance those by Wedemeyer et al. [16], exhibit a dynamic mesh-like pattern of hot filaments from shock waves surrounding cooler post-shock regions. Such 3D simulations are also employed by Wedemeyer-Böhm et al. [17] and Loukitcheva et al. [18] to explore the use of millimetre and submillimetre wavelengths as diagnostic tools for the chromosphere.

    The complex dynamics in the chromosphere have also been reported in numerous observational studies. In particular, small-scale (of the order of 1.5 Mm and smaller) and short-lived (≈100 s or less) bright structures have been observed in the Ca II H and K spectral lines (the so-called H2V/K2V grains), in agreement with the mesh-like patterns seen in the simulations [1922], which are often interpreted as shock signatures [23]. While [24] reported a correlation between these small-scale structures and a weak magnetic field of about 20 G, Lites et al. [25] and Cadavid et al. [26] found no clear relationship between magnetic fields and the K2V grains. Other observational studies of small-scale shock signatures have also found that the magnetic field activity and orientation may play a major role in quiet Sun regions [27,28], where shock waves propagate in both weak (or non-magnetized) and strong field-concentrated regions. Hasan & van Ballegooijen [29] showed that shock waves can produce excess temperatures of about 900 K in small magnetic concentrations in the chromosphere, which is responsible for the excess brightness observed in, e.g. small Ca II H magnetic bright points [30].

    The radiation observed at millimetre wavelengths originates at chromospheric heights from free-free emission in local thermal equilibrium (LTE; see e.g. [31] and references therein). Following the Rayleigh–Jeans Law (i.e. [32]), the measured intensities depend linearly on the local plasma temperature. In order to detect the small-scale dynamics in observational mm-wavelength data, high spatial and temporal resolution is essential. In that regard, the Atacama Large Millimeter/submillimeter Array (ALMA), which offers regular observations of the Sun since 2016, represents a major step forward in terms of resolution, and has already provided insights into the dynamics of mm-wavelength intensities, e.g.[3339]. Modelling has shown that propagating shock waves will cause variation in mm-wavelength observables [13,14]. However, these studies employed 1D models, therefore it is uncertain to what degree the intensity variations are affected by more realistic interactions of shock waves in a 3D environment.

    In this work, we use a 3D model that takes into account essential physical processes such as non-LTE and non-equilibrium hydrogen ionization that have a large impact on the mm-wavelength radiation. With support from more realistic 3D simulations, it is possible to connect the mm-wavelength observables to the underlying physics and, thus, determine and characterize the observable signatures in mm-wavelengths of propagating shock waves.

    The structure of the work is as follows. In §2, the setup of the simulations is explained and in §3 a representative example of a propagating shock wave with the surrounding physical conditions and the resulting signatures in brightness temperature are presented. In §4, we discuss the results and in §5 conclusion is drawn and motivation for future work is given.

    2. Simulation set-up

    A 3D numerical model of the solar atmosphere is created with the radiative magnetohydrodynamic (MHD) code Bifrost [40,41]. The duration of the considered simulation sequence is approximately 1 h with an output cadence of 1 s (matching the highest ALMA cadence in solar mode), so that rapid small-scale events on scales down to a few seconds can be efficiently studied. The simulation box has an extent of 24 × 24 × 17 Mm in horizontal (x, y) and vertical (z) directions, respectively. The number of cells in x and y are both 504, with a constant grid spacing of 48 km (corresponding to approx. 0.066″ at 1 AU). In the z-direction, there are 496 cells with grid spacing varying between 19–100 km, with a spacing of 20 km at chromospheric heights. In the horizontal directions, the boundary conditions are periodic. The bottom boundary, which is located 2.5 Mm under the photosphere, allows flows (e.g. intergranular downdrafts) through, however, the average horizontal pressure is driven towards a constant value at a characteristic time of 100 s. This gives rise to acoustic wave reflection, mimicking the refraction of waves in the solar convection zone and giving rise to p-modes in the simulation. The upper boundary condition is based on characteristic equations and allows for the transmission of magneto-acoustic waves.

    The simulation takes into account non-LTE and non-equilibrium hydrogen ionization, as hydrogen is the major contributor to the number of free electrons. Ionization of other elements are under the assumption of LTE and, thus, given as function of internal energy and total mass density. The simulation represents an ‘enhanced network’ region surrounded by quiet Sun [40]. The magnetic-field strength has an unsigned average value of 50 G (5 mT) at photospheric heights with two opposite polarity regions of magnetic field approximately 8 Mm apart.

    The magnetic field was introduced into a relaxed convection simulation as two patches of opposite polarity at the lower boundary, with the rest of the simulation volume filled with a potential field extrapolation. This initial configuration is quickly modified by the convective flows that sweep the field into the intergranular lanes. The convective flows also cause foot point motions of the loops connecting the two polarities, leading to heating of the chromosphere and corona. The analysis of the simulation is restricted to later times (after 200s), when quasi-equilibrium has been established. For details of the set-up and thermodynamic evolution of the simulation, see [40].

    The same simulation with a cadence of 1 s has been used previously in the study by Wedemeyer et al. [36], to determine the diagnostic potential of solar ALMA observations. Other versions of the same simulation set-up, but with lower cadence have been used in [18,42], where the authors study chromospheric diagnostics at mm and sub-mm wavelengths with focus on the thermal structure and magnetic field.

    The observable intensity at mm-wavelengths is obtained by solving the radiative transfer equation column by column. In LTE, the Rayleigh–Jeans approximation (valid for mm wavelengths) gives

    where TB is the brightness temperature, χν is the opacity, Tg is the gas temperature and τν is the optical depth at height z defined from

    The integrand in equation (2.1) is the contribution function, describing which regions along the line-of-sight that contribute to the observed brightness temperature. We use the Advanced Radiative Transfer code to solve the radiative transfer equation. The code assumes LTE but includes in detail the relevant sources of continuum opacity. For ALMA wavelengths, the opacity is dominated by free–free processes of hydrogen and H [13].

    Two specific wavelengths are used in this study, located at the respective sides of the ALMA spectral band 6 in solar observing mode. The ALMA receiver band 6 is sub-divided into four sub-bands: SB1 (1.298–1.309 mm), SB2 (1.287–1.298 mm), SB3 (1.214–1.224 mm) and SB4 (1.204–1.214 mm). The wavelengths used here are 1.309 mm (229.0 GHz) and 1.204 mm (249.0 GHz), which are at the edges of SB1 and SB4, respectively. The brightness temperature, Tb(ν), is calculated from the radiative intensities, I(ν), through the Rayleigh–Jeans Law approximation.

    Figure 1 shows the brightness temperature for SB1 at a time of 200 s from the start of the simulation. Over the whole approximately 1 h of simulation time, there are many signatures associated with shock waves. An example of a shock wave event, further described in detail in §3, is located in the white circle marked in figure 1, at (x, y) = (21.07, 0.72) Mm. The event is visualized in more detail in figure 2. This example is representative of a typical shock wave in the simulation box, whose characteristics are presented from a qualitative point of view.

    Figure 1.

    Figure 1. Brightness temperature in SB1 (1.309 mm; 229 GHz) at a time of 200 s after the beginning of the simulation. The white circle marks the location of the selected shock wave which is shown in more detail in figure 2. (Online version in colour.)

    Figure 2.

    Figure 2. Gas temperature surrounding the example of shock wave formation for five different time steps, from bottom to top: t = 150, 183, 196, 210 and 230 s. First column: horizontal cuts at a fixed height of 1 Mm above the photosphere. The white cross marks the coordinates of the shock wave, (x, y) = (21.1, 0.72) Mm, which propagates largely in a vertical direction at these coordinates. Second and third column: vertical cuts through this position along the x- and y-axes, respectively, showing heights between approximately 0.4 and 2.3 Mm. The white dotted lines marks the respective coordinate of the shock wave example.The blue dots show the height of optical depth τ = 1.0 for SB1. (Online version in colour.)

    3. Example of shock wave

    The regions surrounding the shock-wave event can be seen in the horizontal and vertical cuts of the gas temperature along the z-, x- and y-axes in figure 2. The vertical cuts are given for five time steps spread out over the time span of the shock event. The bottom row shows the pre-shock phase at t = 150 s, dominated by cooler down-flowing gas. In the second row from the bottom, at t = 182 s, the shock has formed and reached a height of approximately 1 mm. In that moment, the brightness temperature is already at half its maximum value.

    The third row shows the peak phase of brightness temperature at t = 196 s, where the shock front reaches a height of approximately 1.25 Mm. The fourth row from the bottom shows the time t = 210 s when the brightness temperature is at maximum for this location. At this time, the shock wave front has reached above 1.5 Mm. Finally, in the top row at t = 230 s, the shock wave front has reached well above 2 Mm and the cool post-shock medium is evident at formation heights of SB1 and SB4 around approximately 0.8 Mm.

    By comparing the cuts between the time steps, it is possible to see the propagation of the shock wave to a height of z ≈ 2.3 Mm. The shock wave propagation is predominantly vertical along the z-axis, except for a small tilt along the y-axis that can be seen in the rightmost column of figure 2.

    The ambient medium is highly dynamic, with a complex structure that leads to interactions between events. For instance, as the shock wave propagates through the chromosphere, background waves and structures affect the shock front. The pre-shock medium, here seen around t = 150 s (figure 2), is the resulting post-shock medium from preceding shock waves. Thus, the evolution of the shock wave front is dependent on how previous shock events influence the atmosphere. For this reason, it is necessary to use a realistic 3D atmospheric model in place of 1D models in order to make predictions of the mm-wavelength signatures. The mm-wavelength radiation (marked in blue for 1.204 mm in figure 2) is also dependent on the local atmospheric conditions, with a range of formation heights between approximately 0.6 and 2 Mm.

    (a) Contribution function to brightness temperature

    The time-dependent contribution function of the brightness temperature of SB1 (1.309 mm) of the selected shock wave (see the location marked in figure 1) is given in figure 3. The corresponding contribution function for SB4 (1.204 mm), is similar, however, there are small differences resulting in differing heights of optical depth τ = 1.0, marked by the blue and green dots in figure 3. A value for the optical depth of τ = 1.0 is often a good one-point approximation to the full contribution function. As the shock wave front propagates upwards between the heights z1.11.3Mm (at t190210s), the span of the contribution function is very narrow. Here, the formation height that τ = 1.0 corresponds to shows a strong correlation to, and effectively samples, the local gas temperature at a thin layer around the shock front. Immediately before t = 190 s, but mostly evident after t = 210 s, the brightness temperatures are sampled from several components at different heights. As the shock front propagates above z ≈ 1.5 Mm, it is almost completely transparent in mm-wavelengths. There is a small fraction of the total contribution, no more than a few per cent, that comes from the shock front at these heights, visible as a light grey streak in figure 3. In the pre- and post-shock regimes, the contribution function spans a larger extent of heights.

    Figure 3.

    Figure 3. Time evolution of the contribution function for SB1 at 1.309 mm. (The contribution function for SB4, at 1.204 mm, differs slightly but looks nearly visually identical at these scales.) For each time step, the heights of τ = 1.0 for SB1 and SB4 are indicated by the blue and green markers, respectively.The integrated contribution function is normalized to 1.0. (Online version in colour.)

    (b) Gas and brightness temperature

    The time evolution of the gas temperature during the propagation of the selected shock wave is shown in figure 4a for heights between z0.61.8Mm. In addition, the heights where the optical depth, τ, is unity at the wavelengths 1.309 mm (SB1) and 1.204 mm (SB4) are marked as a function of time. The wave steepens into a shock around t = 175 s, close to a height of z = 0.8 Mm and thereafter shows a very rapid increase in height. This height for shock wave steepening is consistent with other studies (e.g. [16]). The formation height varies largely with time as the shock wave propagates through the chromosphere. There are small differences in the formation height between SB1 and SB4, of up to approximately 40 km with a median value of 20 km, although they follow the same trend and keep the same order. As can be seen in figure 4a, the height where τ = 1.0 for both SB1 and SB4 increases from a pre-shock minimum of z ≈ 0.68 Mm to a peak value of approximately 1.38 Mm during the course of 44 s with the propagation of the shock wave. The formation heights (i.e. z(τ = 1)) thus increase by Δz ∼ 0.7 Mm from the low to the middle chromosphere during the upward propagation of the shock wave.

    Figure 4.

    Figure 4. Time-dependent temperatures in the selected shock wave example. (a) Evolution of gas temperature for one column at chromospheric heights between approximately 0.6 and 1.8 Mm. The blue and green dots mark the formation heights (i.e. τ = 1.0) for the wavelengths 1.309 mm (SB1) and 1.204 mm (SB4), respectively. (b) Evolution of the brightness temperatures (Tb) at wavelengths of 1.309 mm (SB1, blue solid) and 1.204 mm (SB4, green solid) and of the gas temperatures (Tg) at heights corresponding to the optical depth of unity at the respective wavelengths (blue/green dotted). The horizontal dashed black line represents the Full Width Half Maximum (FWHM) lifetime of the brightness temperature shock wave signature (with respect of the base temperature, dotted horizontal line).(c) Evolution of the difference between the brightness temperatures of SB1 and SB4. (Online version in colour.)

    From this point on, the formation heights no longer follow the upward propagating shock front. Instead, after about 10 s at the peak height, the formation heights rapidly decrease to z = 0.72 Mm in just 27 s and thus map the post-shock phase. The brightness temperatures map the hot propagating shock wave front up to a certain height where it decouples as a result of the lower opacity at these heights for the mm-wavelengths of SB1 and SB4.

    Figure 4b shows the corresponding evolution of the brightness temperatures of SB1 and SB4 as well as the gas temperature at the specific heights z(τ = 1.0) for SB1 and SB4, respectively. The brightness temperature of SB1 shows a total increase of approximately 4970 K in Δt = 86 s, starting at a pre-shock local minimum of 3660 K at t = 110 s and rising to the peak value of 8780 K at t = 196 s. Thereafter, there is a rapid reduction to 3800 K at the local minimum at t = 256 s in the post-shock phase. The pre- and post-shock temperatures are thus comparable in this example. The local minimum with the highest temperature (here the post-shock minimum) is referred to as the ‘base temperature‘. The time between the two local minima is 146 s. Estimating the lifetime of the observable brightness temperature signature as the Full Width Half Maximum (FWHM) of the peak results in tFWHM = 41 s. The brightness temperature excess and lifetime shown here are in line with values derived from shock wave propagation in 1D simulations [13,14,43]. The resulting strong correlation between the gas temperature at z(τ = 1.0) and the brightness temperature (figure 4b) is expected for mm-wavelength radiation (i.e.[18]).

    Figure 4c shows the time evolution of the difference between the brightness temperatures of SB1 and SB4. In the pre-shock phase around 85–170 s, there is a difference of down to −74 K, which corresponds to a magnitude of about 2% of the total brightness temperature. Later during the shock phase, the difference between SB1 and SB4 increases to a total of 300 K, corresponding to 4% of the total Tb, and finally decreases to approximately zero in the post-shock phase. The propagating shock wave and the pre-shock epoch display two different cases where the temperature gradient between the two sub-bands have a different sign. The brightness temperatures sampled at the propagation of the shock wave front show a positive gradient. That is, SB1 always forms in higher regions than SB4 and thus has a higher temperature. By contrast, during the pre-shock epoch, the gradient is negative and SB1 shows a lower temperature than SB4. The peaks in the time evolution of the brightness temperatures difference (figure 4c) centred around approximately 185 s, 200 s and 230 s originates from the signatures of three distinct wave components with differing propagation speeds. These are seen in the t − z plot of figure 4a as a hot, rapid component, followed by a cooler, slower component going upwards above 1.8 Mm and a third more diffuse component deflecting downwards around 1.2 Mm.

    The rapid and large variations of the gas temperature at z(τ = 1.0) (figure 4b) are clearly connected to the large variations in formation height. The temporal profile of the brightness temperature (figure 4b) is integrated over a span of heights along the specific column and is therefore smoother than the gas temperature, which is sampled at a single height.

    (c) Vertical velocity

    The evolution of the vertical velocity (vz) at the chosen position (cf. figure 1) is shown in figure 5a. There is a bulk downflow of cooler gas (figure 4a) in the pre-shock region, with velocities of up to 20 km s−1. The shock front is met by this downfall and, therefore, experiences a resistance to its motion as it propagates upwards. In the height range from where the mm continuum radiation in SB1 and SB4 originates, the vertical velocity only reaches a maximum velocity vz ∼ 10 km s−1, whereas in the upper chromosphere, there are velocities of up to approximately 20 km s−1. In figure 5a, the markings v1, v2 and v3 point out the shock front at three different heights, approximately 1.0, 1.3 and 1.65 Mm. At these heights, the t − z slope of the sharp transition of the vertical velocity (or the gas temperature in figure 4a) indicates a speed of the vertical propagation of the shock wave of approximately 33, 19 and 83 km s−1, respectively.

    Figure 5.

    Figure 5. Time-dependent vertical velocities for the selected shock wave example. (a) Evolution of vertical velocity for the chosen column (cf. figure 1) at chromospheric heights between approximately 0.6 and 1.8 Mm. Positive velocity (blue colour) indicates outwards motion away from the photosphere. The markings v1, v2 and v3 show the vertical propagation speed of the shock front at different heights. The colour scale is saturated on the positive side (from 23.7 to 20.1 km s−1). The blue and green dots mark the formation heights (i.e. τ = 1.0) for the wavelengths 1.309 mm (SB1) and 1.204 mm (SB4), respectively.(b) Evolution of vertical velocity at heights corresponding to the optical depth unity at the wavelengths for SB1 and SB4. (Online version in colour.)

    The indications of differences in vertical velocity between heights of z(τ = 1.0) for SB1 and SB4 are generally small (i.e. smaller than 1 km s−1). This is a due to the height difference between z(τ = 1.0) for SB1 and SB4 being of the same order as the vertical grid spacing of 20 km.

    (d) Electron density

    The time evolution of the electron density for the column of the shock wave is given in figure 6. During the pre-shock phase, the electron density is slowly decreasing in which the formation height of the brightness temperatures follows the same trend. There is a larger decrease followed by a rapid increase in the electron density approximately around 180 s as a result of the ionization coming from the shock wave. In the post-shock regime, the electron density decreases slowly which combined with the declining gas temperature (figure 4) at some point (t = 225 s) results in a sufficiently low opacity to initiate a sudden jump down in formation height of the mm-wave radiation. The correlation of the gas temperature and electron density in shocks [44] is confirmed by the close relation of these quantities during the shock wave passage between 1.0 and 1.3 Mm.

    Figure 6.

    Figure 6. Time-dependent electron density for the selected shock wave example. (a) Evolution of electron density for one column at chromospheric heights between approximately 0.6 and 1.8 Mm. The blue and green dots mark the formation heights (i.e. τ = 1.0) for the wavelengths 1.309 mm (SB1) and 1.204 mm (SB4), respectively. The dotted and dashed horizontal lines marks the heights of 1.2 and 1.7 Mm. (b) Evolution of electron density at heights corresponding to the optical depth of unity at the respective wavelengths of SB1 and SB4.(c) Evolution of electron density at fixed heights of 1.2 and 1.7 Mm. (Online version in colour.)

    The local atmosphere shows an increased electron density for a significant time span after the shock wave propagates through. [44] find through 1D models that the time scale of relaxation of the local atmosphere through hydrogen recombination after the shock wave has propagated through the chromosphere varies with height. They show a large span of time scales of the order of between approximately 102 and 105 s at chromospheric heights, with a peak value in the mid-chromosphere and rapidly decreasing towards both the photosphere and transition region. In figure 6c, the time evolution of the electron density at two fixed heights, 1.2 and 1.7 Mm are shown. The rate of decrease of electron density is slower at 1.7 Mm than at 1.2 Mm, for this shock wave. The relaxation times are difficult to measure due to the dynamic atmosphere with preceding, as well as succeeding, wave trains at the same position, ensuring that the electron density never reaches a steady state. Estimating the relaxation times by simply extrapolating with the same trend as for the last 30 s between 270 s and 300 s, results in approximately 200 s and 430 s to reach values of previous minima at the heights 1.2 and 1.7 Mm, respectively.

    The difference of the electron density between SB1 and SB4 is on the order of log (Ne) = 0.1 cm−3. The vertical grid spacing is in the same order as the differences in formation height z(τ = 1) of SB1 and SB4. Therefore, as with the velocities in §3c, it is difficult to make use of the differences between electron density of SB1 and SB4 that are sampled at z(τ = 1).

    The opacity of the mm-wavelengths is dominated by free–free processes (e.g. [13,45]) and is reduced by an increase in temperature coupled with a decrease in electron density. This behaviour is seen as the shock wave propagates upwards, resulting in a decrease of the opacity. At some point, the opacity is reduced to the level where the brightness temperature de-couples from the shock wave front, leading to a rapid decrease in the formation height z(τ = 1).

    4. Discussion

    The shock wave example illustrated in this work was found to be representative of a typical shock wave found in the simulation with respect to variations of gas temperature, vertical velocity and electron density. As a result of the complex 3D structure of the atmosphere, the propagation of the shock waves can be intricate, with differing speeds in different directions, alongside changes in the propagation angle of the wavefront. The shock wave example presented in this study exhibits a predominantly vertical propagation, which is to be expected for the region under consideration, with magnetic fields of minimal inclination. The temporal profile of a shock wave event will show deformations depending on the specific propagation properties at different locations (i.e. different columns). Thus, the vertical propagation of the shock front under consideration ensures a simple inference of the brightness temperature, absent of any large 3D components.

    The formation height of mm-waves in the chromosphere is a function of wavelength, with height increasing with wavelength. Accordingly, SB1 forms slightly higher up than SB4. This is seen throughout the entire example of shock wave propagation in the chromosphere (figure 4a), despite the presence of small-scale dynamics in the chromosphere. The difference of 0.105 mm between SB1 and SB4 is shown to be enough to map different layers with brightness temperature differences up to approximately 300 K (figure 4c). There are however notable variations in the differential between SB1 and SB4, sometimes to such a high degree that a reversal of the sub-bands are evident. The strongest reversal is seen in connection to the down-falling cold gas seen in the pre-shock region (approx. t = 140 s) of the illustrated example.

    With larger wavelength separation between the sampled sub-bands, larger Tb differences would be observed as the sampled layers would lie further apart. For that reason, it would be of interest in regards of the temporal domain to track a propagating shock wave from one layer to the other. Observations with ALMA with increased separation between the sub-bands would be favourable. Solar ALMA observations are currently offered in several spectral bands between approximately 0.8 and 3.3 mm. ALMA band 3 (2.83.3mm) offers the largest default separation between the outermost sub-bands of 0.42 mm, which is approximately 4.5 times more than that of band 6 in this study. However, further consideration needs to be made, such as the change in formation height and, thus, potential change of shock wave propagation speed and contrast of the dynamic signatures.

    (a) Contribution function spread over geometric height

    Figure 3 reveals the important relationship between the developing shock and the associated τ = 1.0 region. Importantly, it can be seen that during the initial formation of the shock (approx. 180 s in figure 3), the plasma is both dense and bright, resulting in the observed signatures at the τ = 1.0 location being dominated by the shocked plasma. Initially, this relationship continues to hold as the shock develops and propagates into higher layers of the chromosphere. However, at t ≈ 220 s in figure 3, it can be seen that the contribution function defining the τ = 1.0 surface begins to decouple from the upwardly propagating shock. At this point, the contribution function will be comprised both shocked plasma expanding upwards into more diffuse and optically thin layers of the atmosphere, in addition to cooling plasma beginning to accelerate back towards the solar surface, which is visible in figures 4a and 5a. This results in the contribution function being spread over multiple components spanning a vast assortment of geometric heights. At this point, the signatures extracted at the τ = 1.0 height no longer strictly correspond to the propagating shock, which has important implications for observational studies of such phenomena. For example, recent work surrounding shocks manifesting in sunspot umbrae [46,47] have described the challenges faced when interpreting the spectropolarimetric fluctuations in Stokes I/Q/U/V spectra over the lifetime of a shock. Hence, the opacity response of propagating shocks affects a wide range of observable signals, spanning brightness temperatures in the radio regime through polarimetric intensities across the visible and infrared spectrum. In particular, recent work by Houston et al. [47] interpreted reversals in the Stokes Q/U spectra as evidence for the presence of intermediate shocks, but this interpretation relied upon the observed signals being closely coupled to the developing shock front. As such, future investigations of challenging shock signatures (e.g. [4850]) need to carefully consider the potential effects of the contribution function decoupling from the shock as it propagates into less opaque regions of the solar atmosphere.

    (b) Observations compared to numerical simulations

    Large advances have been made in interferometric observations of the Sun in mm-wavelengths, with ALMA. There are, however, many challenges that come with solar ALMA observations, for example, but not limited to, image reconstruction of interferometric data, limited spatial resolution, absolute temperature measurements, atmospheric noise, etc., that are out of the scope of this work.

    The sophisticated 3D simulations give realistic predictions of how the shock wave signatures would look in brightness temperature as observed at mm-wavelengths. This work considers the brightness temperatures calculated at a horizontal resolution element of 48 km (approx. 0.066 arcsec). Performing actual observations at these wavelengths comes with instrumental resolution limitations that need to be considered. Though further complications arise from the limited spatial resolution of observations, the results of this work point towards the optimistic capability of ALMA with highly resolved data. As a result of limited spatial resolution, the magnitude of the excess Tb of the dynamic profiles will appear less strong due to the pixel filling factor comprising of both the shocked plasma and cooler, quiescent plasma.

    The estimated difference in formation height between the wavelengths 1.309 (SB1) and 1.204 mm (SB4) of up to 40 km with a median value of 20 km is on the limit of the vertical grid spacing of the simulation. The vertical resolution of 20 km in the chromosphere puts a limit on the differences of the small-scale dynamics that can be handled. To study the differences between the sub-bands in more detail, a higher resolved numerical model would be necessary.

    Observations also come with a certain amount of noise. The signal-to-noise ratio needs to be high enough to accurately deal with the magnitude of Tb variations of interest. A few studies of ALMA data have been made where Tb variations of small-scale structures have been reported [3436]. In these studies, the brightness temperatures of the full spectral bands (all four sub-bands combined) were used. Integration over larger spectral or temporal spans can be done in order to increase the signal-to-noise ratio. To accurately map the Tb differences between two sub-bands introduces a larger uncertainty. There are studies where the sub-bands are successfully used separately [51,52], which in this case acts as a proof of concept. Detection of brightness temperature variations as small as 70 K has been reported by Nindos et al. [35], where they use ALMA observations at approximately 3 mm (band 3) with integration over the full band with a cadence of 2 s. The spatial resolution element of their band 3 data (2.5″ × 4.5″), is larger than what currently can be achieved with band 6 data around 1.2041.309mm. This is a direct result of the Fourier sampling (fringe spacing) of the interferometric data scales with wavelength (i.e. [32]). With regards to the ability to spatially resolve small-scale events, the ability to measure precise brightness temperatures should therefore be even more precise in band 6 than in band 3. However, the integration over the full band comes with the inherent loss of sampling different layers as a function of wavelength. There should be an optimal combination of improving upon the signal-to-noise ratio whilst sampling different layers within a spectral band, so that differences of the order of 100 K can be detected, as indicated by the Tb difference between the sub-bands in the simulations (figure 4c).

    Estimating the observable signatures of shock waves in mm-wavelengths with current and potential future modes offered for solar-ALMA observations, including the effect of limited spatial resolution of different spectral bands and the sampling of different physics between the sub-bands will be investigated in a forthcoming paper.

    5. Conclusion

    We use realistic numerical 3D MHD simulations from the Bifrost code, including non-LTE, non-equilibrium hydrogen ionization, of the solar atmosphere to study small-scale dynamics connected to propagating shock waves and how these are perceived in mm-wavelength radiation. An example of a shock wave with nearly vertical propagation and without much interference from neighbouring dynamical features is illustrated. The shock wave propagating upwards in the chromosphere at vertical velocities between approximately 19 and 83 km s−1, and has an associated increase in the local gas temperature of the order of several thousand degrees. We conclude that the brightness temperature at mm-wavelengths corresponding to ALMA band 6 (1.2041.309mm) probes these gas temperatures accurately under the highly dynamical conditions arising from propagating shock waves. The gas temperature at a single height z(τ = 1.0) is quite close to the brightness temperature, which demonstrates the close relationship and the diagnostic potential for determining actual gas temperatures from mm-wavelength observations. The FWHM lifetime of the Tb shock wave signature is 41 s.

    The formation height of the radiation at a certain wavelength is not fixed. The formation height of wavelengths 1.2041.309mm varies of the order of approximately 0.7 Mm, from approximately 0.7 Mm to 1.4 Mm in less than a minute, in the course of the shock wave propagating through the chromosphere. The brightness temperatures at wavelengths corresponding to ALMA band 6 at 1.2041.309mm efficiently maps the shock front while it is propagating from approximately 1.0 Mm up to 1.4 Mm, where the brightness temperatures start to decouple and instead starts to map the post-shock region. The shock wave front continues propagating upwards, unseen by radiation at 1.2041.309mm. In the pre- and post-shock regimes, the radiative contribution function at these wavelengths is more diffuse and spread out over larger span of heights. At some instances, the brightness temperature (at one frequency) has contributions from distinct layers at different heights. This is the scenario right before and after the strong coupling of the brightness temperature with the shock front.

    There is a wavelength dependency of the optical depth which has been explored for wavelengths lying in the furthest apart sub-bands of ALMA spectral band 6. The simulations indicate that the difference in formation height between wavelengths SB1 and SB4 is up to approximately 40 km with a median difference of 20 km. The order of the formation heights with SB1 forming higher up than SB4, is however constant. Because of the correlation between the formation height and the wavelength of radiation, the gradient of brightness temperatures within the spectral band corresponds to a gradient in plasma temperature between the respective formation heights. The brightness temperatures of SB1 and SB4 show differences from about −70 K up to approximately 300 K in the shock wave example. The difference between the sub-bands comes from the local temperature gradient between the mapped layers at the formation heights of the sub-bands. As the brightness temperature is coupled to the shock wave front, SB1 (1.204 mm) has a higher temperature than SB4 (1.309 mm) and there is a positive gradient with increasing temperature with height. In the pre- or post-shock regimes dominated by sampling of cold down flowing gas, the temperature gradient tends to be negative with SB1 colder than SB4.

    The presented simulation results demonstrate that brightness temperatures of wavelengths corresponding to ALMA spectral band 6 (1.2041.309mm) can be used for tracking shock waves from the middle chromosphere and that the gradient of the brightness temperature within the spectral band in principle can be used as a diagnostics tool for probing the small-scale structure of the chromosphere.

    Data accessibility

    The Bifrost simulation with 10 s cadence is publicly available at:

    Authors' contributions

    M.S. and M.C. performed the M.H.D. simulations and radiative transfer computations. H.E. performed scientific analysis, with assistance from S.W., B.S., D.B.J., S.J., S.D.T.G. and M.C. H.E. drafted the manuscript. All authors read and approved the manuscript.

    Competing interests

    We declare we have no competing interests.


    This work is supported by the SolarALMA project, which has received funding from the European Research Council (ERC)under the European Union’s Horizon 2020 research and innovation programme (grant agreement no. 682462) and by the Research Council of Norwaythrough its Centres of Excellence scheme, project number 262622 (Rosseland Centre for Solar Physics). The development of the Advanced Radiative Transfer (ART) code was supported by the PRACE Preparatory Access Type D program(proposal 2010PA3776). Grants of computing time from the Norwegian Programme for Supercomputing are acknowledged. B.S. is supported by STFC researchgrant no. ST/R000891/1. D.B.J. and S.D.T.G. are supported by an Invest NI and Randox Laboratories Ltd. Research & Developmentgrant no. (059RDEN-1).


    D.B.J. and S.D.T.G. are grateful to Invest NI and Randox Laboratories Ltd. for the award of a Research & Development Grant (059RDEN-1). We also acknowledge collaboration with the Solar Simulations for the Atacama Large Millimeter Observatory Network (SSALMON, The support by M. Krotkiewski from USIT, University of Oslo, Norway, for the technical development of ART is gratefully acknowledged. B.S., D.B.J., S.J. and S.D.T.G. wish to acknowledge scientific discussions with the Waves in the Lower Solar Atmosphere (WaLSA; team, which is supported by the Research Council of Norway (project no. 262622) and the Royal Society (award no. Hooke18b/SCTM). The Advanced Radiative Transfer (ART) code was developed by Jaime de la Cruz Rodriguez and Mikolaj Szydlarski and optimised with the help of a PRACE Preparatory Access Type D program (proposal 2010PA3776) supported by M. Krotkiewski from USIT, University of Oslo, Norway.


    One contribution of 16 to a Theo Murphy meeting issue ‘High-resolution wave dynamics in the lower solar atmosphere’.

    Published by the Royal Society under the terms of the Creative Commons Attribution License, which permits unrestricted use, provided the original author and source are credited.