The Astrophysical Journal, 862:59 (19pp), 2018 July 20 https://doi.org/10.

3847/1538-4357/aacc29
© 2018. The American Astronomical Society.

Modeling of the Hydrogen Lyman Lines in Solar Flares

Stephen A. Brown1, Lyndsay Fletcher1, Graham S. Kerr2 , Nicolas Labrosse1 , Adam F. Kowalski3 , and
Jaime De La Cruz Rodríguez4
1
SUPA School of Physics & Astronomy, University of Glasgow, Glasgow, G12 8QQ, UK; s.brown.6@research.gla.ac.uk
2
NASA Goddard Space Flight Center, Heliophysics Sciences Division, Code 671, 8800 Greenbelt Rd., Greenbelt, MD 20771, USA
3
Astrophysical & Planetary Sciences, University of Colorado, Boulder, CO 80309, USA
4
Department of Astronomy, Stockholm University, AlbaNova University Center, Stockholm, SE-106 91, Sweden
Received 2018 February 21; revised 2018 May 28; accepted 2018 June 10; published 2018 July 23

Abstract
The hydrogen Lyman lines (91.2 nm < λ<121.6 nm) are significant contributors to the radiative losses of the
solar chromosphere, and they are enhanced during flares. We have shown previously that the Lyman lines observed
by the Extreme Ultraviolet Variability instrument onboard the Solar Dynamics Observatory exhibit Doppler
motions equivalent to speeds on the order of 30 km s−1. However, contrary to expectations, both redshifts and
blueshifts were present and no dominant flow direction was observed. To understand the formation of the Lyman
lines, particularly their Doppler motions, we have used the radiative hydrodynamic code, RADYN, along with the
radiative transfer code, RH, to simulate the evolution of the flaring chromosphere and the response of the Lyman
lines during solar flares. We find that upflows in the simulated atmospheres lead to blueshifts in the line cores,
which exhibit central reversals. We then model the effects of the instrument on the profiles, using the Extreme
Ultraviolet Variability Experiment (EVE) instrumentʼs properties. What may be interpreted as downflows
(redshifted emission) in the lines, after they have been convolved with the instrumental line profile, may not
necessarily correspond to actual downflows. Dynamic features in the atmosphere can introduce complex features in
the line profiles that will not be detected by instruments with the spectral resolution of EVE, but which leave more
of a signature at the resolution of the Spectral Investigation of the Coronal Environment instrument onboard the
Solar Orbiter.
Key words: line: formation – radiation: dynamics – radiative transfer – Sun: chromosphere – Sun: flares – Sun: UV
radiation

1. Introduction flows are typically believed to constitute a high-velocity upflow

(“chromospheric evaporation”) detectable in high-temperature
The solar chromosphere emits increased levels of radiation
species, such as Fe XIX, with an accompanying low-velocity
during flares, which can be exceptionally energetic (∼1032 erg)
downflow (“chromospheric condensation”) in cooler species
events, initiated by a reconfiguration of the solar magnetic field.
such as He II and O V (Fisher 1989; Milligan et al. 2006;
The magnetic reconnection that triggers flares typically
Taroyan & Bradshaw 2014). Observations of these flows are
originates in the Sunʼs tenuous corona, but results in energy
important in testing flare models, as the speed, direction, and
being transported through the atmosphere toward deeper layers
duration of these flows are tied to the transport and deposition
(Benz 2008; Fletcher et al. 2011). The prevailing theory of flare
of flare energy (Fisher et al. 1985; Allred et al. 2005).
energy propagation is that electrons accelerate along the
The solar EUV output is monitored by the Extreme
reconnected field lines and deposit their energy in the
Ultraviolet Variability Experiment (EVE) instrument onboard
chromosphere via Coulomb collisions (Korchak 1967;
the Solar Dynamics Observatory. EVE consists of two Multiple
Brown 1971). Additional energy transport and dissipation
EUV Grating Spectrographs (MEGS); MEGS-A observes the
mechanisms, such as by high-frequency Alfvén waves, have
wavelength range from 5–37 nm with an FWHM of around
also been proposed (Emslie & Sturrock 1982; Fletcher &
0.1 nm. It has an additional pinhole camera (MEGS-SAM)
Hudson 2008; Kerr et al. 2016; Reep & Russell 2016).
The majority of flare emission is generated in the lower capable of measuring the region from 0.1–5 nm with a
atmosphere and emitted in the visible and ultraviolet (Fletcher resolution of around 1 nm. MEGS-B observes the region from
et al. 2011; Kretzschmar 2011; Milligan et al. 2014). The 35–105 nm, with an additional photodiode (MEGS-P) operat-
extreme ultraviolet (EUV) region is particularly interesting ing at 121.6 nm. As with MEGS-A, the MEGS-B detector has a
during these events due to its variability, which can be as high resolution of roughly 0.1 nm and a wavelength sampling of
as several orders of magnitude. This variability, observable in 0.02 nm. A MEGS-B flare spectrum is shown in Figure 1.
the Ly-α line, directly affects the Earthʼs atmosphere, resulting Unfortunately, MEGS-A is no longer operational as of May
in detrimental effects on satellites and communication systems 2014, and MEGS-B functions on a reduced duty cycle of 3 hr
(Woods et al. 2012; Kretzschmar et al. 2013). It has also been per day. (Woods et al. 2012; Brown et al. 2016)
established that flares drive flows in the chromosphere. These Emission line spectroscopy of the hydrogen Lyman lines,
which are well-observed in the MEGS-B spectrum (91.2 nm
<λ<121.6 nm), was performed by Brown et al. (2016). In a
Original content from this work may be used under the terms
of the Creative Commons Attribution 3.0 licence. Any further sample of six M and X class flares, they found red- and
distribution of this work must maintain attribution to the author(s) and the title blueshifted Lyman lines, suggesting plasma flows on the order
of the work, journal citation and DOI. of several tens of km s−1, but with some flares showing upflows

1
The Astrophysical Journal, 862:59 (19pp), 2018 July 20 Brown et al.

Figure 1. A MEGS-B spectrum of the hydrogen Lyman lines (after preflare subtraction), shortly after the onset of the X5.4 SOL2012-03-07T00:07 flare. This
spectrum was obtained from Version 6 of the Level 2 EVE data. The C III and O VI lines are also prominent in this wavelength region.

and others showing downflows. It could be that upflow 2.1. RADYN

signatures originate from eruptive features, but the observed
The RADYN code is a powerful and versatile tool for
Doppler speeds are much lower than the observed projected
probing an input atmosphereʼs response to the injection of
speeds of accompanying ejecta. energy. The code was developed by Carlsson & Stein
In this paper, we explore the potential causes of red- and (1992, 1997) to study acoustic waves in the chromosphere; it
blueshifted emission in the Lyman lines. We first use was subsequently extended by Abbett & Hawley (1999) to
simulations computed by the RADYN code (Carlsson & calculate the response to heating by a beam of non-thermal
Stein 1992, 1997; Allred et al. 2015) to simulate energy electrons (thus simulating a flare). Additional modifications to
deposition into a model atmosphere via the injection of electron the code include a Fokker–Planck beam description, which
beams of varying properties. Two values for the beam spectral more accurately models the diffusion of beam particles by
index (δ) are used, allowing us to observe the response of a pitch-angle scattering, and improved soft X-ray and EUV
“hard” beam (δ=3) and a “soft” beam (δ=8). We then backwarming (Allred et al. 2005, 2015).
extract atmospheric snapshots from these simulations and use RADYN solves the nonlinear, nonlocal equations of
the RH code (Uitenbroek 2001) to calculate model line profiles, radiation hydrodynamics, which couple the hydrodynamic
in which the important effects of partial redistribution (PRD) in equations to the non-LTE radiative transfer equation and the
the Lyman lines can be included in the radiative transfer. Both nonequilibrium time-dependent atomic-level population
codes provide detailed predictions for the Lyman line profiles, equations. The non-LTE formalism is important in cases where
which we then degrade by convolving with an approximation the radiative rates contribute significantly to the level
of the EVE instrumental profile. The comparison between these population density. Obtaining the correct level populations,
synthetic profiles and those observed by EVE reveals that the and thus the correct emission and absorption coefficients, for
loss of detailed features due to instrumental convolution can the low-density chromosphere requires a non-LTE treatment.
lead to apparent Doppler shifts that mask the true flow direction The plane-parallel equations of radiative hydrodynamics are
of the Lyman lines. solved simultaneously, accounting for the conservation of
mass, momentum, charge, and internal energy density, along
with the level population equation and the radiative transfer
2. Numerical Tools equation. These coupled equations are solved on a spatially
adaptive grid that dynamically adjusts to resolve strong
The RADYN and RH codes have been employed in many gradients and shocks (Dorfi & Drury 1987). It should be noted
other studies; they are state-of-the-art resources to obtain model that we used a modified version of RADYN that has 300 grid
predictions for the shapes and intensities of spectral lines (and cells rather than the typical 191.
continua) emitted from solar and stellar flare atmospheres. RADYN solves the level populations in our models for three
Recently, Kuridze et al. (2015) probed the origin of the elemental species that are of paramount importance in
asymmetries in Hα lines during flares. Synthetic Hα and Ca II describing the radiation field in the chromosphere and
8542 Å lines from RADYN were compared to IBIS observa- transition region. These consist of a six-level plus continuum
tions of an M class flare by Rubio da Costa et al. (2015). IBIS hydrogen atom, a nine-level plus continuum helium atom, and
observations were again compared to synthetic Na I D1 profiles a six-level plus continuum Ca II ion. Inclusion of the continua
by Kuridze et al. (2016). Simões et al. (2016) studied the allows for the calculation of bound–free transitions in addition
formation of the He lines, and Kerr et al. (2016) investigated to the bound–bound transitions. Transitions are computed for
the formation of the Mg II h and k and Ca II 8542 Å lines. five different viewing angles with up to 201 frequency points;
RADYN was also used to compare the ratio of Hα to Hβ line this is carried out under the assumption of complete
intensities in two simulations versus IBIS observations of a C redistribution (CRD), which posits that a photon undergoing
class flare, as reported by Capparelli et al. (2017). Modeling of a scattering or absorption is re-emitted with a wavelength that
the spectra of M-dwarf flares has also been performed (Allred is uncorrelated to its original wavelength. In RADYN, we
et al. 2006; Kowalski et al. 2015, 2017). Here, we outline some mitigate the effects of CRD by modeling the Lyman lines as
salient details of each code, but see Allred et al. (2015) and Doppler profiles (Leenaarts et al. 2012), but a more accurate
Uitenbroek (2001) for full descriptions. approach would be to use PRD.

2
The Astrophysical Journal, 862:59 (19pp), 2018 July 20 Brown et al.

2.2. RH beam contains more high-energy electrons, which can penetrate

and heat deeper into the chromosphere.
The radiative transfer code RH, developed by Uitenbroek
We focus on four particular simulations from RADYN: two
(2001), allows the computation of spectral line profiles with the
that deposit energy rather deep into the chromosphere with
inclusion of PRD effects. Based on the accelerated lambda
moderate and high beam fluxes, and two that deposit a greater
iteration (ALI) method for multilevel atoms (Rybicki &
fraction of their energy at a higher altitude in the atmosphere. All
Hummer 1991), it solves the equations of statistical equilibrium
simulations were obtained from the online grid of RADYN
and radiative transfer.
models, all of which are freely available, courtesy of the European
While the treatment of chromospheric lines with the
Commission-funded F-CHROMA collaboration (https://star.pst.
assumption of CRD is computationally advantageous for
qub.ac.uk/wiki/doku.php/public/solarmodels/start). The proper-
RADYN, it has been demonstrated this formalism is not
ties of each simulation are described in Table 1. Each electron
accurate for the Lyman lines (Vernazza et al. 1973; Hubeny &
beam was injected into a loop of half-length 10 Mm.
Mihalas 2014). As mentioned above, CRD assumes a lack of
The initial pre-flare atmosphere in each of the simulations
coherence between the absorbed and emitted photon (due to
was a VAL3C-like atmosphere (Vernazza et al. 1981). Our
collisions). However, in low-density environments—such as
starting atmosphere was constructed from the VAL3C temp-
parts of the solar chromosphere—there may be insufficient
erature structure, from which the heating required to sustain
collisions before a photon is re-emitted. In these cases, the
that temperature structure was computed. This atmosphere plus
coherence is not destroyed, and the wavelength of the emitted
heating function was then allowed to relax to an equilibrium
photon is correlated with that of the absorbed photon.
state, so the upper chromosphere differs somewhat from
Assuming CRD may lead to overestimation of the line wing
VAL3C (M. Carlsson 2018, private communication).
intensity, as a photon absorbed in the line core can be re-
Next, we explore line formation by following the approach
emitted in the wing—whereas in PRD, a photon absorbed in
of Carlsson & Stein (1997), examining the contribution
the core will more likely be re-emitted with a wavelength close
function to the emergent intensity (CI), which is given by the
to that of the core. CRD may also lead to an inaccurate number
integrand in Equation (1). Integrating CI over height in the
of wing photons being scattered into the core; by contrast, in
atmosphere results in the emergent intensity (Iν), meaning that
PRD, they will more readily escape because the wing is more
CI effectively tells us the locations that contribute to line
optically thin.
formation in the atmosphere. We use CI to describe features
Vernazza et al. (1973) could not replicate observations of the
observed in each simulation at key times further on in this
quiet-Sun Ly-α line without computing it by using a mixture of
paper.
CRD and PRD. Rubio da Costa et al. (2015) also note that the
application of PRD is desirable, not only for correct z1 cn
computation of the Lyman lines, but also for Hα, because this In = òz 0
Sn tn e-tn
tn
dz. (1 )
line shares an upper level with Ly-β. Similarly, Uitenbroek
(2002) found that the CRD approximation led to overestima- Due to the assumption of CRD in RADYN, Sν is constant
tions of the radiative rates in the solar Ca II K line. For these across a line profile, but does vary as a function of height. In
reasons, we use RH to obtain Lyman line profiles under the RH, we employ PRD, so Sν is a function of frequency. The
assumption of PRD. tn e-tn component, which is large when τ=1, describes how
Because RH is a time-independent code, we input atmo- radiation is attenuated as a function of height for a given
spheric snapshots obtained from RADYN simulations to the 1D c
frequency. The final component, tn , is large when there are
version of RH, building a dynamic picture from these individual n
many emitting particles at a low optical depth (where χν is the
times. We note that, in this approach, we are effectively
monochromatic opacity per unit volume), which highlights
neglecting the “history” of the atmosphere because RH re-solves
flows in the atmosphere (Carlsson & Stein 1997).
the hydrogen populations in statistical equilibrium. We input the
The RADYN outputs are then prepared for input to RH. This
non-equilibrium electron density, as computed from RADYN,
is done by decomposing the time-resolved atmospheric arrays
which mitigates this to some extent. Additionally, non-thermal
into multiple atmospheric “snapshots,” which span the duration
collisions with the beam are neglected. The non-LTE calcula-
of each simulation. These snapshots define the temperature (T),
tions are performed for any atomic species included as “active,”
electron density (ne), macroscopic velocity (Vz), and micro-
and an LTE assumption is used for species included as
turbulent parameter (Vturb=2 km s−1) on a column mass depth
background (“passive”).
scale. Note that we input the non-equilibrium electron density
into RH, which somewhat mitigates using statistical equili-
brium to obtain the hydrogen level populations. RH is then run
3. Description of Simulations and Methods
with each of these atmospheres in sequence, using a six-level
Simulations from the RADYN code were used to emulate a hydrogen atom with a continuum level. Each of the Lyman
variety of flare types, which are differentiated by the lines are treated with PRD effects.
parameters of the electron beam. The beam is characterized The final step in analyzing the output line profiles from
by its low-energy cutoff (Ec), spectral index (δ), and non- RADYN and RH is to simulate instrumental effects and
thermal energy flux. The flux varies in time, with a 20 s measure the Doppler shifts that an instrument like EVE would
triangular profile, peaking at t=10 s in all simulations. The δ observe (Brown et al. 2016). This is done by first rebinning the
and Ec values effectively determine the altitude of the energy non-constant RH data (between 30 and 130 nm) to a constant
deposition. A high δ means that the number of beam electrons wavelength spacing of 0.005 nm, before convolving the
drops off sharply as a function of energy. Therefore, the RADYN and RH profiles with the instrumental profile of
distribution is weighted to the lower-energy electrons, which EVE, given by a Gaussian of FWHM 0.085 nm (see Crotser
are stopped higher in the atmosphere. Conversely, a low δ et al. 2007), and then resampling the data to replicate EVEʼs

3
The Astrophysical Journal, 862:59 (19pp), 2018 July 20 Brown et al.

Figure 2. The atmospheric variables for the F10D3 simulation, with the pre-flare atmosphere plotted in black. The evolution of the atmosphere is represented by the
varying line colors. Negative velocities correspond to upflows. Quantities are plotted at 1.5 s intervals.

wavelength bin size of 0.02 nm. Measurements of the line shifted to the blue and the emergent profiles exhibit a
centroid variations are performed by both Gaussian fitting and blueshifted core.
intensity-weighted means. This is done for both the RADYN The line profiles appear similar to each other, with both Ly-α
output and the RH output, providing two sets of velocity results and Ly-β exhibiting central reversals. The central reversals are
for each simulation. a consequence of the line source functions having local
maxima closer to the wing formation height (around
4. Results z=1.5–1.8 Mm) than the core formation height, resulting in
stronger emission at the wing frequencies than those of the line
4.1. The F10D3 Simulation core, which forms at an altitude where Sν is relatively weaker.
The F10D3 simulation (δ=3, Ec=25 keV) describes a The blueshift in the τν=1 surface causes the central
moderate amount of heating spanning the first 20 s, followed by reversal to also be shifted to the blue. A distinct asymmetry in
30 s of cooling and relaxation. The evolution of the atmosphere the lineʼs central reversal can be seen in both Ly-α and Ly-β,
is shown in Figure 2, which shows that the atmospheric which effectively acts to reduce the amount of emission in the
temperature is quick to respond to the beam injection and blue wing relative to the red wing.
increases at all heights around the chromosphere and transition By t=26 s, the transition region has moved upward through
region (which rises to a higher altitude). The atmosphere the atmosphere to 2.6Mm, and the τν=1 surfaces for both
expands for the duration of the flare, with the transition region Ly-α and Ly-β show that the line cores are now formed at the
settling at an altitude >3 Mm by t=50 s. The chromosphere top of the chromosphere, with the line formation region now
now extends over a greater height range than does the pre-flare spanning a greater range in height. The atmosphere is still
atmosphere. In this case, ne remains elevated by almost two upflowing at the core formation height, so the τν=1 surface
orders of magnitude between z=2–3 Mm above the pre-flare consequently maintains a blue asymmetry while the line core is
value, which results in significant cooling via the high amount of still blueshifted. The line profiles are weaker in intensity
radiative losses. because Sν has decreased in the region encompassing the line
The temperature increase at early times is accompanied by formation. They still retain a central reversal, due to Sν in the
an atmospheric upflow that attains a speed of 80 km s−1 at line formation region being weaker at the core formation height
z=3Mm. The beam injection also causes a net increase in the than at lower altitudes. As a result of the asymmetric τν=1
overall electron density through a combination of the rise in surface, the central reversals are still strongly shifted to
temperature and a significant number of non-thermal collisions the blue.
with the beam itself. To understand how the detailed line profiles from the
We present the line contribution functions from RADYN for F10D3 simulation would be recorded by the EVE instrument,
Ly-α and Ly-β in Figure 3. At t=20 s, the electron beam has we perform the rebinning and Gaussian convolution
just stopped heating the atmosphere. At this time, panels (a) procedures described in Section 3. This allows us to emulate
and (c) of Figure 3 show that both Ly-α and Ly-β have the EVE instrumental profile (Crotser et al. 2007) and
developed a blue asymmetry in their τν=1 surfaces, peaking compare the simulated velocities—as they would be observed
blueward of the rest wavelength at a height around 2.3Mm. by EVE—to those reported in Brown et al. (2016). Examples
The lower right sections of these panels indicate that the line of degraded line profiles for Ly-α are shown in Figure 4; they
contribution function forms around the τν=1 height, clearly display a stark reduction in detail, compared to the
confirming that the lines are optically thick during the flare. raw data.
Here, we define the core of the line as that part of the line that
forms the highest in the atmosphere (similar to Rathore &
4.1.1. Velocities from RADYN
Carlsson 2015). The line wings form lower in the atmosphere
because wing photons are not as readily absorbed and can Velocities are obtained by measuring the deviation of the
escape more easily. line centroid positions from their rest wavelengths. These
Because the line core forms in a region of upflowing plasma measurements are achieved by two methods; the first is fitting a
(with velocities of 40–50 km s−1), the opacity structure is also four-parameter Gaussian to the line profile and obtaining the

4
The Astrophysical Journal, 862:59 (19pp), 2018 July 20 Brown et al.

Figure 3. Line contribution functions for Ly-α (a and b) and Ly-β (c and d) at two times during the F10D3 simulation that correspond to the end of the energy
deposition (t=20 s) and atmospheric relaxation (t=26.0 s). Each subfigure images a separate constituent of the overall contribution function, which is stated in the
lower right section of each panel. Darker colors indicate higher values of their respective quantities. The dashed red and green lines show the atmospheric velocity and
τν=1 surface, respectively. The dot-dashed blue and yellow lines respectively map the Planck and source functions as a function of height, in units of radiation
temperature. The emergent line intensities are overplotted on the lower right panels in red and purple.

measured line centroid; the second is by averaging the on the time-integrated line profiles (for 10 s) in order to
wavelength values weighted by the intensity at each wave- facilitate a comparison to EVE observations (Brown et al.
length bin. These methods are visualized in Figure 5. Velocity 2016). The time-integrated velocities demonstrate an additional
profiles for Ly-α through Ly-δ are plotted as a function of time loss of information due to the instrumentation when compared
in Figure 6(a). We also include a panel of velocities calculated to the results with high temporal resolution.

5
The Astrophysical Journal, 862:59 (19pp), 2018 July 20 Brown et al.

Figure 4. The Ly-α line at various times during the F10D3 simulation before (in black) and after (in green) instrumental convolution. This involves smoothing of the
line by a Gaussian of width 0.85 Å before rebinning to a wavelength sampling of 0.2 Å.

Figure 5. Gaussian fits (in red) to the degraded Ly-α line profiles (in black), as in Figure 4, at various times throughout the F10D3 simulation. Line centroid positions
derived from the Gaussian fit are indicated by the dashed red lines, and those derived from intensity weighting by the broken green lines.

Figure 6. Doppler velocities of the F10D3 line profiles obtained by emulating the instrumental effects and simulating observations from the EVE instrument.
Velocities are displayed for the RADYN (a) and RH (b) profiles, with positive velocities indicating downflows. Circular data points were obtained via the intensity-
weighting method, whereas those plotted with a diamond symbol were obtained via Gaussian fitting. The lower panels also show Doppler velocities, but with the line
irradiances time-integrated for 10 s before Gaussian fitting is performed, in order to fully emulate an EVE observation.

Each of the lines plotted in Figure 6(a) exhibits redshifts of present, although this could have been due to the variability in
varying magnitude at t=20 s, with Ly-α and Ly-β displaying the measured velocities of the higher-order EVE lines also being
the most prominent shifts during the beam deposition. Once the greater due to their weaker irradiances.
beam stops heating the atmosphere (after t=20 s), the redshift It is illuminating that the velocity signatures from the
signatures quickly peak, suggesting downflows of 30–35 km s−1 degraded profiles suggest downflows throughout the duration
in the Ly-α and Ly-β lines when the Gaussian fitting method is of this simulation. From Figures 3(a) and (c), it is clear that the
used. Ly-γ and Ly-δ exhibit redshifts corresponding to flows of beam-heating stage results in blueshifts in the line cores.
20–25 km s−1. There is a general ordering to the derived However, because the line profiles are centrally reversed and
velocities throughout the simulation, with lower-order lines shifted to the blue, the red peak becomes dominant. The red
suggesting higher speeds. This is particularly interesting because peak drags the derived line centroids further redward than the
the opposite effect was observed by Brown et al. (2016) in line core. Similarly, during relaxation, it can be seen from
velocity profiles with no preflare-subtraction when ordering was Figures 3(b) and (d) that the persistence of the blueshift in the

6
The Astrophysical Journal, 862:59 (19pp), 2018 July 20 Brown et al.

line core maintains an overall red asymmetry in each of the profiles after the beam heating stops. Further, the CRD and
lines, meaning that the velocity profiles continue to exhibit PRD results computed by RH were surprisingly similar during
redshifts until the end of the simulation. the main phase of the flare, suggesting an alternate cause for the
It should also be noted that the Doppler shifts observed difference between synthetic spectra from RH and RADYN.
between t=0–5 s are an artifact of the instrumentally RH computes the atomic level populations using the
convolved line profiles transitioning from absorption profiles equations of statistical equilibrium, whereas RADYN employs
to emission profiles, which skews the Gaussian fit as the line non-equilibrium excitation and ionization, both of which will
profiles briefly flatten. The initial absorption profiles do not affect the resulting line and continuum emission, with the latter
persist for as long in the higher flux simulations. being important for emission far from the line core. It is
The cause of the redshifts observed in Figure 6(a) is less therefore possible that the majority of the differences between
obvious once the profiles undergo degradation by the the RADYN and RH profiles arise from the alternative methods
instrument, and it would be easy for the line profiles to be used for obtaining the level populations, and not because of the
misinterpreted as emitting an excess in the red wing as opposed application of PRD.
to being highly absorbing in the blue wing. While the beam is being injected, the additional heating and
direct excitation of the plasma by the beam act to increase the
amount of collisional excitation, which enhances the popula-
4.1.2. Velocities from RH
tions of the upper levels. This also increases the recombination
The atmosphere snapshots from the RADYN simulation rate as ne increases. As a result, the conditions in RADYN
were also used as input to RH. These were run with a six-level- closely approximate statistical equilibrium. However, once the
with-continuum hydrogen atom with each of the Lyman lines beam switches off, the dynamics of the atmosphere continue to
computed via PRD. The output RH profiles were then degraded evolve rapidly, while the level populations may take tens of
as before, and the simulated Doppler velocities were calculated. seconds to return to equilibrium levels (Carlsson & Stein 2002).
These velocities are displayed in Figure 6(b). This means that statistical equilibrium becomes a poor
The RH velocity profiles agree rather well with those approximation once the beam is switched off, leading to
obtained from RADYN throughout the beam-heating stage differences in the computed RH profiles with respect to
(t=0–20 s), with redshifts found in each of the lines and Ly-α RADYN. Further investigation is required, but it seems likely
and Ly-β again displaying the more prominent signatures. At that the enhanced electron density facilitates more collisions,
t=20 s, Ly-α suggests a peak downflow speed of 20 km s−1. such that the PRD solution approaches the CRD solution.
Significant differences arise between RADYN and RH after the Across both sets of results in Figure 6, the maximum speeds
electron beam is switched off. In RH, the shift in Ly-α decays of the apparent flows lie roughly between 20–40 km s−1. This
to zero, while the other Lyman lines abruptly transition into is in line with previous velocity measurements for chromo-
exhibiting strongly blueshifted signals. The blueshifted signa- spheric lines (Milligan & Dennis 2009; Brown et al. 2016).
tures peak at t=26 s, with Ly-γ and Ly-δ now exhibiting the However, the apparent flow direction does not reflect that of the
strongest flows; the former reaches speeds of around 50 km s−1. atmospheric velocity, because the line profiles are centrally
These signatures then decay over the remainder of the reversed. From Figure 3, the core formation height is upflowing
simulation. with a velocity of 50 km s−1 at t=20 s, while Figure 6 shows
To understand why the velocity profiles deviate so downflows at this time. Velocities obtained from intensity-
significantly from each other after t=20 s, Figure 7 shows a weighting are not as large as those obtained from Gaussian
comparison between the RADYN and RH profiles for Ly-α fitting, but do verify the direction of the observed shifts.
and Ly-γ. We include Ly-γ because it exhibits the strongest
blueshifted signatures in Figure 6(b). In addition to the RH
4.2. The F10D8 Simulation
profiles computed with PRD, we include those obtained from
RH when CRD is applied. It can be seen that, after t=20 s, To assess the impact of injecting energy at higher altitude,
both Ly-α and Ly-γ as computed from RH have much weaker we increased the spectral index to a value of 8 while keeping
wing intensities than the RADYN profiles, while the core other parameters fixed. Consequently, this electron beam has a
intensities show better agreement. sharper drop-off of electron number as a function of energy,
The drop in the wing intensities in Ly-α from RH is not as and therefore is constituted of a greater number of lower-
pronounced as in Ly-γ. The consequence of this is that the energy electrons, which are stopped at higher altitudes in the
RADYN profiles remain centrally reversed as a result of the model atmosphere.
high wing intensities, while the RH profiles, although still The evolution of the atmosphere in this simulation is shown
centrally reversed, are more prominently peaked due to the in Figure 8. Compared to the F10D3 simulation (Figure 2), the
lower intensities in the far wings (more so in the higher-order lower atmosphere (below z=1.7Mm) develops a much
lines). This means that, while RADYN continues to produce steeper temperature gradient throughout the first 20 s, which
redshifted profiles as a result of blueshifts acting in the dissipates as the transition region moves upwards. The steep
centrally reversed line cores, the blueshifts in the RH profiles temperature gradients are cospatial with localized enhance-
for Ly-β through Ly-δ act on more emissive features. ments in the electron density. Upflows are again initiated by the
Therefore, they strengthen the blue wing, which in turn means beam injection, with velocities in the lower chromosphere
that upflows are obtained in the velocity profiles. slightly lower (around 60 km s−1) than those of the F10D3
While the application of PRD was expected to produce line simulation.
profiles that differed from those computed by RADYN, it is We examine the formation of the Ly-α line in Figure 9 at
evident from Figure 7 that even the RH profiles computed times similar to those of the F10D3 simulation. At t=20 s, the
using CRD can deviate from their RADYN counterparts, with electron beam has just stopped heating the atmosphere and has
significant differences between the RADYN and RH (CRD) resulted in an upflow through much of the chromosphere. The

7
The Astrophysical Journal, 862:59 (19pp), 2018 July 20 Brown et al.

Figure 7. Ly-α and Ly-γ line profiles obtained from the F10D3 RADYN simulation, with RH output calculated via both CRD and PRD formalisms as a function of
time. To help illustrate the features causing the differences between Figures 6(a) and (b), we do not show these profiles with instrumental effects.

Figure 8. Atmospheric evolution during the F10D8 simulation. The pre-flare atmosphere is indicated by the thick black line. Quantities are plotted at 1.5 s intervals.

atmosphere has a peak velocity of almost 60 km s−1, but the positioned entirely within the blue wing, as a result of the core
velocity gradient above the transition region is steeper than in being formed in the presence of a 50 km s−1 upflow.
the F10D3 simulation, because the coronal layers are almost Curiously, the Ly-α line at t=20 s in the F10D3 simulation
stationary. As a result of the upflow, the τν=1 surface is also forms in the presence of an upflow of this speed, but the
noticeably asymmetric and the opacity structure of the moving blueshift in its central reversal is not as pronounced as in this
plasma is shifted to the blue as it is carried upward. The height model. From Figure 3(a), it can be seen that the line core forms
of this surface indicates that the line core is formed at over a larger vertical extent in the F10D3 simulation than in the
z=2.05Mm, whereas it was formed higher at the same time F10D8 simulation. In the latter, the line core forms in a very
in the previous simulation narrow layer. This perhaps suggests that the blueshift in the
The line contribution function shows that almost all of the line core of Ly-α in Figure 9(a) is more indicative of the dynamics,
is optically thick, with only small contributions from altitudes as it samples less of the atmosphere.
above the τν=1 surface. Once again, Sν exhibits a maximum at At t=25 s (Figure 9(b)), there is no longer any energy
an altitude below the formation height of the line core, meaning being injected and the atmosphere is relaxing. Plasma
that the resulting line wings contain a greater amount of emission continues to be carried upward by an atmospheric flow, which
than the line core and the line profile is centrally reversed. The still exhibits a prominent velocity gradient at coronal heights
central reversal in Ly-α is heavily blueshifted; in fact, it is and maintains a blue asymmetry in the τν=1 surface. The

8
The Astrophysical Journal, 862:59 (19pp), 2018 July 20 Brown et al.

Figure 9. Line contribution functions for Ly-α at two times during the F10D8 simulation. All lines and panels retain their meanings from Figure 3.

Figure 10. Doppler velocities of the F10D8 line profiles from RADYN (a) and from passing the RADYN atmospheres through RH (b). Some of the atmospheres
toward the end of the beam-heating stage did not converge in RH.

height of this surface has also increased, indicating the Ly-α profiles, obtained from both the RADYN profiles and those
line core now forms at a height of 2.3Mm. after calculation in RH, are shown in Figure 10. The velocities
As before, Sν peaks at an intermediate altitude between the from RADYN (Figure 10(a)) are dominated by persistent
core and wing formation heights, and it has undergone an redshifted signatures in all of the Lyman lines throughout the
overall decrease since t=20 s. As a result, the emergent Ly-α entire duration of the simulation.
line is centrally reversed and less intense than at t=20 s. The These redshifts suggest downflows of between 20–30 km s−1.
line contribution function at t=25 s indicates a slightly larger As in the F10D3 simulation, the source of these redshifts can
fraction of optically thin emission, relative to t=20 s, that is be easily understood by visually inspecting the line profiles.
mainly concentrated in the red wing. The upflow velocity at the Figure 9 makes it clear that, while the line cores are heavily
core-formation height has not diminished, so the central blueshifted, they are also centrally reversed—so when the
reversal still persists exclusively in the blue wing. profiles undergo convolution, it is the red wing that becomes
accentuated.
4.2.1. Velocities from RADYN
4.2.2. Velocities from RH
The Lyman lines from the F10D8 simulation are again
convolved with EVEʼs instrumental profile, and Doppler shifts The velocity profiles obtained from the RH solutions
are calculated in the Ly-α through δ lines. Resulting velocity (Figure 10(b)) are more complex. RH fails to converge

9
The Astrophysical Journal, 862:59 (19pp), 2018 July 20 Brown et al.

Figure 11. Atmospheric evolution during the 3F10D8 simulation. The pre-flare atmosphere is indicated by a thick black line. Quantities are plotted at 1.5 s intervals.

between t=10–20 s, likely as a result of the steep gradients in more prominent spikes in electron density in the upper
the atmosphere at these times. The velocity profiles throughout chromosphere, compared to the F10D8 simulation.
the first 10 s generally agree with those found from RADYN, We can see that the plasma reaches much higher speeds than
with redshift signatures suggesting downflows of around in both F10 simulations. We also note a more complex,
20 km s−1, although the initial RH solutions for Ly-γ and Ly- structured temperature profile, with a maximum of around
δ (t=1–2 s) are temporarily skewed because only half of each T=106 K being attained early in the simulation. Figure 11
line profile is modeled correctly. shows that the altitude of the transition region briefly decreases
As observed in the F10D3 simulation, after the beam-heating during the beam injection, but moves upward at later times.
stops, evolution of the line profiles becomes significantly The spikes in the electron density between t=10–20 s are
different with respect to the profiles computed from RADYN. concurrent and cospatial with relatively narrow but deep
As before, the higher-order lines eventually transition into troughs in the atmospheric temperature profile and enhance-
displaying blueshifts as the wing intensities in the line profiles ments in the total population of neutral hydrogen (Figure 11).
diminish, allowing the blueshifts in the line core to be retained These features move upward with time and could indicate a
after the instrumental convolution. However, Ly-α exhibits front of cool, dense material being pushed upward as
curious behavior, with velocities from Gaussian fitting briefly chromospheric evaporation takes place. The simultaneous
suggesting downflow speeds of almost 100 km s−1. enhancement in the total population of neutral hydrogen
Shortly after the beam stops heating the atmosphere indicates that this is a propagation of mass, and not just an
(t=20–25 s), the Ly-α line computed from RH continues to ionization front.
display a central reversal (which is still blueshifted), but We investigate the properties of the line contribution
crucially does not have comparable intensities in the blue and function for Ly-α in Figure 12. Because the higher-order
red wings. In RADYN, the wings on either side of the central Lyman lines behave similarly to Ly-α (and for conciseness),
reversal are typically similar in intensity. Conversely, in RH, we display only Ly-α in order to be able to present more
the emergent blue wing of Ly-α shortly after the beam injection timesteps in the simulation.
is notably less intense than the red wing. The combination of In Figure 12(a), we again observe line characteristics similar
the blueshifted central reversal and lack of an appreciably to those seen in Figure 3, with an asymmetric τν=1 surface
intense blue wing acts to produce a line profile that has a very and Sν peaking at an intermediate height between the core and
strong red asymmetry. This leads to the strongly redshifted wing formation heights. The plasma is slightly upflowing at
signatures in Figure 10(b). this point, and there appears to be complex structure within the
The F10D8 simulation shares similarities with its low-δ central reversal. Much of the wing emission is optically thin, as
counterpart, namely the presence of clear blueshifts in the indicated by the contributions to the intensity that appear above
centrally reversed line cores. As with the F10D3 simulation, the τν=1 surface.
these features are lost when the line profiles undergo Between the beamʼs peak-time and the beam switching off,
convolution and produce profiles with strong red asymmetries. the line formation becomes much more complex. At t=16.5 s
Velocities at the core formation height are similar to those in (Figure 12(b)), the line source function still has a maximum
the F10D3 simulation, but the altitude of this region is below the core formation height, so the profile is still centrally
generally lower in the δ=3 case. The velocity structure in the reversed. A secondary feature of the τν=1 surface has
atmosphere is comparably more complex, consisting of a sharp emerged. This feature originated from the original asymmetric
gradient at higher altitudes that is not observed in the F10D3 surface close to the line core, and it propagates through the blue
simulation. As in the F10D3 simulation, neither the direction or
wing, peaking at higher altitudes as a function of time.
the extent of the atmospheric velocity is recovered by the line
Furthermore, the maximum height of the secondary τν=1
profiles after convolution.
surface feature lies at the same height as the maximum velocity
in the atmosphere, which has grown considerably since t=6 s,
4.3. The 3F10D8 Simulation to almost 150 km s−1.
The overall beam flux in the 3F10D8 simulation (δ=8, This secondary feature in the τν=1 surface acts to
Ec=25 keV) is slightly enhanced above the previous two introduce an additional, highly blueshifted source of line
simulations. As with the F10D8 simulation, the high δ leads to emission that is linked in formation height to the peak speed of
the electrons being stopped higher in the chromosphere than the atmosphere. This secondary component is predominantly
the F10D3 simulation. The evolution of the atmospheric optically thick, mainly contributing emission at the same height
variables is shown in Figure 11, in which we can see much as the τν=1 surface. This additional source of emission is

10
The Astrophysical Journal, 862:59 (19pp), 2018 July 20 Brown et al.

Figure 12. Line contribution functions for Ly-α throughout the duration of the 3F10D8 simulation. All lines and panels retain their meanings from Figure 3.

likely related to the cool, dense front revealed by the features undergone an overall decrease, leading to a general reduction
seen in the temperature and density profiles in Figure 11. in emission from the line profile.
As the feature moves higher through the atmosphere By t=45 s (Figure 12(d)), the atmospheric velocity has
(Figure 12(c)), it propagates further through the blue wing of weakened in the region of the chromosphere and is extended over
the line as the atmospheric velocity increases. The highly a large range in height. The intensity of the line profile has
blueshifted “core” remains optically thick, while some diminished considerably, due to the low magnitude of the source
additional optically thin wing emission is also produced by function over the line formation height, but the secondary
the dense upflow. By t=20 s, the source function has blueshifted component persists and acts to strengthen the blue wing.

11
The Astrophysical Journal, 862:59 (19pp), 2018 July 20 Brown et al.

Figure 13. Doppler velocities of the 3F10D8 line profiles from RADYN (a) and from passing the RADYN atmospheres through RH (b). Note that many atmospheres
failed to converge in RH, so only the initial part of the simulation is solved in RH.

4.3.1. Velocities from RADYN structure of the atmosphere in this simulation is computation-
ally difficult for RH to solve, so only a subset of snapshots
Synthetic Doppler velocities are again calculated after
during the beam-heating stage converged. The first 10 s are
simulating the EVE instrumental response and measuring the
well-sampled, and display redshifts in the lines on the order of
line centroid shifts. The velocities obtained from analysis of the
10–20 km s−1, which agrees rather well with the RADYN
RADYN and RH profiles obtained from the 3F10D8 simulation
velocity profile (Figure 13(a)).
are plotted in Figure 13. Due to the very steep gradients in the
The presence of the high-density upflow evidently causes
atmospheric structure in this simulation, many of the snapshots
problems for convergence, as only two snapshots at around
could not converge in the RH code; ergo, the results in 13(b) do
t=17.5 s converge. These snapshots, however, suggest very
not constitute the entire length of the simulation as in 13(a).
weak redshifts in the lines. In Figure 13(a), we see very weak
The initial 10 s of the simulation display redshifts
blueshifts at this time, so it is likely that we are again seeing the
corresponding to downflows of around 15 km s−1 in Ly-α
effect of the secondary line component acting to temporarily
and Ly-β, and 5–10 km s−1 in the higher-order lines. As with
increase the amount of emission in the blue wing.
the F10D3 simulation, these perceived downflows are not due
It is interesting to note that, despite the brief (around
to redshifted emission, but rather are a consequence of the EVE
t=12–18 s) excursion to blueshifts due to the secondary line
instrumental profile masking a blueshifted core with an
component, the velocity profiles again register redshifts
asymmetric red wing. The profiles after instrumental convolu-
throughout the first 10 s as a result of the centrally reversed
tion are initially in absorption between t=0–3 s, transitioning
line cores being blueshifted. Maximum flow velocities of
into emission earlier than in the F10D3 simulation.
around 60 km s−1 are observed at late times in the higher-order
Between t=10–20 s, the effects of the dense upflow
Ly-δ and Ly-γ lines, with speeds across the lower-order Lyman
become apparent in the velocity profiles. While the secondary
series remaining around the 20 km s−1 level. The diminishing
line component is self-reversed, it also has narrow wing-like
of flows between t=10–20 s and the high upflow velocities at
enhancements apparent in Figure 12(b) that act to increase the
late times are both caused by the appearance and persistence of
overall amount of emission in the blue wing relative to the red
the blue-wing feature observed in Figure 12.
wing. This feature counteracts the reduction in blue-wing
It is worth noting that the production of the secondary line
emission from the primary central reversal and leads to the
component in this simulation is reliant on both the increased F
observation of faint blueshifts in the line between t=10–20 s.
and δ values, as similar features are not found in lower-flux
At t=20 s, the secondary line component has developed a
simulations with the same δ (F10D8 model) or in higher-flux
deeper reversal. The primary central reversal also remains
simulations with a lower δ (F11D3 model).
slightly blueshifted, and the combination of these two features
leads to a greater excavation of emission in the blue-wing. This
leads to the simulated observations once again obtaining 4.4. The F11D3 Simulation
redshifts of the order 10–20 km s−1. The final simulation in our study is a high-flux, hard beam
After t=20 s, the beam heating has stopped and the (δ=3, Ec=25 keV). Like the F10D3 simulation, this beam
atmosphere is in the process of relaxing. The line profiles are deposits a larger fraction of its energy in the lower chromo-
now faint, but the secondary blue-wing component does persist sphere, compared to the δ=8 simulations. The atmospheric
to an extent, eventually contributing an overall enhancement in evolution is shown in Figure 14, in which we see a prominent
emission to the blue wing. Therefore, the velocity profiles spike in the lower chromospheric temperature shortly after the
register blueshifts during this time. beam is injected, along with a brief compression of the
chromosphere as the transition region makes a small excursion
to lower altitudes. The transition region does not return to its
4.3.2. Velocities from RH
initial altitude by the end of our simulation. The atmospheric
In Figure 13(b), the velocities obtained from passing the temperature during beam heating is higher than that in
RADYN snapshots through the RH code are plotted. The the F10D3 simulation, attaining temperatures of almost

12
The Astrophysical Journal, 862:59 (19pp), 2018 July 20 Brown et al.

Figure 14. Atmospheric evolution during the F11D3 simulation. The pre-flare atmosphere is indicated by the thick black line. Quantities are plotted at 1.5 s intervals.
A complex velocity structure can be seen close to the transition region at t=46 s (in red).

T=106.5 K, as compared to T=105.9 K in the F10D3 now-upflowing plasma, and a slight blue-wing enhancement
c
simulation. can be seen in tn . The line contribution function indicates that
n
The flows initiated in this simulation are of a much greater the emission is produced in an extremely thin region cospatial
magnitude than in any of the other simulations, with upflows with the altitude from which the atmospheric flow changed
attaining speeds of almost 400 km s−1 at z=3Mm. The sharp direction. The resulting line profile is complex, with a dominant
temperature boundary at t=50 s indicates that the transition red-wing enhancement as a result of the core-formation height
region settles at an altitude of 1.4Mm. undergoing downflow. This is accompanied by a lesser
Near the end of the simulation (t=45–50 s), a sharp feature contribution in the blue wing, likely produced by the
of enhanced electron density appears close to the transition rebounding upflow. The line is clearly more intense than at
region. At this time, we also observe a downflow in the t=45 s, indicating that the complex dynamics at this time act
atmospheric velocity of around 100 km s−1. It is possible that to produce heightened levels of emission.
the presence of downflowing material results in a compression The level populations for hydrogen are plotted in Figure 16
of this region, which would enhance the electron density. for a series of timesteps around this feature, and a marked
The contribution functions for Ly-α at varying simulation increase is seen in all levels at the same time as the appearance
times are plotted in Figure 15. Again, because the other Lyman of the wing emission.
lines exhibit similar behavior in this simulation, we forgo Given the population enhancements and the rapid change in
plotting higher-order lines in favor of presenting more direction of the atmospheric velocity just above the τν=1
timesteps. surface, it is possible that a fast downward flow results in an
Panels (a) and (b) of Figure 15 reveal that, while the increase in the local plasma density, allowing for an increase in
transition region and corona are subject to a very fast the amount of collisional excitation at this height and
(v>200 km s −1) upflow, the height of core formation population of the upper levels. This then leads to emission
(z=1.5Mm) is stationary. This was also found at the core- by de-excitation. A fraction of this emission may be blueshifted
formation height of the Na I D1 in an F11 simulation by as the material is carried upward by the “rebounding” upflow.
Kuridze et al. (2016). Because of this, the line core is not
Doppler-shifted; instead, it is centrally reversed as Sν peaks
deeper in the atmosphere. At t=9 s, upflows reaching 4.4.1. Velocities from RADYN
20 km s−1 are present below the core-formation height, with The simulated Doppler velocity profiles for the F11D3
c
a very steep velocity gradient above. An enhancement in tn can simulation are plotted in Figure 17, again for both the RADYN
n
be seen at this time in the blue wing, as a result of these and RH output. The first 20 s in the RADYN velocity profile
upflows leading to an asymmetry in the τν=1 surface. (Figure 17(a)) are primarily dominated by blueshifts. The line
Curiously, a slight redshift is observed in the centrally reversed contribution functions for Ly-α at t=9 and t=20 s
line core at t=20 s, despite the lack of any downflow (Figures 15(a) and (b)) show profiles with shallow central
signature in the atmospheric velocity. At this time, the blue reversals, but with heightened amounts of emission in the blue
wing is also slightly more intense than the red wing, as a result wing relative to the red wing. These are directly linked to the
c
of an enhancement in tn . upflow in the atmosphere. While the blue wing is less
n
At t=45 s, the atmosphere is in the process of relaxation pronounced at t=20 s, it is accompanied by a slight redshift
(Figure 15(c)). The τν=1 surface has regained its symmetry in the centrally reversed core, which acts to further accentuate
and all parts of the line are in emission due to the source the emission in the blue wing. These factors lead to blueshifts
function peaking slightly above the core-formation height, being registered in the velocity profiles throughout the beam-
which is very close to the transition region. Additionally, the heating stage, which continue throughout the majority of the
line is largely symmetric, as it is not formed in the presence of simulation. As before, Ly-α and Ly-β provide more
any appreciable flows in the atmosphere. At this time, a pronounced signatures, suggesting upflow speeds of 23 and
downflow with a velocity of around 100 km s−1 can be seen 10 km s−1, respectively, when obtained from Gaussian fitting.
propagating down from the corona. Around t=45 s, all Lyman lines abruptly transition into
This downflow reaches the core-formation height, then displaying redshifts, with Ly-α indicating downflows of
abruptly changes direction as if rebounding. This can be seen in 20 km s−1. As a result of the change in direction of the
Figure 15(d), and it clearly has an effect on the Ly-α line. The atmospheric flow described earlier, it can be seen from
line source function is largely concentrated in the wake of the Figure 15(d) that the Ly-α line is strongly redshifted because

13
The Astrophysical Journal, 862:59 (19pp), 2018 July 20 Brown et al.

Figure 15. Line contribution functions for Ly-α throughout the duration of the F11D3 simulation.

of the core-formation height undergoing downflow. As a result of 4.4.2. Velocities from RH

this, the velocity profiles produce redshifts for the remainder of
the simulation. For this simulation, the flows are mostly linked to The velocities from RH (Figure 17(b)) closely match those
emitting features, and therefore the synthetic Doppler velocity obtained from RADYN. Blueshifts are again observed
profiles generally match the genuine flow direction in the throughout the beam-heating stage, with Ly-α and Ly-β
atmosphere. It should be noted that the redshift in the centrally exhibiting the fastest speeds of 25 and 12 km s−1, respectively.
reversed line cores visible at t=20 s will partially contribute to Ly-γ and Ly-δ suggest only weak upflows of around 5 km s−1.
the observed blueshifts in the velocity profiles at this time. Throughout the intermediate time (t=20–45 s), little to no

14
The Astrophysical Journal, 862:59 (19pp), 2018 July 20 Brown et al.

Figure 16. Level populations for hydrogen in the chromosphere around the time corresponding to the abrupt appearance of emission in the line wings (t=45 s)
during the F11D3 simulation.

Figure 17. Doppler velocities of the F11D3 line profiles from RADYN (a) and from passing the RADYN atmospheres through RH (b).

Doppler shifting is measured in the lines, whereas the RADYN red asymmetry in Hα during the early stage of their F11
profiles show sustained blueshifted signatures. This is not simulation was a result of the blueshift of the centrally reversed
dissimilar to the post-beam phase in the F10D3 simulation, line core. These results show that caution must be taken in
where the wing intensities for the RADYN and RH profiles associating line asymmetries to flows in the same direction.
evolved differently and led to differences in the velocity
profiles for the two codes (Section 4.1.2). 5. Simulation Of SPICE Profiles
From t=45–50 s, the RH velocity profiles again echo those
observed from the RADYN line profiles, with each of the The Solar Orbiter satellite will accommodate the Spectral
Lyman lines transitioning into exhibiting redshifts. In RH, Investigation of the Coronal Environment (SPICE) instrument
however, the Ly-α response to the change in direction of the (Fludra et al. 2013). SPICE covers two EUV wavelength
bands, one of which will include the Ly-β line. While we do
atmospheric flow is less pronounced than in the higher-order
not attempt to perform an exhaustive analysis of SPICEʼs
lines, which display redshifts corresponding to downflows of
future capability regarding Ly-β observations, we do explore
20 km s−1.
some of the basic concepts.
In summary, the F11D3 simulation remains dominated by
Much of this paper has focused on emulating the degradation
blueshifts in the Lyman lines while the beam is injected, which of detailed model line profiles via the EVE instrument. Here,
primarily originate from emitting features in the line profiles— we perform a similar analysis, but with SPICEʼs design
although this is bolstered by a slight redshift in the centrally parameters. The long-wavelength (LW) band will cover the
reversed line cores close to t=20 s. These blueshifted range of 97.3 nm<λ<104.9 nm. The spectrograph will
signatures transition into redshift when the downflowing disperse sunlight onto the detector at a resolution of 0.0083
plasma interacts with the core-formation height, which nm per pixel at 101 nm, and the line spread function will have
produces elevated levels of emission with distinctly redshifted an extent of four pixels, corresponding to a FWHM of around
profiles. 0.04 nm about the Ly-β line (Fludra et al. 2013).
In all simulations, the Lyman lines are capable of indicating Given these parameters, we convolve the Ly-β profiles from
the motion of plasma upflowing through the atmosphere. RADYN with a Gaussian that has an FWHM of 0.04 nm, and
Crucially, this behavior is not recovered in simulated observa- then rebin the resulting convolved profiles to a wavelength
tions if the line profiles are centrally reversed (as in the F10 spacing of 0.0083 nm per bin. We believe that this should
simulations), as the absorbing nature of the blueshifted line reasonably approximate the effects of SPICEʼs instrumentation
core acts to reduce the amount of emission in the blue wing on the line. Because the exposure time of the instrument can
rather than enhance it. Kuridze et al. (2015) also found that the vary, we do not include the effects of time integration. In

15
The Astrophysical Journal, 862:59 (19pp), 2018 July 20 Brown et al.

Figure 18. Snapshots of the Ly-β line during the 3F10D8 simulation with SPICE instrumentation effects (a) and the Doppler velocities obtained from measuring the
I
line centroid variations of the degraded line profile, plotted alongside a running measurement of the line asymmetry. The line asymmetry is quantified by IR , and is
B
shown for both the RADYN output (blue) and for the convolved profiles (green).

Figure 18, we show the results of SPICEʼs instrumentation on that from post-convolution, and it can be seen that SPICE
the Ly-β line profile and the resultant Doppler velocities for the performs reasonably well at detecting the asymmetries in this
3F10D8 simulation. We choose this simulation because the line.
resulting Ly-β line is particularly complex.
l + 0.58Å
In Figure 18(a), it is clear that the SPICE instrumentation ål =
0
l0 Il
allows the strengthened red wing in the line profile to be A= . (2 )
l - 0.58Å
detected during the first 10 s. It is not capable of resolving the ål =
0
l0 Il
secondary blue peak in the profile (around −0.5 Å) at t=16 s
SPICE is capable of resolving the central reversal of Ly-β at
(or any detailed feature at any time). However, it does result in
certain times during this simulation, and this does lead to
a relatively symmetric profile being produced around
deviations in velocity results at some points because the single
t=13–25 s, after which the red wing again becomes dominant
Gaussian becomes a very poor fit. Nonetheless, it is
as the secondary component contributes less emission. At
encouraging to find that these interesting central reversals
t=2 s, we note that the primary central reversal is also subtly
could be retained in SPICE data, although such features may
suggested by a slight dip in emission at the line core in the
appear weak and exhibit only shallow dips in the core intensity.
degraded profile. We find that, across all simulations, SPICE
In Figure 19, snapshots of the Ly-β line are shown at
may be capable of observing the central reversals at certain
t=18 s for the 3F10D8 model after convolution with
times in the evolution of the Ly-β line.
Gaussians of increasingly narrow FWHM values. The profiles
We find that the Doppler velocities obtained in Figure 18(b)
are then rebinned to SPICE’S wavelength spacing (0.083 Å per
are not very different from those in Figure 13(a), but while the
pixel). It can be seen that clear detection of the secondary line
shape of the velocity profiles remain similar, the results from
component and the central reversal in the line core would
SPICE allow greater maximum speeds to be measured. This is
require an instrument with a profile around 0.2 Å wide, while
unsurprising, as the more detailed profiles from the instrument
the central reversal is visible from FWHM values smaller
allow for Doppler shifts to be detected at a finer resolution than
than ∼0.25 Å.
from EVE. As in Section 4.3.1, strong blueshifts are found at
late times as a result of the persistence of the secondary line
6. Discussion and Conclusions
component, while the stationary line component has dimin-
ished in intensity. The modeling of these chromospheric lines has led to some
We define the line asymmetry (A) as the ratio of the emission interesting results across a variety of different beam injection
in the red wing to that in the blue wing (Equation (2)). In parameters. It should, however, be stressed that the RADYN
Figure 18(b), a running measurement of the asymmetry in Ly-β simulations are not without approximation. Flares are not one-
is shown for both the emergent line profile from RADYN and dimensional, single-loop structures. However, as a numerical

16
The Astrophysical Journal, 862:59 (19pp), 2018 July 20 Brown et al.

Figure 19. A snapshot of the Ly-β line at t=18 s from the 3F10D8 simulation, convolved with Gaussians of increasingly narrow FWHM values. The FWHM of the
Gaussian is shown in blue (in angstroms), while the resulting profile is rebinned to SPICEʼs wavelength spacing (0.083 Å per pixel).

experiment, the output from RADYN can be illuminating and and strongly blueshifted, leading to perceived downflows in the
used to obtain predictions. Further understanding could be synthetic velocity profiles. Compared to the δ=3 case, the line
gained from attempting a multithreaded approach, via the core appears to form in a narrower region in altitude, and it
addition of multiple beams over a spread of time with varying δ better samples the flow structure of the atmosphere.
and Ec values. The 3F10D8 simulation presents an interesting case of a
The simulations used in this paper cover four different beam- secondary source of line emission. An additional line
injection schemes, ranging from moderate to high flux, with a component is linked to a high-velocity atmospheric upflow.
variety of deposition heights. In this section, we outline the key While the secondary feature is self-reversed, its overall effect is
results before briefly revisiting each simulation in more detail. to enhance the blue wing and counteract the blueshifted
The key results are: absorption in the centrally reversed core. From Figure 3, it can
be seen that the velocity profiles throughout the heating stage
1. According to the simulations, the Lyman lines can be are initially not too dissimilar from those in the F10D3
blueshifted. simulation, with differences arising due to the initiation of the
2. The Lyman lines in the simulations often have a centrally secondary line source. It is interesting to note that a fast upflow
reversed core, as a result of Sν being large between the in the line contribution functions is clearly present between
core and wing formation heights. t=10–20 s (Figures 12(b) and (c)), but its presence is not
3. EVE—and instruments with similar properties—may not clearly detected in the Doppler velocity profiles because the
be able to detect blueshifts in the Lyman lines if their weak emissive features in the lines are lost when degraded to
cores are centrally reversed. In this case, the blueshifted EVEʼs resolution.
core is not resolved and the line will present an overall The F11D3 simulation presents an emitting feature during
red asymmetry. the relaxation process that appears to be linked to the flow
4. Accounting for non-equilibrium effects appears to be a structure in the atmosphere. Around 25 s after the energy
more important factor in obtaining consistent line profiles deposition stops, an atmospheric flow propagates down from
between RADYN and RH, with the effects of PRD being the corona and rebounds from a height very close to the
less significant (although not negligible). transition region, transitioning back to flowing upward. This
height is also where the majority of the line core emission is
The F10D3 simulation typically reveals blueshifted, cen-
being formed (Figures 15(c) and (d)), and Figure 16 indicates
trally reversed cores, as a result of Sν being greater at altitudes that the rebounding flow causes a rise in density that results in
below the core formation height. The overall effect of this is to collisional excitation of the levels in hydrogen. This leads to
remove blue-wing irradiance from the line profiles, such that line emission, which is affected by the flow structure, meaning
the line could be observationally interpreted as being red- that downflows are registered in the velocity profile.
shifted. Figure 4 shows that degradation of the lines by an These simulations suggest that the Lyman lines are sensitive
EVE-like instrument removes any delicate features, so it is to atmospheric upflows and that optically thin emission can be
clear that care should be taken when trying to interpret the flow produced in either red or blue wings, depending on the
direction obtained from Doppler-shifted EVE lines. In this circumstances, while the core emission is always optically
simulation, despite the line cores clearly being blueshifted, thick. More importantly, the simulations have illustrated the
convolution with the instrumental profile smears the asym- complications that may arise when interpreting the flow
metric red peaks and blueshifted central reversals, giving the direction from genuine observations of these lines.
overall appearance of a redshifted line core. This leads to By emulating the effects of instrumental response from EVE
velocity profiles (Figure 6) that could easily be interpreted as (used to measure Doppler shifts in the lines in Brown
suggesting downflows during the beam injection. et al. 2016), it is clear that one must be careful when it comes
The F10D8 simulation alters the deposition altitude of the to assigning a direction to a flow observed by a detector with a
electron beam, leading to differences in the resulting velocity wide instrumental profile. Our simulations show that an
structure of the atmosphere as compared to when a low-δ beam instrument like EVE is unable to observe flows in the centrally
is used. The atmosphere is upflowing, with a prominent reversed part of a line, as any self-reversal will not be retained
velocity gradient at the leading edge of the upflow. As with the after being affected by the instrumentation. Similarly to
F10D3 simulation, the line cores are both centrally reversed Kuridze et al. (2015), we find that a blueshift in the centrally

17
The Astrophysical Journal, 862:59 (19pp), 2018 July 20 Brown et al.

Table 1 However, care should be taken when interpreting observations of

List of RADYN Simulations Used these lines, as instrumental effects can misrepresent the detailed
Beama Fpeak (erg cm−2 s−1) δ Ec (keV) Ftot (erg cm−2) structure of the line profiles and the Doppler velocities may not
10
fully reveal the maximum speed or direction of an atmospheric
F10D3 1×10 3 25 1×1011
flow. Furthermore, time integration performed by the instru-
F10D8 1×1010 8 25 1×1011
3F10D8 3×1010 8 25 3×1011
mentation must be considered, as this represents an additional
F11D3 1×1011 3 25 1×1012 loss of information.

Note. S.A.B. is grateful for the support from an STFC research

a
The simulations used are publicly available on the grid of RADYN models at studentship. L.F. and N.L. are grateful for STFC funding under
https://star.pst.qub.ac.uk/wiki/doku.php/public/solarmodels/start. The F10D3, grant numbers ST/L000741/1 and ST/P000533/1. G.S.K. is
F10D8, 3F10D8, and F11D3 simulations have model numbers 55, 60, 66, and 67, supported by an appointment to the NASA Postdoctoral Program
respectively. at Goddard Space Flight Center, administered by USRA through
a contract with NASA. J.dlC.R. is supported by grants from the
Swedish Research Council (2015-03994), the Swedish National
reversed core of the Ly-α line (and higher-order Lyman lines)
Space Board (128/15), and the Swedish Civil Contingencies
can lead to an observed red asymmetry, which masks the true
Agency (MSB). The authors are grateful to M. Carlsson and the
direction of the flow. Similarly, the 3F10D8 simulation shows
F-CHROMA collaboration for the production and availability of a
that, in cases where there are flows indicated by the line profile,
grid of RADYN simulations. The research leading to these results
but the overall profile is symmetric, the flows may be difficult
has received funding from the European Communityʼs Seventh
to detect with the EVE instrument.
Framework Programme (FP7/2007-2013) under grant agreement
It also seems that the effects of assuming CRD on the line
No. 606862 (F-CHROMA), and from the Research Council of
profiles may not be too severe while the electron beam heats
Norway through the Programme for Supercomputing.
the atmosphere. Figure 7 shows that, during the relaxation
stage of the F10D3 simulation, the ratio of core to wing
ORCID iDs
emission is not comparable between RADYN and RH, and
even the RH profiles computed with CRD are significantly Graham S. Kerr https://orcid.org/0000-0002-5632-2039
different from those in RADYN. In order to produce more Nicolas Labrosse https://orcid.org/0000-0002-4638-157X
accurate model profiles for the Lyman lines, it would be Adam F. Kowalski https://orcid.org/0000-0001-
preferable to take non-equilibrium effects into account while 7458-1176
also applying the PRD formalism. Non-equilibrium effects
appear to be the most significant factor in causing the RH References
profiles to deviate from those in RADYN, but PRD remains
important during the decay phase. PRD should also be heeded Abbett, W. P., & Hawley, S. L. 1999, ApJ, 521, 906
Allred, J. C., Hawley, S. L., Abbett, W. P., & Carlsson, M. 2005, ApJ, 630,
during the heating phase, as it is likely to be highly dependent 573
on the electron density. Allred, J. C., Hawley, S. L., Abbett, W. P., & Carlsson, M. 2006, ApJ, 644, 484
These simulations have revealed a variety of interesting Allred, J. C., Kowalski, A. F., & Carlsson, M. 2015, ApJ, 809, 104
atmospheric features. From Figures 2, 8, 11, and 14, it is clear Benz, A. O. 2008, LRSP, 5, 1
Brown, J. 1971, SoPh, 18, 489
that the deposition of energy alters the altitude of the transition Brown, S. A., Fletcher, L., & Labrosse, N. 2016, A&A, 596, A51
region and leads to atmospheric upflows in all three Capparelli, V., Zuccarello, F., Romano, P., et al. 2017, ApJ, 850, 36
simulations. The flow speed generally increases with the Carlsson, M., & Stein, R. F. 1992, ApJL, 397, L59
amount of energy deposited, and the F10 simulations display Carlsson, M., & Stein, R. F. 1997, ApJ, 481, 500
Carlsson, M., & Stein, R. F. 2002, ApJ, 572, 626
the weakest flows. Crotser, D. A., Woods, T. N., Eparvier, F. G., Triplett, M. A., &
Furthermore, it is clear that the magnitude of the atmospheric Woodraska, D. L. 2007, Proc. SPIE, 6689, 66890M
flows is not fully represented by the apparent Doppler Dorfi, E. A., & Drury, L. O. 1987, JCoPh, 69, 175
velocities in the Lyman lines after convolution with the Emslie, A. G., & Sturrock, P. A. 1982, SoPh, 80, 99
instrumental profile; they provide only a sample of the Fisher, G. H. 1989, ApJ, 346, 1019
Fisher, G. H., Canfield, R. C., & McClymont, A. N. 1985, ApJ, 289, 414
atmospheric velocity. The flow speed of the Lyman lines is Fletcher, L., Dennis, B. R., Hudson, H. S., et al. 2011, SSRv, 159, 19
typically found to be several tens of km s−1 in these Fletcher, L., & Hudson, H. S. 2008, ApJ, 675, 1645
simulations, which is not implausible for these lines (Brown Fludra, A., Griffin, D., Caldwell, M., et al. 2013, Proc. SPIE, 8862, 88620F
et al. 2016). However, the apparent Doppler velocities rarely Hubeny, I., & Mihalas, D. 2014, Theory of Stellar Atmospheres (Princeton, NJ:
Princeton Univ. Press)
match the atmospheric speed (from RADYN) over the line Kerr, G. S., Fletcher, L., Russell, A. J. B., & Allred, J. C. 2016, ApJ, 827,
formation region and generally do not return the correct flow 101
direction when the line profiles are being centrally reversed. Korchak, A. A. 1967, SvA, 11, 258
Interpreting the correct flow direction from observations Kowalski, A. F., Allred, J. C., Uitenbroek, H., et al. 2017, ApJ, 837, 125
Kowalski, A. F., Hawley, S. L., Carlsson, M., et al. 2015, SoPh, 290, 3487
requires knowledge of whether the line profile is centrally Kretzschmar, M. 2011, A&A, 530, A84
reversed. While both red and blue line shifts are observed in Kretzschmar, M., Dominique, M., & Dammasch, I. 2013, SoPh, 286, 221
each of the Doppler velocity profiles, the initial atmospheric Kuridze, D., Mathioudakis, M., Christian, D. J., et al. 2016, ApJ, 832, 147
flows initiated in RADYN are typically upflowing. Kuridze, D., Mathioudakis, M., Simões, P. J. A., et al. 2015, ApJ, 813, 125
Leenaarts, J., Carlsson, M., & Rouppe van der Voort, L. 2012, ApJ, 749, 136
This work has shown that, similar to the results in Brown et al. Milligan, R. O., & Dennis, B. R. 2009, ApJ, 699, 968
(2016), the Doppler velocities of the Lyman lines during flare- Milligan, R. O., Gallagher, P. T., Mathioudakis, M., et al. 2006, ApJL,
type simulations suggest flow speeds of around 20–40 km s−1. 638, L117

18
The Astrophysical Journal, 862:59 (19pp), 2018 July 20 Brown et al.

Milligan, R. O., Kerr, G. S., Dennis, B. R., et al. 2014, ApJ, 793, 70 in the Space Instrumentation Era, ed. I. Dorotovic, C. E. Fischer, &
Rathore, B., & Carlsson, M. 2015, ApJ, 811, 80 M. Temmer (San Francisco, CA: ASP), 197
Reep, J. W., & Russell, A. J. B. 2016, ApJL, 818, L20 Taroyan, Y., & Bradshaw, S. J. 2014, SoPh, 289, 1959
Rubio da Costa, F., Kleint, L., Petrosian, V., Sainz Dalda, A., & Liu, W. 2015, Uitenbroek, H. 2001, ApJ, 557, 389
ApJ, 804, 56 Uitenbroek, H. 2002, ApJ, 565, 1312
Rybicki, G. B., & Hummer, D. G. 1991, A&A, 245, 171 Vernazza, J. E., Avrett, E. H., & Loeser, R. 1973, ApJ, 184, 605
Simões, P. J. A., Fletcher, L., Labrosse, N., & Kerr, G. S. 2016, in ASP Conf. Vernazza, J. E., Avrett, E. H., & Loeser, R. 1981, ApJS, 45, 635
Ser. 504, Coimbra Solar Physics Meeting: Ground-based Solar Observations Woods, T., Eparvier, F., Hock, R., et al. 2012, SoPh, 275, 115

19