Formation of the thermal infrared continuum in solar flares

Observations of the Sun with the Atacama Large Millimeter Array have now started, and the thermal infrared will regularly be accessible from the NSF's Daniel K. Inouye Solar Telescope. Motivated by the prospect of these new data, and by recent flare observations in the mid infrared, we set out here to model and understand the source of the infrared continuum in flares, and to explore its diagnostic capability for the physical conditions in the flare atmosphere. We use the 1D radiation hydrodynamics code RADYN to calculate mid-infrared continuum emission from model atmospheres undergoing sudden deposition of energy by non-thermal electrons. We identify and characterise the main continuum thermal emission processes relevant to flare intensity enhancement in the mid- to far-infrared (2-200 $\mu$m) spectral range as free-free emission on neutrals and ions. We find that the infrared intensity evolution tracks the energy input to within a second, albeit with a lingering intensity enhancement, and provides a very direct indication of the evolution of the atmospheric ionization. The prediction of highly impulsive emission means that, on these timescales, the atmospheric hydrodynamics need not be considered in analysing the mid-IR signatures.


Introduction
Solar flares are the product of sudden energy release in compact regions in the solar atmosphere. The energy transferred from the magnetic field in active regions to the local plasma causes heating and particle acceleration, leading to the emission of electromagnetic radiation. The solar flare phenomenon was first discovered independently by Carrington (1859) and Hodgson (1859), in visible light, and since then, the broadband nature of the electromagnetic emission of flares became clear, following the first detection of flare emission in several spectral ranges, radio (Hey 1983), hard (Peterson & Winckler 1959) and soft X-rays (White 1964), γ-ray lines (Chupp et al. 1973;Talon et al. 1975). Sporadic observations in the infrared (IR) and sub-millimetric (sub-mm) wavelengths have been performed over the years, but not aimed at flares studies (see Loukitcheva et al. 2004, for a summary of these observations). Only recently, however, observations in the mid-infrared range become feasible for flare studies (Kaufmann et al. 2013;Penn et al. 2016). As we will demonstrate in this paper, flare IR continuum observations are of particular value as they give us a very direct diagnostic of atmospheric electron content, with a temporal and spatial resolution adequate to discriminate between different flare energy transport models.
The first tentative description of IR flare emission was proposed by Ohki & Hudson (1975) as a part of early unsuccessful searches for flare signatures (cf. Hudson 1975). They suggested two main possible emission mechanisms: optically thin thermal free-free radiation from the chromosphere or blackbody radiation (dominated by H − opacity) from the heated photosphere, following the ideas proposed to explain the formation of the white-light continuum in flares.
More recently, Heinzel & Avrett (2012) performed calculations of thermal IR flare emission using existing semiempirical model atmospheres of flares (Machado et al. 1980;Mauas et al. 1990). Their results strongly indicated the dominance of ion free-free emission from the chromosphere. 1 A similar method was employed by Trottet et al. (2015) to investigate the origin of the 30 THz (10 µm) emission in the M-class flare SOL2012-03-13 (Kaufmann et al. 2013), and led to similar conclusions. Penn et al. (2016) have shown that the observed IR spectrum during SOL2014-09-24 is consistent with free-free emission from a chromospheric region with an upper limit to the electron density of ≈ 4 × 10 13 cm −3 . Kašparová et al. (2009a) used non-LTE radiative hydrodynamic models to investigate the effects of non-thermal excitation and ionisation of hydrogen in the formation of sub-mm emission. Cheng et al. (2010) used radiative hydrodynamic simulations performed with RADYN (Carlsson & Stein 1995, 1997 to investigate the formation of the continuum at 1.56 µm (Xu et al. 2006). Kleint et al. (2016), using observations in the near IR, optical and near ultraviolet (NUV), have concluded that the NUV contiuum excess intensity is more consistent with H free-bound continuum, while the optical and near-infrared contiuum excess intensity could be explained by a blackbody spectrum.
In this work, we employ 1D radiative hydrodynamic simulations using RADYN to investigate the formation of the IR continuum in the dynamically evolving flare atmosphere. We study very short episodes of energy injection, to be able to isolate the most important timescales namely ionization and recombination. This important limit has previously been studied by Kašparová et al. (2009a), who showed the possibility of sub-second emission signatures in the sub-mm ranges.We find that the continuum IR signature, particularly at longer wavelengths, tracks energy injection with sub-second fidelity, and that the only parameter of importance is forming the radiation is atmospheric ionization. This means that details of the hydrodynamic development of the atmosphere are irrelevant for interpreting this measurement on short timescales, consistent with the simplifying "constant density" limit of Shmeleva & Syrovatskii (1973), in which the atmospheric structure does not change while the radiation is being formed. We present a detailed analysis of the formation height of the IR continuum, main sources of opacity, synthetic spectra and light curves. We also compare our results with IR data from Penn et al. (2016).

Calculations of the infrared continuum in solar flares
The formation of the infrared continuum is significantly simpler than the formation of recombination continua or spectral lines because the continuum source function is given by the local Planck function S ν = B ν (T ), with T being the kinetic (or electron) temperature. Nonetheless, the chromospheric opacities depend on electron, proton and neutral hydrogen densities that must be treated under fully non-LTE conditions. In this work, we used RADYN (Allred et al. 2005 simulations to obtain the dynamic flaring chromosphere, subjected to collisional heating by a 1 We refer to free-free emission mechanism on ions, predominantly H ii, as "ion free-free", and on neutrals as "neutral freefree" continua. The latter is often also called "H − free-free emission" for hydrogen. Generally free-bound contributions are relatively small in the mid-IR. non-thermal distribution of electrons, as described in Sect. 3 below. For the calculations of the IR emission we follow a similar approach to previous works (e.g. Heinzel & Avrett 2012) where a 1D model atmosphere is used to calculate the relevant opacities and the radiative transfer for a range of IR wavelengths. In this work we calculate the radiation emerging from the top of the 1D atmosphere, i.e. with a heliocentric angle cosine µ = 1. Center-to-limb effects are likely to be important; these cannot be assessed with 1D models, as one needs to account for the intrinsic anisotropy of the flare emission, and also its propagation through the quiescent atmosphere surrounding the flare site.
In the chromosphere, the main source of opacity for the infrared continuum is ion free-free continuum, while in the temperature minimum region and below the opacity is dominated by neutral free-free opacity.
The ion free-free opacity κ ν (H) (in cm −1 ) is calculated by (Rybicki & Lightman 1986) κ ν (H) = 3.7 × 10 8 T −1/2 n e n p ν −3 g ff (1) where n e and n i are the electron and proton densities respectively, and T is the kinetic (or electron) temperature. Note that n e includes contributions from hydrogen, helium and metals, as calculated by the employed model, see Section 3. Accurate values for the Gaunt factor g ff were obtained by interpolating the table of numerically calculated g ff provided by van Hoof et al. (2014). Lower in the atmosphere, around the temperature minimum region and below, H − free-free opacity dominates. We use the expression for the H − free-free opacity given by Kurucz (1970).
where n H is the neutral hydrogen density, and the numerical coefficients are A 1 = 1.3727 × 10 −25 , A 2 = 4.3748 × 10 −10 , and A 3 = 2.5993 × 10 −7 (Kurucz 1970). A similar expression given by John (1988) is in good agreement with Eq. 2, within 1-3%. Using Eq. 1 and 2, the total opacity is then with the term (1 − e −hν/k b T ) being the correction for stimulated emission, where h and k b are the Planck and Boltzmann constants.
Here we define the contribution function CF (given in erg s −1 cm −3 Hz −1 sr −1 throughout this paper), following e.g. Carlsson (1998): where j ν = κ ν B ν (T ) is the emission coefficient, and B ν (T ) is the Planck function. We define CF as a function of height h where the height is related to the optical depth scale via dτ ν = −κ ν dh, the optical depth. The contribution function CF indicates the formation height of emission. Finally, the specific intensity is calculated by integrating the CF along the line of sight: which we use to calculate the brightness temperature T b using the Rayleight-Jeans law: 3. Radiative hydrodynamic flare modelling RADYN is a well-established numerical tool for modelling solar flares. Abbett & Hawley (1999) and Allred et al. (2005) modified the code of Carlsson & Stein (1995, 1997 to simulate flares, and since then RADYN has been used by several authors to better investigate the hydrodynamic and radiative response of the solar atmosphere to the injection of flare energy (e.g. Kuridze et al. 2015;Kennedy et al. 2015;Rubio da Costa et al. 2016;Kerr et al. 2016) RADYN solves the coupled, non-local equations of hydrodynamics, radiation transport and atomic level populations (including time-dependent effects) on a 1D adaptive grid. Backwarming by soft X-ray, extreme ultraviolet and ultraviolet (XEUV) radiation is included self-consistently in the code as an additional heating term. For a recent description of RADYN, along with a list of atomic transitions solved, the reader is encouraged to consult Allred et al. (2015). In part because of the restrictions of 1D, the results presented here are numerical experiments and should not be interpreted as direct modelling of specific flares. Instead, in this work we have tried to understand the ways in which the IR emission shows the response of the chromospheric plasma activated by flare energy input. In our use of the RADYN code we explicitly assume electron energy transport, as in the "thick-target" model (Brown 1971;Hudson 1972), in which an electron beam formed in some coronal region precipitates into the pre-flare atmosphere. We further explore this model only in the restrictive case of vertical, unidirectional beam propagation (but including RADYN's Fokker-Planck treatment of electron propagation), and ignore other forms of energy transport such as the Poynting flux (Emslie & Sturrock 1982;Fletcher & Hudson 2008).
The electron beam is described as a power law in energy with spectral index δ and low-energy cutoff E c . We have run a set of RADYN simulations with spectral index δ = 5, 6, 7, and 8, for each value of the low energy cutoff E c = 20, 25, and 30 keV.
is the Heaviside function. In all 13 models the energy input followed a Gaussian time profile with a FWHM of 0.553 s and with a maximum of F = Ef (E)dE = 10 11 erg s −1 cm −2 at t = 2.0 s. Hereafter, we define our energy input function as F11t, following the typical notation in previous works e.g. F11 for F = 10 11 erg s −1 cm −2 in Allred et al. (2005). Note that this beam power per unit area roughly matches the solar luminosity per unit area. In Figure 1 we summarize the properties of the atmosphere at the time of peak energy deposition. Different electron beams will heat the atmosphere in slightly different ways (see detailed discussion in Allred et al. 2015). The lowenergy cutoff E c effectively controls the height of the energy deposition, i.e. how deep the lowest-energy electrons (which carry most of the energy) can penetrate the atmosphere. A beam with a harder spectral index (lower δ value) will spread the energy deposition deeper in the atmosphere than a softer (higher δ value).
The pre-flare atmosphere that we use for all runs, as shown in many of our figures, is a radiative equilibrium model (Abbett & Hawley 1999;Allred et al. 2015), with a maximum temperature of 1 MK at the top of the loop, which has a total length of 10 Mm. The pre-flare transitionregion pressure is about 0.22 dyne cm −2 . Its details may have only minor relevance to the flare effects. For comparison the standard VAL-C model (Vernazza et al. 1981) has a transition-region pressure of about 0.2 dyne cm −2 . This short time scale that we have adopted allows us to follow the impulsive development of ionization and mid-IR emission in the constant-density limit (Shmeleva & Syrovatskii 1973), in which the flows resulting from the energy input cannot play a major role for lack of time in which to develop. For each model, we used the output atmosphere at time intervals ∆t = 0.1s to calculate the opacities κ ν (H) and κ ν (H − ), the contribution function CF, and brightness temperature T b (Eq. 1 to 6) for wavelengths in the infrared: λ = 2, 10, 50, 200 µm (ν = 150, 30, 6 and 1.5 THz, respectively).

Formation height and optical depth
In Figure 1 we show the flaring atmospheric parameters from RADYN simulations and IR calculations at λ = 10 µm for all the models, at the time of maximum energy input (t = 2.0 s). Different beam models deposit the energy in slightly different ways (see panel Fig. 1e), causing different atmospheric responses e.g. temperature and pressure (especially above 1 Mm). The chromospheric temperatures ( Fig. 1a) rise from a few 10 3 K to a few 10 4 -10 5 K between 0.7 < h < 1.5 Mm, along with a large increase of the electron density ( Fig. 1c), from about 10 10 cm −3 to a few 10 13 cm −3 , as hydrogen becomes fully ionised ( Fig.  1f), at and above the locations where the peak heating occurs. There are no significant effects below h ≈ 0.5 Mm. The short heating timescales in these models is not sufficient to increase the temperature in the mid-chromosphere enough to result in strong backwarming by UV radiation and subsequent photospheric temperature enhancement (e.g. Machado et al. 1978).
Qualitatively, the production of mid-IR (10 µm) emission is very similar for all of our models. During the simulated flares, the IR emission has an optically thin contribution from the flaring chromosphere superposed on an optically thick contribution from the undisturbed photosphere (Section 4.1), with the flare IR emission increase due almost entirely to changes in the optically thin contribution. The contribution function CF (Figure 1b) shows an enhancement in the chromosphere, with a maximum just above 1 Mm, varying within 0.1 Mm depending on the beam model. The optical depth (Figure 1d) increases by many orders of magnitude, but remains optically thin, reaching 0.1 < τ ν < 0.01, depending on the beam model. The enhancement of τ ν , and consequently the CF is a direct result of the increase of the electron density rather than an increase in local temperature.
The main differences in CF and τ ν for the various beam models appear roughly at heights 0.6 < h < 1 Mm. These differences can be traced to the effect of the beam on the hydrogen ionisation and consequent increase of the electron density. Beam models with higher low-energy cutoff E c values and harder spectral indices δ will have more electrons with higher energies, that can penetrate and deposit en-  ergy deeper in the chromosphere (e.g. Emslie 1978), where the higher density of hydrogen provide a larger supply of electrons.
The beam model that gives the largest CF and τ ν has E c = 30 keV and δ = 5, which is able to heat deeper in the atmosphere, ionising more hydrogen than the other beam models and thus producing more free electrons. The weakest cases are E c = 20 keV, with both δ = 7 and δ = 8, meaning that, proportionally, most of the electrons collisionally stop higher in the chromosphere (about 0.1 Mm higher than the latter case), and hence produce fewer free electrons.
We now characterize the IR emission at wavelengths λ = 2, 50 and 200 µm, comparing with the results at 10 µm ( Figure 1). Figure 2 shows the contribution function CF and the optical depth τ ν for all the beam models, at the four selected wavelengths. The CF at 2 µm and 10 µm ( Figure  2a and Figure 2c) are very similar for all the models, and τ ν at 2 µm (Figure 2b) is smaller than at 10 µm (Figure 2d).
At 50 µm the opacity regime in the chromosphere starts to change, with the optical depth reaching τ ν > 1 ( Figure  2f). This affects the CF (Figure 2e): as the optical depth increases in the chromosphere past τ ν > 1, the contribution to CF from the undisturbed photosphere is strongly reduced. As the hydrogen recombines following the cessation of flare energy deposition, the opacity regime reverts to being optically thin (this happens around t = 6 s in our simulations). At 50 µm, during the energy deposition phase, the flare emission is predominantly optically thick from the chromosphere.
At 200 µm, the chromospheric emission becomes very optically thick (τ ν > 10) below ≈ 1 Mm at t = 2.0 s. In fact, as soon as H starts to ionise (around t = 1.6 s), the electron density reaches n e ≈ 10 12 cm −3 , which is sufficient to make τ ν > 1. After the energy input, the chromosphere remains highly ionised, supplying enough free electrons to maintain τ ν > 1 until the end of the simulated flares, at t = 60 s. These differences in CF and τ ν for a range of IR  wavelengths are a direct consequence of the ion free-free opacity, Eq. 1, which increases for lower frequencies (longer wavelengths).
In our simulations, the change in the opacity regime, from optically thin to thick, occurs around 50 µm. However, the wavelength at which this change happens will depend slightly on the initial atmosphere and strongly on the intensity of the energy deposition. The density profile of the atmosphere defines the column depth for penetration of the electron beam, and hence the heating. It also dictates the amount of hydrogen that can be ionised, and therefore, sets a limit to the largest electron density available.

Comparison with semi-empirical atmospheric models
Semi-empirical atmospheric models are often used to assess the formation height of the radiation under investigation or properties of the atmosphere (e.g. Xu et al. 2004;Trottet et al. 2015). Here we compare our calculations of the time-averaged (1 < t < 3 s) contribution function CF and τ ν at 10 µm from our model with δ = 5 and E c = 30 keV, with the same properties obtained from the semiempirical models for the quiet-Sun VAL-C (Vernazza et al. 1981) and C7 (Avrett & Loeser 2008), and flare models F1 and F2 (Machado et al. 1980). CF and τ ν are presented in Figure 3, along with the temperature T and electron density n e for each model. The reason for using time-averaged quantities from our calculations is simply to approach the time-and space-average nature of the semi-empirical models. Qualitatively, all models are in good agreement: all models present an optically thick (τ ν > 1) contribution formed in the photosphere (h < 0.1 Mm), with an optically thin (τ ν ≪ 1) chromospheric contribution (h ≈ 1 Mm) in the flare models. The location and intensity of the chromospheric CF of a given model is a direct consequence of its electron density profile.
The F2 model has a stronger CF in the chromosphere than in the photosphere, as is the case for F1 and our calculations. This effect comes directly from the model's larger chromospheric n e values, which in turn, comes from the larger hydrogen density between 0.5-1 Mm (Machado et al. 1980). The larger n e produces a larger optical depth, that approaches τ ν ≈ 1 in the chromosphere, essentially blocking the photospheric contribution. Note that this does not occur for the 'weaker' F1 model, and its photospheric contribution is in very good agreement with the quiet-Sun models VAL-C and C7, and also with our RADYN calculations.
Note the excellent agreement of our CF and τ ν calculations for the models C7 and F2 with the calculations by Trottet et al. (2015), who used the PAKAL radiative transfer code (De la Luz et al. 2010. We have also calculated FLA and FLB flare models from Mauas et al. (1990), and note that they do not produce a strong enhancement of the IR emission against the quiescent photospheric emission, given that both models have a lower n e in the chromosphere compared to F1 and F2.

Sources of IR opacity enhancement
The near-and mid-IR emission in the quiet Sun is dominated by the H − opacity near disk center. The physical height of the τ ν = 1 layer falls just longward of the H − absorption edge, slightly below the height of the visible photosphere (the "opacity minimum"), but increases monotonically for longer wavelengths. It reaches the height the visible photosphere again roughly at the temperature minimum region, near 200 µm.
Under flaring conditions, the flare energy deposition ionises hydrogen, helium, and the metals, which supplies free electrons that enhance the ion free-free opacity and creates the IR flare excess emission. The interplay of these two sources of opacity as a function of wavelength shown A&A proofs: manuscript no. radyn_ir_arxiv in Figure 4 evidence the dominance of ion free-free opacity (blue) in the chromosphere over the contribution from the undisturbed photospheric H − opacity (red).
Note that RADYN explicitly calculates the radiative transfer necessary to estimate the contributions from hydrogen and helium, including continua and many lines (e.g. Allred et al. 2015). Helium does not play an important role as a source of electrons, but its ionization to He ii and then He iii contributes to defining the electron temperature. Thus, hydrogen ionization generally dominates the supply of free electrons.

Synthetic IR flare spectrum
Qualitatively, we can describe the development of the flare emission as follows. Initially hydrogen ionizes at the depths the incident electrons can reach, with the height of full ionization moving downward as the energy input continues. As this happens, the chromospheric ionization contributes an increasing opacity. Initially, throughout the mid-IR the chromosphere remains optically thin and the photosphere remains visible through the chromospheric contribution. However, ultimately the electron column depth becomes sufficient to produce a larger optical depth, starting with the longest wavelengths. Figure 5 shows the evolution of the IR flare T b spectrum for the model with δ = 5 and E c = 30 keV. The shape of the T b spectrum deviates from the typical form for an isothermal source (Dulk 1985). For λ > 50µm the brightness temperature tends to reflect the actual electron temperature in the non-uniform, optically-thick upper chromosphere. Note that the spectrum remains enhanced after the energy input has finished around t = 3 s, especially at longer wavelengths. The short-wavelength branch can be identified with the Rayleigh-Jeans law expected from an optically-thin slab, T b (ν) ∝ ν 2 , more easily seen in the inset panel of Figure 5, which shows the flare emission after subtracting the pre-flare spectrum. Figure 6 shows the evolution of the brightness temperature T b and contrast excess CE = ∆I/I λ (0) for 2, 10, 50, and 200 µm. For all models, T b shows an impulsive rise, quickly decreasing after the energy deposition finishes. Neither T b nor CE are very significant at 2 µm (CE < 1%), but they increases substantially for longer wavelengths and effectively increase the observability of the process. At the shorter wavelengths (2 and 10 µm in Figure 6), where the chromospheric emission is optically thin, there is a small time difference between the peak of the energy deposition and the peak of the emission. This is because the ionisation timescale is shorter than the recombination timescale, allowing the accumulation of free electrons in the chromosphere. As the energy deposition diminishes, the ionisation rate decreases and the recombination rate starts to dominate, so that the IR emission fades.
The IR emission T b directly depends on the electron density n e in the chromosphere. Figure 7 shows the T b (10µm) from the models that produced the maximum and minimum T b (10µm) values, respectively δ = 5 and E c = 30 keV (red curve) and δ = 8 and E c = 20 keV (blue curve), compared with the electron density n e at the height corresponding to the peak of CF(10 µm). Although the thickness of the ionised chromosphere and its temperature structure evolve with time, the short-wavelength IR emission is a clear proxy for the chromospheric electron density, and thus ionisation state of hydrogen. The direct implication of this finding is that mid-IR observations can be used to track the region where energy is being deposited, without the extra complications of formation and optical depth of spectral lines or recombination continua, often associated with typical flare observations, e. g. Hα (Kuridze et al. 2015), Mg ii (Kerr et al. 2015), or white-light continuum (Kowalski et al. 2017).
At longer wavelengths (50 and 200 µm in Figure 6), the emission is optically thick, and T b saturates at the local electron temperature T . As the heating increases, hydrogen ionises higher in the chromosphere, shifting the τ ν = 1 layer upwards, where the electron temperature is progressively higher, hence T b increases. When the energy deposition finishes, hydrogen recombines, the τ ν = 1 layer moves back to lower altitudes, and T b decreases. Note that at these longer wavelengths, the decay time scale following the end of energy deposition is longer than for the shorter wavelengths.
The impulsive IR emission is directly linked to the evolution of the H ionisation which, in our simulations, is determined by the heating function used. Penn et al. (2016) observations suggested impulsive heating timescales shorter than ∼4 s, with an exponential cooling with a e-folding time of 10 s. It is likely that the IR sources are not resolved in the  spectrum for the model with δ = 5 and Ec = 30 keV. For λ > 50µm the brightness temperature tends to reflect the actual electron temperature in the optically-thick upper chromosphere. The inset planel shows the flare excess emission, after subtracting the pre-flare spectrum. This highligths the positive slope of the short-wavelength branch, T b (ν) ∝ ν 2 , that can be identified with the Rayleigh-Jeans law expected from an optically-thin slab.
observations, and that multi-thread sources on fine scales must be considered (Hori et al. 1998;. The presence of multiple, unresolved, fine-scale events might help to explain the extended temporal emission that is observed by Penn et al. (2016).

Comparison with observations
In a recent observation Penn et al. (2016) obtained the first two-frequency spectrophotometry in the mid-IR, in broad wavelength bands centered at about 5.2 and 8.2 µm. We will consider here how well our impulsive modelling provides a qualitative match to these observations. Figure 8 summarizes these observations, obtained at the NSF's McMath-Pierce (McMP) Solar Facility telescope with quantum well infrared photodetector (QWIP) imaging arrays, capable of sub-second time resolution. These unique observations have clearly defined properties, i.e. double "footpoint" source, close to the white-light continuum and hard X-ray sources, and rapid variability. Fig. 8 shows maps of the contrast excess at the peak of IR emission (17:49:21 UT) during the event SOL2014-09-24 (C7.0; S10E31). These contrast maps were obtained by subtracting a time average of 6 images before the flare from the image at the time of interest. Two flare sources are clearly visible with highly significant T b and contrast excesses CE of a few percent, greater at the longer wavelength. Light curves of a single pixel in the northern source are also shown in Fig. 8, along with brightness temperature T b values for the flare excess for both wavelengths.
At 8.2 µm, at the time of the maximum emission the excess flare brightness T b is observed to be ≈ 600 K, for a contrast of ≈ 13%. At 5.2 µm the signal is weaker relative to the background, at T b ≈ 270 K, and contrast ≈ 5.7%. These imaging observations have well-defined absolute calibration from the preflare observations of the quiet Sun (e.g. Allen 1973). As noted by Penn et al. (2016), these values roughly match the expectation from a Rayleigh-Jeans spec-Article number, page 7 of 11 A&A proofs: manuscript no. radyn_ir_arxiv  trum, consistent with a warm, optically thin, chromospheric slab source.
RHESSI detected hard X-rays, and HMI observed white-light emission from this event. From the hard X-ray observations and assuming a standard thick-target for the spectrum fitting, Penn et al. (2016) report a spectral index of δ = 3.8, and estimate an energy deposition rate of 1.1 × 10 11 erg s −1 cm −2 , for electrons above E c = 50 keV. These values motivated our choices for the RADYN models presented in this paper, and calculations of T b and contrast for the Penn et al. (2016) wavebands 5.2 and 8.2 µm are shown in Figure 9. These capture the main features of the observation, specifically the contrasts, the spectral ratio, and the time scales (considering that the IR sources are not spatially resolved in the observations, see Sect. 4.4). In particular the Rayleigh-Jeans-like spectral ratio appears as a natural consequence of the small chromospheric optical depth at these wavelengths. This straightforward comparison indicates that our simulations capture the main mechanisms involved in the production of the IR emission in solar flares. This strongly suggests that the IR emission is due to thermal ion free-free emission generated by the ion free-free mechanism, as proposed by Ohki & Hudson (1975), Kašparová et al. (2009a), andTrottet et al. (2015), and tentatively confirmed by the Penn et al. (2016) observations. The necessary content of free electrons are mainly supplied by the ionisation of hydrogen in the chromosphere. For these wavelengths, the observations are consistent with optically thin ion free-free emission, as described by our simulations.
These consistencies are heartening and confirm our interest in the constant-density limit, and in the parameter space we have sampled. Based upon these models, we can anticipate much more diagnostic power from mid-IR to submm ranges (especially for sources nearer the limb, where greater heights in the chromosphere can contribute). Although the duration of our simulated flares is shorter than the overall duration of the observed event, the general properties are also in excellent agreement: a fast, impulsive rise of the mid-IR emission and also a quick, but slightly slower, decrease. The second, slower decrease of the emission we have discovered in these models, due most probably to hydrogen ionisation left over in the upper chromosphere, could not readily be isolated in the Penn et al. observations.

Discussion
Our modeling has assumed that the energy transport is entirely by electron beams, which may only be part of the story (e.g. Fletcher & Hudson 2008). It is also possible to heat the chromosphere during flares via the dissipation of Alfvenic waves, mainly via ion-neutral friction (Emslie & Sturrock 1982;Reep & Russell 2016;Kerr et al. 2016   electrons present in the chromosphere though whether in a the form of a beam from the corona or not is unknown, and of course, both mechanisms may contribute. Therefore the presence of free electrons during solar flares to enhance the IR opacity in the chromosphere is not particularly modeldependent. The non-thermal electrons play an important role in enhancing the excitation and ionisation of H in the chromosphere. In general, non-thermal collisional rates dominate the H excitation and ionisation from level n = 1 (see e.g. Karlický et al. 2004;Kašparová et al. 2009b), but the H ionisation from level n = 2 and above are still controlled by thermal rates. Since the presence of non-thermal electrons enhances H ionisation and excitation in the chromosphere, and the IR flare emission is essentially dependent on high electron density values (n e > 1 × 10 12 cm −3 ), one would expect HXR signatures from non-thermal electrons and the IR emission to be very well associated, as it was in fact observed by Penn et al. (2016).
The height distribution of the energy deposition does matter, in that a large IR enhancement requires increased ionization in a dense atmospheric layer but the short mid-IR wavelengths (2-10µm) being optically thin, track the chromospheric 'total electron content' and do not provide information on the chromospheric electron density distribution.
The departures from the Rayleigh-Jeans law seen in Figure 5 suggest that future flare observations with the Atacama Large Millimeter/sub-millimeter Array (ALMA, Wedemeyer et al. 2016) will help to find the temperature structure of the flaring chromosphere, providing valuable information to identify the energy transport mechanisms. Note, however, that at mm-waves the contribution from synchrotron emission from non-thermal electrons may become important or even dominant (e.g. Trottet et al. 2002;Lüthi et al. 2004b,a;Trottet et al. 2011;Giménez de Castro et al. 2013;Tsap et al. 2016;Fernandes et al. 2017).
We have shown that the mid-IR (2-10 µm) continuum intensity can vary on time scales close to those of the energy input, in agreement with the findings of Kašparová et al. (2009a) for the sub-mm/mm range. However, it persists on a longer time scale after the end of the energy input, and we attribute this to the accumulation of free electrons in Article number, page 9 of 11 A&A proofs: manuscript no. radyn_ir_arxiv  the chromosphere due to the non-equilibrium ionization.. The mid-IR continuum intensity is then a good proxy for the evolution of (primarily) hydrogen ionization, at around 1 Mm in these electron-beam-driven models, and also for the chromospheric electron content.

Conclusions
Using radiative hydrodynamic simulations, we have confirmed that the enhanced mid-infrared continuum in solar flares forms mainly in the chromosphere, via the ion freefree mechanism. We confirm with these models the rapid ionization and recombination of hydrogen needed to allow the mid-IR to exhibit time scales close to those of an impulsive energy input (the Shmeleva-Syrovatskii "constant density" approximation). The rapid ionisation of neutral hydrogen due to the energy input supplies the free electrons for the enhancement. Ionisation of He i to He ii can also contribute free electrons in the upper chromosphere. The brightness temperature evolution of the mid-IR continuum emission is a good proxy for the chromospheric 'total electron content', but its peak slightly lags the peak of the energy input. At higher altitudes in the chromosphere, seen in the simulations in the continuum longward of 200 µ, we see a substantially different decay profile.
The IR continuum emission is optically thin from the shorter wavelengths (2 µm) to about 50 µm, formed in the chromosphere. The τ ν = 1 layer throughout the mid-IR remains in the undisturbed photosphere for our parameter choices. More detailed spectral observations in this (and longer) wavelengths could define the development of electron content with better detail, possibly sufficient to clarify the energy-transport mechanism. In the meantime a simple uniform slab model can be used to assess the electron density where the IR continuum is formed: T b = τ ν T = κ f f ν LT ∝ n 2 e LT −1/2 , for a slab with thickness L.
On the other hand, for λ > 50µm, the emission becomes optically thick for our parameter ranges. The τ ν = 1 layer moves to the upper chromosphere, and observations should give a direct measurement of the kinetic temperature of the electrons, since T b = T for τ ν >> 1. This temperature is highly dependent upon model parameters, specifically for the treatment of energy transport.
Observationally, we expect great progress from continuum observations via continuing mid-IR observations, from observations in the near-IR with new large-aperture telescopes such as the 4-m DKIST, and from sub-mm/mm observations with ALMA.