RADIATIVE HYDRODYNAMIC SIMULATIONS OF ACOUSTIC WAVES IN SUNSPOTS

Periodic, short-lived brightening of small regions of sunspot umbrae, known as umbral flashes (Beckers & Tallant 1969), are easily observable periodic phenomenon present in most, if not all, sunspot umbrae. Photometrically, umbral flashes manifest themselves in a rapid increase in the intensity followed by a more gradual decrease seen as asymmetries in the light curve. Spectroscopically, short-lived periodic features in the H and K lines of singly ionized calcium are observable in the chromosphere above sunspot umbrae, even with moderate resolution. The dramatic increase in brightness of the H and K lines is a consequence of increased emission from the local solar material. With the current resolution of the best ground-based telescopes, such as the Swedish 1 m Solar Telescope (Scharmer et al. 2003) with a resolution of 0.1 arcsec, these features can be studied in considerable detail. Thus, temporally and spatially high-resolution spectra of the Ca ii H and K lines provide diagnostic tools to probe the chromospheres above sunspot umbrae.

The observational properties of umbral flashes have been widely reported (Beckers & Tallant 1969; Wittmann 1969; Firstova 1980; Turova et al. 1983; Thomas 1985; Bard & Carlsson 1997; Socas-Navarro et al. 2000; López Ariste et al. 2001; Rouppe van der Voort & Krijger 2003; Rouppe van der Voort et al. 2003; Nagashima et al. 2007). The commonly reported values are brightening relative to the umbra ≈2, lifetimes ≈50 s, diameters ≈3 arcsec, periods ≈140–160 s, and Doppler blueshifts ≈7 km s⁻¹. Also, Beckers & Tallant (1969) noted that the light curves were asymmetric with rapid rise in brightness followed by a more gradual fall and that a rapid redshift occurs just before the blueshifted flash is observed.

In parallel with this continued observational effort, a large body of theoretical investigations has grown to ascertain the nature of the physical processes that cause umbral flashes (Havnes 1970; Moore 1973; Zhugzhda et al. 1983; Lee & Yun 1985; Lites 1992; Zhugzhda 2002) and there is some consensus that they are manifestations of wave motions in sunspots. Further, there is observational evidence (Thomas 1985) in comparing chromospheric lines from different heights that these waves are driven from deeper layers.

The work presented here is motivated by the observational evidence and theoretical efforts outlined above. In our effort to gain a better understanding of the dynamics of the chromospheric regions above sunspot umbrae, we examine the spectroscopic consequences of umbral flashes for the H line of singly ionized calcium.

It is commonly assumed that the mechanism responsible for umbral flashes is sound waves originating in the photosphere and propagating outward, in an increasingly rarefied atmosphere and developing into shock fronts in the chromosphere, similar to the mechanism producing bright grains in the quiet Sun as shown by Carlsson & Stein (1997). In the present article, we test this hypothesis by making detailed radiation hydrodynamic simulations in an umbral atmosphere and investigating how the Ca ii H-line core emission depends on the photospheric power spectrum of sound waves, the coronal incident radiation field, and possible extra mechanical heating in the chromosphere.

The structure of this paper is as follows. In Section 2, we describe the methods we use to solve the relevant equations. In Section 3, we analyze how umbral flashes are formed, and in Section 4 we present our grid of simulations where we have varied input parameters. We present the results of this parameter study in Section 5 and finish the paper with a discussion and conclusions in Section 6.

We assume that the magnetic field in the umbra is strong and vertical; the main effect of the magnetic field is to restrict the plasma motions to the vertical direction. We therefore implicitly solve the one-dimensional hydrodynamic equations of mass, momentum, and energy conservation (ignoring magnetic field effects) together with the non-LTE radiative transfer and population rate equations, on an adaptive mesh. We use the same code and the same methods as in Carlsson & Stein (2002): the radiative transfer is treated using a linearization of the rate equations (Scharmer & Carlsson 1985) but with a banded approximate lambda operator (Rybicki & Hummer 1991) instead of the global Scharmer operator. The advection terms are treated using van Leer's (1977) second-order upwind scheme to ensure stability and monotonicity in the presence of shocks. An adaptive mesh is used (Dorfi & Drury 1987) in order to resolve the regions where the fluid properties are changing rapidly (such as in shock fronts). The equations are solved implicitly to ensure stability in the presence of radiative energy transfer and to have the time steps controlled by the rate of change of the variables and not by the Courant time for the smallest zones.

The effects of non-equilibrium ionization, excitation, and radiative energy exchange from several atomic species (H, He, and Ca) on fluid motions and the effect of motion on the emitted radiation from these species are calculated. We model hydrogen and singly ionized calcium by six level atoms and helium with a nine level atom. For helium, we collapse terms to collective levels and include the 1s², 2s, and 2p terms in the singlet system and the 2s and 2p terms in the triplet system of neutral helium, and the 1s, 2s, and 2p terms of singly ionized helium. In addition, we include doubly ionized helium. We include in detail all the transitions between these levels. For singly ionized calcium, they are the H and K resonance lines, the infrared triplet, and the photoionization continua from the five lowest levels. We use 31–101 frequency points in each line and 4–23 frequency points in each continuum, a total of 1424 frequency points. All lines are treated assuming complete frequency redistribution (CRD) but with truncated profiles for the hydrogen Lyman lines to mimic the lower intensity in the wings in partial redistribution (PRD). For the studied calcium H resonance line, effects of PRD could potentially be important but comparisons at a number of representative times in the simulations show no difference in derived quantities between CRD and PRD; see Section 6. Continua from elements other than H, He, and Ca are treated as background continua in LTE, using the Uppsala atmospheres program (Gustafsson 1973).

Our initial atmosphere, used throughout our model parameter space, is a radiative equilibrium model with T_eff = 4000 K which is typical of sunspot atmospheres. log g is set to 4.44 (cgs units). The microturbulent-velocity parameter is set to 2.0 kms⁻¹ throughout the atmosphere, but experiments with decreasing this parameter give essentially unchanged results. A corona with a temperature of 436,000 K at 5 Mm is also added.

The upper chromosphere is modified by the incoming coronal radiation. The effect of the coronal radiation field is included in our calculations by transforming, as done in Wahlstrøm & Carlsson (1994), the incoming photon flux values given in the revised solar EUV flux model of Tobiska (1991) to the incoming intensity at the coronal boundary. The incoming intensity thus calculated is an average value over the solar disk. However, measurements from the Transition Region And Coronal Explorer (TRACE) satellite (Handy et al. 1999) of intensity in the 195 Å band over our active region at the time of our observations are five times this value. We thus scale the incoming intensity values from Wahlstrøm & Carlsson (1994) by a factor of 5 in our reference model.

The temperature, sound speed, pressure scale height, and acoustic cutoff frequency in our reference model are shown in Figure 1 for the height range of relevance for umbral flash formation. Of special relevance is the acoustic cutoff frequency that peaks at 6.5 mHz in the lower photosphere and stays at about 6 mHz for about 0.6 Mm. Waves with frequencies below 6–6.5 mHz are thus not expected to propagate into the chromosphere.

Zoom In Zoom Out Reset image size

The top boundary is transmitting but with the temperature set at 436,000 K at 5 Mm. The bottom boundary is situated where the radiation at all frequencies is optically thick, 90 km below τ₅₀₀ = 1. The bottom boundary is a non-transmitting piston with the velocity prescribed. The velocity spectrum is taken to mimic the observed photospheric velocities.

2.1. Piston Power Spectrum

Although the aim of the current work is to investigate many points in a model parameter space, we nevertheless wish to operate with piston power spectra that are representative of the photospheric wave power spectrum above solar active regions (sunspots). We have therefore constructed piston power spectra based on Doppler velocities derived from the Fe i 630.25 nm spectral line (Lites et al. 1998). This line is formed 190 km above where τ₅₀₀ = 1 (Bruls et al. 1991).

2.1.1. Power Spectrum at Fe i Formation Height

The power spectrum in the observations by Lites et al. (1998; see their Figure 6) is represented with power increasing linearly from zero at 2 mHz to a maximum value at 3 mHz, then decreasing linearly back down to zero at 4.5 mHz. From 4.5 mHz to the Nyquist frequency (500 mHz), we adopt a power law that decreases linearly on a log scale (1/e per 2 mHz). We scale the power in the two bands separately to obtain given root-mean-square (rms) velocities. As a reference run, we use an rms velocity of 41 m s⁻¹ in the frequency band of 2–4.5 mHz and 11.5 m s⁻¹ in the frequency band 4.5–7 mHz. This is very close to their values of 41 and 11 m s⁻¹, respectively (see their Table 1). Note that the power is not set to zero above 7 mHz—it is just the rms velocity measurement that is confined to a band that ends at 7 mHz. The adopted power continues above this limit with the same exponential falloff of 1/e per 2 mHz. The adopted power spectrum for the reference simulation is shown in Figure 2.

Zoom In Zoom Out Reset image size

2.1.2. Power Spectrum at Piston

The piston is located at a height where all transitions are optically thick in order to have a well-behaved lower boundary condition for the radiative transfer. This is 90 km below where the optical depth at 500 nm is unity. We parameterize our piston power spectrum by the power spectrum at the Fe i 630.25 nm formation height, but we then need the velocity 280 km further down in the atmosphere. To construct a velocity power spectrum at the piston that reproduces the chosen power spectrum at the Fe i formation height, we must take into account the changes in phase and amplitude as a function of frequency of waves propagating up from the piston through a stratified atmosphere. These effects are accounted for by calculating a velocity transfer function that relates the velocity spectra at the two heights. The desired power spectrum at the formation height of the Fe i 630.25 nm line is multiplied with this transfer function to obtain the needed piston power spectrum; see Carlsson & Stein (1997) for details.

The simulated Ca ii H-line intensity as a function of wavelength and time for our reference simulation is shown in Figure 3 together with the velocity at two different heights in the atmosphere. Many of the observed characteristics are seen in the simulated spectrum: there is an emission core that displays a sawtooth behavior with a slowly increasing redshift followed by a rapid shift to the blue with a strong brightening—an umbral flash. We see clear umbral flashes at t = 1040 s, 1190 s, and 1320 s followed by weaker activity and again two strong flashes at t = 2870 s and 3000 s. The emission intensity power spectrum has a peak between 6.5 and 7.5 mHz, corresponding to periods between 154 and 133 s, respectively. The umbral flashes have an intensity about twice that of the average intensity outside flashes with the emission blueshifted by 4.2–5.6 km s⁻¹. The Doppler shift of the emission has a striking resemblance with the actual velocity in the simulation at heights of 1–1.2 Mm and the umbral flashes coincide with when the crest of a steepened sound wave passes this height interval. Note, however, that the blueshift of the umbral flash emission is larger than the velocity at 1 Mm height.

Zoom In Zoom Out Reset image size

We proceed with studying the formation of the umbral flashes in more detail by splitting the contribution function to intensity in the same way as was done in Carlsson & Stein (1994, 1997).

The outgoing intensity for a vertical ray is

$$\begin{equation*} I_\nu =\int _{-\infty }^\infty S_\nu (z)e^{-\tau _\nu }\chi _\nu dz, \end{equation*}$$

where S_ν(z) is the source function at height z, τ_ν is the monochromatic optical depth, and χ_ν is the monochromatic opacity. The integrand is the contribution function to intensity and it can be rewritten as

$$\begin{equation*} C_I=S_\nu (z) \tau _\nu e^{-\tau _\nu }{\chi _\nu \over \tau _\nu }. \end{equation*}$$

This way of splitting the contribution function highlights several different aspects. The source function term highlights the local contribution to the radiation field, the $\tau _\nu e^{-\tau _\nu }$ term has a maximum at τ_ν = 1 emphasizing that the emergent intensity is approximately given by the source function at that depth (the Eddington–Barbier relation) while the χ_ν/τ_ν factor gives weight to those depths where we have a velocity gradient such that we have large opacity at a given frequency but a small optical depth.

In Figures 4–7, we show the full contribution function to intensity (lower right panels) together with the three factors above. Figure 4 shows the situation at the time of the umbral flash at t = 1190 s. The intensity profile (lower right panel) shows a pronounced peak blueshifted to 5.8 km s⁻¹. The contribution function to intensity has a peak at a height of 1.1 Mm (lower right panel). This peak comes from a maximum in the source function because of increased temperature in the compression phase of the wave (upper right panel) sampled by the optical depth being unity at that wavelength (lower left panel) and amplified by the χ_ν/τ_ν term (upper left panel) being large because of large opacity and small optical depth (because of the velocity gradient). The corresponding red peak where the optical depth is unity crosses the maximum of the source function is weak because of a small value of the χ_ν/τ_ν term on the red side: small opacity (large distance between this wavelength and the Doppler-shifted line center) and large optical depth (small distance between the line-center position higher up in the atmosphere and the red wavelength). The emission peak is not formed at the local line-center wavelength (which is Doppler-shifted by the local velocity and is thus situated at the position of the velocity curve in the figures) but on the blue side because this is when we can get the maximum combined effect of the three terms making up the contribution function. This explains why the emission peak has a larger blueshift than the local velocity at the height of formation of the intensity. Ten seconds later (Figure 5), the wave has steepened further and moved higher up in the atmosphere. The size of the Doppler shift is larger but the source function is decoupled from the Planck function so the emission goes down. This continues at t = 1210 s (Figure 6): the wave moves higher, gets stronger, and starts to shock. The source function at the position of the wave crest gets smaller and the position in the line profile where optical depth is unity at the increasing height of the wave gets closer and closer to local line center. The increase of the wave amplitude delays the decrease of the emission intensity and even sometimes causes the Doppler shift of the emission to increase after the flash, causing curved paths in the intensity versus time plots (Figure 3). At t = 1280 s (Figure 7), we have a maximum redshifted emission that is formed rather high in the atmosphere at z = 1.4 Mm.

Zoom In Zoom Out Reset image size

Zoom In Zoom Out Reset image size

Zoom In Zoom Out Reset image size

Zoom In Zoom Out Reset image size

The umbral flashes are thus caused by steepened acoustic waves in the same way as Ca ii H_2V bright grains in the quiet Sun. The next question we ask is whether we can trace the waves back to the photosphere.

Figure 8 shows the velocity scaled with the square root of the density (for a propagating acoustic wave with constant speed that would be proportional to the acoustic energy flux) for our reference simulation. Umbral flashes are shown at their time of occurrence and height of formation. It is obvious that the umbral flashes occur during times of the passage of steepened waves at a height of 1.0–1.3 Mm. In the photosphere, the full velocity field (left panel) is dominated by standing waves (corresponding to the peak at 3 mHz in the input piston spectrum; see Figure 2) and it is difficult to follow the umbral flash causing wave back to the photosphere. When only frequencies above the cutoff frequency of 6.5 mHz are retained (right panel), it is easy to see that the wave that causes the umbral flash is also present in the photosphere. Downwardly propagating waves (streaks going from lower right to upper left) are waves reflected off the transition region steep temperature increase.

Zoom In Zoom Out Reset image size

The umbral flashes are thus caused by propagating acoustic waves (frequency above 6.5 mHz) that steepen in the chromosphere. We now proceed to investigate how the characteristics of the umbral flashes depend on the input velocity power spectrum, the incoming coronal radiation field, and possible extra heating in the chromosphere.

We set up a grid of experiments with a view to determining how observable quantities depend on the chosen parameters in the grid. The four axes in our parameter space are piston power spectrum below and above 4.5 mHz, incoming coronal radiation, and additional mechanical heating in the chromosphere. In the following sections, we describe how we have set these parameters in our grid.

4.1. Piston Power Spectrum

4.2. Incoming Coronal Radiation

In the reference model, we have scaled the incoming radiation field used by Wahlstrøm & Carlsson (1994) by a factor of 5 to take into account the increased coronal radiation field over a sunspot compared with the average Sun. To test the influence of the incoming coronal radiation on the behavior of the Ca ii H line, we have constructed models with four different values for the incoming coronal radiation. The scaling factors used are 0, 1, 5, and 10. The initial temperature structure for these different scaling factors is shown in Figure 9.

Zoom In Zoom Out Reset image size

4.3. Heat Input

4.4. Grid Summary

We have discussed our parameter space and the methods we have employed for exploring this space in previous sections. In this section we examine, first qualitatively, the effects of varying each parameter in turn keeping all other parameters fixed to their reference values from Table 1. This will obviate the need for defining quantitative measures which will admit a statistical examination of these effects. Before we begin examining the effects of varying the physical parameters, we need to gain a measure of the changes produced by different random seeds for the piston power spectrum. The bounds of this deviation then define the lower threshold of what can be considered to represent a physical change in the simulation. We discuss these random seed effects next.

5.1. Effects of Random Seed Variation

To gain a measure of the effects that different random seeds might have on the spectra, we calculate four models, each with a different random seed, but the same set of standard values (see values in bold type in Table 1). In Figure 10, we plot the Ca ii H-line intensity as a function of Doppler shift about the center frequency and time. In Section 5.6, we describe our procedures for defining a statistical measure to be used for the comparison of the key parameters between models. We find that there are very small differences in the overall behavior displayed in the four panels of Figure 10. The periods, the intensity increases, and excursions of the line core from the center frequency are statistically very similar in the different realizations, the main difference being a shift in time. This is the expected effect of a change in the random seed (giving differences in the phase of waves/shocks).

Zoom In Zoom Out Reset image size

5.2. Effects of I⁻_ν Variation

Figure 11 shows the time evolution of the line core profile and its Doppler excursions shown in velocity units as a function of increasing values of incoming coronal radiation, I⁻_ν. The scaling of I⁻_ν is as factors of the quiet Sun values used in Wahlstrøm & Carlsson (1994, WC). These are, from left to right, 0, 1, 5, and 10 times WC. With no incoming coronal radiation field (panel (a) of Figure 11), the flashes are very weak and there is no emission between the flashes. Even with quiet Sun values for the incoming coronal radiation (scaling factors of 1, panel (b)) the emission between times of flashes is very weak. Only by increasing the incoming coronal radiation to five times WC do flashes appear with emission at all intermediate times. Panels (c) and (d) in Figure 11 show increasingly clear flashes for 5 and 10 times WC, respectively, of incoming coronal radiation. With increasing incoming coronal radiation, we get increasing flash brightness and increasing emission between flashes while the Doppler excursions do not seem to change significantly.

Zoom In Zoom Out Reset image size

5.3. Effects of Power Variation in the 2.0–4.5 mHz Range

Figure 12 shows the time evolution of the intensity (gray scale) as a function of Doppler shifts, given in velocity units, for successively increasing V_rms values in the 2.0–4.5 mHz range (see Table 1). The periodic rapid shift of the line core to blue followed by its gradual return to the red is very clear as is the onset of flashes in all panels. However, there is no discernible change in this behavior for different rms values of the velocity in this frequency range in spite of an eightfold increase in V_rms from panel (a) to panel (c), indicating that power above 4.5 mHz must be the determinant in flash formation. We discuss this in the next section.

Zoom In Zoom Out Reset image size

5.4. Effects of Power Variation in the 4.5–7.0 mHz Range

Figure 13 is a clear demonstration of the effects of varying power in the 4.5–7.0 mHz frequency range. Removing piston power in this range produces a nearly symmetric line core profile with no Doppler shifts (panel (a)) while increasing power produces stronger and stronger flashes with increasing excursions to the blue and to the red. With increased power, the transition from a redshifted emission to the umbral flash is also more sudden with the increase of emission almost over the whole profile (see especially panel (d)). This is because the shock formation height is lower when the amplitude is higher so the wave has shocked before reaching the formation height of the umbral flash.

Zoom In Zoom Out Reset image size

5.5. Effects of Volume Heating

Figure 14 shows the effects of volumetric heat input in the simulated atmosphere. Our reference model (values in bold face in Table 1) includes no heat input (panel (a)). With increasing additional volumetric heat input, the overall brightness increases (both the flash brightness and the emission at times between the flashes). The Doppler excursions do not seem to depend on the heat input.

Zoom In Zoom Out Reset image size

5.6. Statistical Measures

As we have seen above, the main characteristics that change in the simulated Ca ii H-line intensity profiles are the flash brightness, the emission at times between flashes, and the size of the Doppler excursions. In order to facilitate a quantitative comparison of the different simulations, we proceed with defining three measurable quantities: the mean flash intensity, the mean minimum emission, and the mean of the Doppler excursion of the flashes. In the following, we will first define and illustrate each of the these statistical measures for a given subset of our models. We then proceed to examine the whole of our model grid in light of these measures.

5.6.1. Flash Intensity and Mean Minimum Emission

As a measure of the intensity evolution, we take the maximum of the core part of the intensity as a function of time. The core part is defined to be the part closer to line center than the minimum intensity of the profile.

As a measure of the flash brightness, we use the average flash brightness. A flash is defined as a maximum intensity that is larger by at least 0.6% of the photospheric continuum intensity than the intensity 30 s before or after the flash.

As a measure of the emission between flashes, we take the mean of the minima in the time evolution of the maximum of the core part intensity.

Figure 15 demonstrates these ideas and shows how flashes are defined in the following. Maximum emission core intensity, as a fraction of the continuum intensity, is plotted here against the time for different values of incoming coronal radiation. Curves correspond to panels (a) (solid), (b) (dotted), (c) (dashed), and (d) (dot dashed) of Figure 11.

Zoom In Zoom Out Reset image size

5.6.2. Flash Doppler Velocity

For each of the flashes fulfilling the criteria given above, we measure the Doppler shift of the maximum intensity point. The mean of these Doppler shifts is taken as a measure of the Doppler excursions of the emission core. Inspection of Figures 10–14 suggests that the Doppler excursions may vary with the intensity of the flash and thus our measure might be strongly dependent on the exact definition of a flash. Figure 16 shows the Doppler shift versus the intensity of flashes in the simulation with V_rms = 46.0 m s⁻¹ in the frequency range 4.5–7.0 mHz. There is a clear trend of larger Doppler shift with larger flash intensity, thus confirming the suspicion that the measurement is dependent on how a flash is defined. There is also a rather large spread in both the flash intensities and Doppler shifts.

Zoom In Zoom Out Reset image size

5.6.3. Statistical Summary

Figure 17 shows how the statistical measures defined above vary across our grid. When one parameter is changed, all other parameters have the value from the reference model. When the velocity amplitude in the 2.0–4.5 mHz frequency range is varied (leftmost column in Figure 17) there is no significant variation in any of the three measures. These low-frequency oscillations are below the cutoff frequency of the atmosphere and have no effect on the Ca ii H-line emission even when the power in this frequency range is much higher at the formation height of the Fe i 630.25 nm line than in the higher frequency range.

Zoom In Zoom Out Reset image size

The power above 4.5 mHz (see second column of Figure 17) has a very strong influence on the Ca ii H-line emission behavior with increases in both the Doppler excursions (top panel) and flash intensity (middle panel). The mean minimum of the emission intensity goes slightly down with increased velocity power because increased shock strength also means lower temperature behind the shock.

The incoming coronal radiation (third column of Figure 17) also plays an important role for the Ca ii H-line emission. The effect on the Doppler excursions is small but the effect on both the flash intensity and minimum emission intensity is large.

The extra heating (rightmost column in Figure 17) has a small effect on the Doppler excursions. Both the flash intensity and the mean minimum emission intensity increase with increased heating.

The Ca ii H profiles have been calculated assuming CRD although both the inner wings and the peak intensity are expected to be lower in PRD (see, e.g., Shine et al. 1975). To estimate the effect on our results, we calculated the emergent intensity in PRD from the four snapshots shown in Figures 4–7. Both effects are present but at a level of only a few percent and the derived quantities in Figure 17 change less than the variation caused by different seeds (shown with error bars in the figure). Our conclusions are thus not affected by our assumption of CRD in the profile calculations.

We neglect magnetic effects in all our simulations. The justification is that the main effect of a strong, vertical magnetic field is to restrict the plasma motions to the vertical direction, thus lending support to the one-dimensional approximation. In the limit of a strong magnetic field (Alfvén speed much greater than the sound speed), the slow mode wave will propagate with the sound speed along the magnetic field. The acoustic waves in our hydrodynamic simulations would, in the presence of the strong magnetic field in the umbra, be more correctly called slow mode waves, but their properties would be very similar.

We have seen that umbral flashes are caused by slow mode waves that propagate to chromospheric heights and steepen into shocks. The maximum umbral flash emission is formed when the wave is between 1.0 and 1.3 Mm above where optical depth at 500 nm is unity (in the umbra). This is well below the shock formation height in the reference model. With increased velocity power in the photosphere above 6.5 mHz, the shock formation height is lower in the atmosphere and the waves have shocked before they reach the formation height of the umbral flashes. This formation height is also lowered because the strong shock produces an increase in the source function. In the simulation with V_rms = 46 m s⁻¹ in the 4.5–7.0 mHz range (Figure 13, panel (d)) the second umbral flash comes 20 s earlier than in the reference run (Figure 13, panel (b)) at t = 1170 s and the contribution function to intensity has its maximum at a height of 0.9 Mm compared with 1.1 Mm in the reference simulation.

Waves with a frequency below the cutoff frequency of the atmosphere (6.5 mHz) do not have a significant effect on our measures.

How do these statistical measures compare with observations? The observed flash emission is often relative to the quiet umbra and a factor of 2 brightening is reported. We have chosen to relate the flash brightness to the photospheric continuum instead (since the photospheric continuum is a more robust measure that is not affected by stray-light issues) and it is therefore difficult to compare. If we take the minimum emission as a measure of the quiet umbra, we get a flash brightening of a factor of 2 in our reference calculation. Varying the incoming coronal radiation affects both the flash brightness and the minimum brightness in the same way such that the brightening factor remains about 2. Increasing the power in the 4.5–7 mHz band increases the ratio from two for the reference model to five in the model with largest power. The effect of scattered light in observations would decrease the observed brightening so we conclude that the values in our simulations are similar to the observed values.

The reported blueshifts of the flashes (7 km s⁻¹) are higher than what we get in our reference run and only reproduced when we increase the power above 4.5 mHz by a factor of 4 compared with the value reported by Lites et al. (1998). The power in the 4.5–7.0 mHz band is dominated by the non-propagating part in the 4.5–6.5 mHz band so it is quite possible that a factor of 4 increase in power above 6.5 mHz is still compatible with their observations. It is also possible that their measured values are lower than the typical values for umbrae with clear flashes. A detailed comparison should be made where the photospheric velocity field and the Ca ii H-line spectra have been recorded simultaneously. Such a study is underway.