{'fraction_non_alphanumeric': 0.030733366126997026, 'fraction_numerical': 0.01800799880231837, 'mean_word_length': 4.600431189363996, 'pattern_counts': {'":': 0, '<': 10, '<?xml version=': 0, '>': 3, 'https://': 0, 'lorem ipsum': 0, 'www.': 0, 'xml': 0}, 'pii_count': 0, 'substrings_counts': 0, 'word_list_counts': {'cursed_substrings.json': 0, 'profanity_word_list.json': 0, 'sexual_word_list.json': 0, 'zh_pornsignals.json': 0}}
Highly Efficient Modeling of Dynamic Coronal Loops

Observational and theoretical evidence suggests that coronal heating is impulsive and occurs on very small cross-field spatial scales. A single coronal loop could contain a hundred or more individual strands that are heated quasi-independently by nanoflares. It is therefore an enormous undertaking to model an entire active region or the global corona. Three-dimensional MHD codes have inadequate spatial resolution, and 1D hydro codes are too slow to simulate the many thousands of elemental strands that must be treated in a reasonable representation. Fortunately, thermal conduction and flows tend to smooth out plasma gradients along the magnetic field, so"0D models"are an acceptable alternative. We have developed a highly efficient model called Enthalpy-Based Thermal Evolution of Loops (EBTEL) that accurately describes the evolution of the average temperature, pressure, and density along a coronal strand. It improves significantly upon earlier models of this type--in accuracy, flexibility, and capability. It treats both slowly varying and highly impulsive coronal heating; it provides the differential emission measure distribution, DEM(T), at the transition region footpoints; and there are options for heat flux saturation and nonthermal electron beam heating. EBTEL gives excellent agreement with far more sophisticated 1D hydro simulations despite using four orders of magnitude less computing time. It promises to be a powerful new tool for solar and stellar studies.

transition region footpoints; and there are options for heat flux saturation and nonthermal electron beam heating. EBTEL gives excellent agreement with far more sophisticated 1D hydro simulations despite using four orders of magnitude less computing time. It promises to be a powerful new tool for solar and stellar studies.

Introduction
An abundance of observational and theoretical evidence indicates that much of the corona is highly dynamic and evolves in response to heating that is strongly time-dependent. The evidence further suggests that the cross-field spatial scale of the heating is very small, so that unresolved structure is ubiquitous. In particular, many if not all coronal loops are bundles of thin strands that are heated impulsively and quasi-randomly by nanoflares. It is estimated that a single loop contains several tens to several hundreds of such strands. See Klimchuk (2006) for a detailed justification of these ideas and references to relevant work.
Three-dimensional (3D) magnetohydrodynamic simulations are extremely useful for studying the source of coronal heating (instabilities of electric current sheets, reconnection, turbulence, etc.), but they cannot adequately address the complexity that is present in a single coronal loop, much less an entire active region. A more feasible approach is to treat the magnetic field as static and to solve the one-dimensional (1D) hydrodynamic equations along many representative flux strands using an assumed heating rate. The individual strands must be treated separately. It is not valid to approximate a loop as a monolithic structure with uniform heating corresponding to the average for the component strands. This gives a completely different and incorrect result.
There is reason to believe that the diffuse corona that lies between distinct bright loops is also comprised of elemental strands (e.g., Aschwanden et al. 2007). If roughly 100 strands are present in a single loop, then the numbers present in active regions and the global Sun are truly staggering. While it is possible to construct a limited number of model active regions with time-dependent 1D simulations (Warren & Winebarger 2007), it is not possible to investigate a wide range of values for the coronal heating parameters that must be assumed, such as the dependence on magnetic field strength, loop length, etc. (Mandrini, Démoulin, & Klimchuk 2000). This is a major limitation, since we are still struggling to identify the properties and physical origin of the heating mechanism. Progress in the foreseeable future must therefore rely on simplified solutions to the hydro equations that treat field-aligned averages and are much less computationally intensive. These are sometimes called "0D models" because there is only one value of temperature, pressure, and density at any given time in the simulation.
0D models have been developed previously by Fisher & Hawley (1990) and Kopp & Poletto (1993), but the best known is that of Cargill (1994). It has been used to study a variety of topics, including coronal loops (Cargill & Klimchuk 1997, 2006Klimchuk & Cargill 2001;Parenti et al. 2006), flares (Reeves & Warren 2002;Patsourakos, Antiochos, & Klimchuk 2002), post-eruption arcades (Reeves & Forbes 2005), and active stellar coronae . We have learned a great deal with the Cargill model, and our understanding has now advanced to the point where a more accurate and flexible model is required. This article presents an improved 0D model called Enthalpy-Based Thermal Evolution of Loops (EBTEL). As the name suggests, a key aspect of the model is an explicit recognition of the important role that enthalpy plays in the energy budget.
EBTEL improves upon the Cargill model in several important ways. First, whereas the Cargill model is limited to an instantaneous heat pulse, EBTEL accommodates any timedependent heating profile and can include a low-level background heating if desired. Second, EBTEL accounts for thermal conduction cooling and radiation cooling at all times during the evolution. The Cargill model assumes that only one or the other operates at any given time. Third, EBTEL has options for heat flux saturation and nonthermal electron beam heating. Finally, EBTEL is unique among 0D models in that it provides the time-dependent differential emission measure distribution of the transition region footpoints. Emission from the transition region plays a critical role in spatially unresolved observations, such as stellar observations and observations of the solar spectral irradiance, which is important for space weather (Lean 1997). Note that footpoint emission is not limited to the cooler (< 1 MK) plasma traditionally associated with the transition region. It can also include hot emissions that originate from the base of very hot loops. The so-called moss seen in the "coronal" channels of the Transition Region and Coronal Explorer (TRACE) is an example Martens, Kankelborg, & Berger 2000).

Coronal Model
Both 0D and 1D models of the corona are traditionally referred to as "loop" models. However, as we have discussed, there is good reason to believe that the loop-like intensity features seen in coronal images are actually comprised of many individual strands. EBTEL treats these individual strands, which are mini flux tubes in which the plasma is approximately uniform within the cross section. Multiple strand models can be combined to form a loop bundle.
Under static equilibrium conditions, the coronal portion of a strand is characterized by an exact balance between energy input (coronal heating) and energy losses by radiation and thermal conduction (Rosner, Tucker, & Vaiana 1978;Craig, McClymont, & Underwood 1978;Vesecky, Antiochos, & Underwood 1979). Some of the coronal heating energy-less than half-is radiated directly to space, and a heat flux carries the remainder to the transition region, from where it is more efficiently radiated. Temporal variations in the heating rate produce a well-defined response involving the transfer of mass between the chromosphere and corona. Heating increases cause the coronal temperature to rise and the downward heat flux to intensify. The transition region is unable to radiate the extra energy, and heated plasma flows into the strand in response to enhanced pressure gradients. This is the wellknow process of chromospheric "evaporation." The inverse process ("condensation") occurs when the heating rate decreases. As the coronal temperature declines, the reduced heat flux is insufficient to power the transition region radiation. The plasma cools, pressure gradients drop to sub-hydrostatic values, and material drains from the strand.
The basic idea behind EBTEL is to equate an enthalpy flux of evaporating or condensing plasma with any excess or deficit in the heat flux relative the transition region radiation loss rate. An excess heat flux drives evaporative upflows, while a deficient heat flux is compensated by condensation downflows. The key assumption of the model is that the radiative losses from the transition region and corona maintain a fixed ratio at all times. This ratio is the same one that applies during static equilibrium conditions. We defer justification of this assumption until later and now derive the equations that define the model.
We begin with the 1D time-dependent equation for energy conservation: where is the combined thermal and kinetic energy density; s is the spatial coordinate along the magnetic field; n, T, P, and v are the electron number density, temperature, total pressure, and bulk velocity, respectively; F is the heat flux; Q is the volumetric heating rate; g is the component of gravity along the magnetic field; and Λ(T ) is the optically thin radiative loss function, for which we use a piecewise continuous form based on the atomic physics calculations of J. Raymond (1994, private communication) and twice the coronal elemental abundances of Meyer (1985): 1.09 × 10 −31 T 2 , T ≤ 10 4.97 8.87 × 10 −17 T −1 , 10 4.97 < T ≤ 10 5.67 1.90 × 10 −22 , 10 5.67 < T ≤ 10 6.18 3.53 × 10 −13 T −3/2 , 10 6.18 < T ≤ 10 6.55 3.46 × 10 −25 T 1/3 , 10 6.55 < T ≤ 10 6.90 5.49 × 10 −16 T −1 , 10 6.90 < T ≤ 10 7.63 1.96 × 10 −27 T 1/2 , 10 7.63 < T. (3) The highest temperature range of the loss function is dominated by thermal bremsstrahlung (Cox & Tucker 1969). Equation (1) assumes a constant cross-sectional area, which is appropriate for distinct coronal loops and their constituent strands (Klimchuk 2000;Watko and Klimchuk 2000;López Fuentes, Klimchuk, & Démoulin 2006), but probably not for the diffuse corona. We also assume that the loop is symmetric, so only one half need be considered. We define s to increase from footpoint to apex. The downward heat flux is therefore a negative quantity. To simplify the discussion, we do not at this time include the energy and particle fluxes of a possible nonthermal electron beam. These will be added later.
If the flow is subsonic (v < C s = 1.5×10 4 T 1/2 = 2.6×10 7 cm s −1 at T = 3 MK) and the loop is shorter than a gravitational scale height (z apex < H g = 5.0×10 3 T = 1.5×10 10 cm at T = 3 MK), then the kinetic energy and gravity terms in equation (1) can be neglected, leaving 3 2 We now define the base of the corona, designated by subscript "0," to be the location where thermal conduction changes from being a cooling term above to a heating term below. This occurs at the top of a thin transition region, close to the chromospheric footpoint. Integrating equation (4) over the coronal portion of the strand and noting that the velocity and heat flux both vanish at the apex due to symmetry, we obtain The over bars indicate spatial averages along the coronal section, which has length L from the coronal base to apex. The first two terms on the right-hand-side of equation (5) are the enthalpy flux and heat flux at the coronal base, and R c is the radiative cooling rate per unit cross sectional area in the corona (erg cm −2 s −1 ). Since temperature, pressure, and density typically vary by less than a factor of 2 along the coronal section, the averages are quite characteristic of the entire section. We sometimes refer to them as simply the coronal values.
If we instead integrate equation (4) over the transition region, spanning from the top of the chromosphere to the base of the corona, we obtain a similar result: except that the spatial averages are now along the transition region, which has thickness l and radiative cooling rate R tr . In deriving this result, we made use of the fact the heat flux and enthalpy flux are both ignorable at the top of the chromosphere. During evaporation, a very small heat flux does in fact reach the top of the chromosphere, but most of the heat flux is dissipated throughout the transition region, heating each layer to the next higher temperature. No heat flux reaches the chromosphere during condensation.
Concerning the enthalpy flux at the top of the chromosphere, we note that conservation of mass requires that the electron flux be nearly constant through the transition region during both evaporation and condensation: Together with the ideal gas law, where k is Boltzmann's constant and we have assumed a fully ionized hydrogen plasma, equation (7) implies that the enthalpy flux is proportional to temperature. The enthalpy flux is therefore much smaller at the top of the chromosphere than at the base of the corona and can be safely ignored in equation (6).
Because the transition region is so thin, we can neglect the terms involving l in equation (6) and are left with When |F 0 | > R tr , there is an excess heat flux that drives a positive enthalpy flux (evaporation). When |F 0 | < R tr , there is a negative enthalpy flux (condensation) that combines with the heat flux to power the radiation. Static equilibrium corresponds to an exact balance |F 0 | = R tr .
Combining equations (5) and (9), we obtain the following equation for the evolution of the coronal pressure: We note that the same equation is obtained if we keep the terms involving l in equation (6) and interpretP andQ as the spatial averages along the entire strand, including both the transition region and coronal sections, with L then being the total length. Equation (10) reflects the energetics of the combined corona-transition region system. Energy enters the system only through coronal heating, and energy leaves the system only through radiation. Thermal conduction and flows transport energy between the corona and transition region, but they do not add or remove energy from the system.
It has been suggested that most "coronal" heating occurs in the transition region, in which case lQ tr >> LQ and equation (9) is not a good approximation. We do not believe this is a likely possibility, however. For one thing, the transition region is very thin. For another, it moves up and down a significant distance in response to changes in the spatially integrated heating rate within the strand (Klimchuk 2006). One might expect the positional dependence of the heating to be more closely related to the magnetic field than to the variable location of the transition region plasma. It is nonetheless possible to use EBTEL to study direct heating of the transition region. Its effect is very similar to that of a nonthermal electron beam, which is discussed in Sections 2.1 and 4.5. A minor difference is that an electron beam will slightly decrease the coronal mass.
We have defined the transition region to be the section of the strand where the heat flux is an energy source term. By this definition, its thickness is roughly 10% of the strand half length (l/L ≈ 0.1). This is not exceptionally small. However, most of the coronal heat flux is deposited within the extreme lower part of transition region, which is also where most of the radiation is emitted. We could therefore redefine the transition region to be much thinner and our model would be substantially unchanged. As an example, consider an equilibrium strand with an apex temperature of 2 MK and half length L = 7.5×10 9 cm. 80% of the heat flux is deposited over a distance of only l/L = 0.013, and 50% of the heat flux is deposited over an even shorter distance of l/L = 0.00069.
We wish to express the basic pressure equation, equation (10), in terms of the timedependent variablesP ,T , andn. We therefore approximate the radiative loss rate from the corona as R c ≈n 2 Λ(T )L.
This would be exact if the coronal density and temperature were perfectly uniform instead of approximately so.
Next, we assume that the radiative loss rates of the transition region and corona maintain a fixed ratio at all times: Since we want the model to apply during slow evolution as well as fast, c 1 should be equal to the static equilibrium value. One difficulty is that R tr /R c is different for different equilibrium conditions. In particular, it depends on the apex temperature of the strand, T a . Table 1 lists R tr /R c determined from exact equilibrium solutions in a semi-circular strand of half length L = 2.5×10 9 cm. Six apex temperatures ranging from 0.8 to 10.4 MK correspond to six different spatially-uniform heating rates. Except for lowest temperature case, R tr /R c increases monotonically with T a from 1.8 to 20.7 MK. In one implementation of EBTEL, we let c 1 vary according to a third order polynomial fit to these data. However, after some experimentation, we found that a constant value c 1 = 4.0 provides the best overall agreement with 1D simulations, especially in cases of impulsive heating. Table 2 lists R tr /R c for a longer equilibrium strand with L = 7.5×10 9 cm. The ratio is reasonably close to 4 for apex temperatures ranging from 1-4 MK. All of the results presented in this paper use a constant value c 1 = 4.0.
We do not yet have a compelling physical argument for why R tr /R c should be constant even when the strand is far from equilibrium. Fortunately, this does not appear to be an important assumption, at least not for the simulations presented in this paper. It is certainly not important during times of strong evaporation, when the evolution is essentially a balance between the downward heat flux and upward enthalpy flux, and radiation plays no significant role. The radiative losses are only 10 −3 of the heat flux during the strong evaporation phase of the nanoflare simulation of Section 4.1 (example 1). Radiation is very important during times of strong condensation, on the other hand. The assumption that R tr /R c equals the equilibrium value could then potentially cause problems. It turns out that condensation is fairly mild in all of our example simulations in the sense that the radiative losses never greatly exceed the heat flux. They are at most a factor of 3.6 larger in example 1.
We have considered a wide variety of heating scenarios, some discussed in Section 4 and others not reported in this paper. The fact that EBTEL is able to reproduce exact 1D solutions as well as it does gives us considerable confidence in the approximations of the model, including the assumption that R tr /R c = 4.0.
With equations (11) and (12), we can now express equation (10) for the evolution of coronal pressure in terms of the fundamental variablesP ,T , andn.
We next move on to an equation for the coronal density. The total mass contained in the coronal section of the strand changes as material evaporates and condenses. Specifically, the time derivative of the electron column densitynL (electrons per unit cross sectional area) is equal to the flux of electrons through the coronal base: This can be derived trivially by integrating the 1D equation of mass conservation from the base of the corona to the apex. Combining equations (7), (8), and (9), we can write the electron flux as Substituting into equation (13), we get where we have introduced c 2 for the ratio between the average coronal temperature and apex temperature, and c 3 for the ratio between coronal base temperature and apex temperature, Tables 1 and 2 list the values of these two ratios for the exact equilibrium solutions in the short and long strands, respectively.T /T a is very close to 0.87 in all cases, while T 0 /T a varies over a fairly narrow range from 0.22 to 0.61. In the implementation of EBTEL with variable c 1 , we also let c 3 vary based on a polynomial fit to the data in Table 1. As already indicated, constant values give better overall agreement with the 1D simulations; c 3 = 0.5 seems to work best and is the value used in the examples presented here. Note that c 2 = 0.87 is not far from 7/9, which corresponds to a constant heat flux solution.
Using equations (11) and (12), we can express R tr in terms of our fundamental variables n andT , but we still need an expression for F 0 . The classical expression for the heat flux is where κ 0 = 1.0×10 −6 in cgs units. Noting that we can approximate the heat flux at the base as where T a =T /c 2 . The precise value of the coefficient depends on the details of the temperature profile; 2/7 corresponds to a constant heat flux, while 4/7 corresponds to a constant heat flux divergence.
The classical heat flux is, however, unphysically large during times of exceptionally high temperature and/or exceptionally low density, such as during the earliest phase of an impulsive heating event. Under these conditions, the heat flux saturates at approximately where m e is the electron mass and β is a flux limiter constant that we set to 1/6 (Luciani, Mora, & Virmont 1983;Karpen & DeVore 1987). We consider two possibilities in our simulations. First, we set F 0 = F c at all times, regardless of the temperature and density. Second, we use the form which reduces to F c when |F c | ≪ |F s | and to F s when |F c | ≫ |F s |. We can now express equation (15) for the coronal density evolution in terms of the fundamental variables.
The last governing equation, for the coronal temperature evolution, follows straightforwardly from the ideal gas law: Note that this is not exact because the ideal gas law is not exact when average values of P , T , and n are used.
Summarizing, the coronal part of EBTEL is defined by: evolutionary equations (10), (15), and (23); the assumption given by equation (12); the approximations given by equations (11) and (20) or (22); and parameters c 1 = 4.0, c 2 = 0.87, and c 3 = 0.5. Note that c 2 and c 3 appear only together as a ratio in equation (15), so there are really only two parameters in the model.
The plasma velocity at the base the corona can be obtained straightforwardly from the electron flux, equation (14), according to Using T in place ofT gives the velocity at that temperature in the transition region.

Nonthermal Electron Beam
The mechanism that directly heats the coronal plasma may also produce energetic particles. It is thought, for example, that a sizable fraction of the total energy of a flare goes into nonthermal electrons (Saint-Hilaire & Benz 2002;Emslie et al. 2005). We therefore have incorporated a nonthermal electron beam into EBTEL. We assume that the electrons originate from the existing strand plasma and stream freely along the magnetic field to the coronal base. We further assume that all of the beam energy goes into the enthalpy of evaporating plasma and that any chromospheric radiation that may be produced is negligible. Because we do not consider the details of the energy deposition (i.e., how it depends on column depth), our calculation of the differential emission measure of the transition region (Section 3) is not reliable when nonthermal electrons are included.
The effect of the electron beam on the coronal energy budget is straightforward. The corona gains energy from the enthalpy of the evaporated plasma, but it loses energy because electrons must be removed from the thermal pool to supply the seed particles for the beam. In general, the gain far exceeds the loss because the mean energy of the accelerated electrons, E, is much greater than their thermal energy, (3/2)kT .
If F and J are the energy flux and particle flux of the beam, respectively, so that then we must modify our equations by subtracting F from the right of equation (9), adding (3/2)kT J to the right side of equation (5), and adding J to the right side of (13). Note that F and J are both negative quantities. The evolutionary equations for pressure and density are then and We have avoided the tricky issue of electron return currents. Quasi neutrality of the plasma requires either that protons are accelerated with the electrons, which is thought to be unlikely, or that an electron return current replenishes the electrons lost to the beam. With a return current, the unity term inside the last set of parentheses disappears from equation (27), and the temperature in the kT /E term of equation (26) must be replaced by the temperature difference between the strand plasma and the replenishing electrons. Note that both of these terms are negligible as long as E ≫ kT , even without a return current.

Differential Emission Measure
Most observed plasmas are multi thermal, even within a single observational pixel. An important quantity is therefore the differential emission measure, DEM(T ), which describes how the plasma is distributed in temperature. Spatial variations in the coronal temperature tend to be greatest across the magnetic field. In a multi-stranded loop bundle, for example, the different strands will have different temperatures if the heating is steady but unequal or if it is impulsive but out of phase. There is also some temperature variation along the field. We here consider the differential emission measure of a single stand of unit cross sectional area: The transition region is treated carefully in the next section. For the corona, we make the crude approximation that the total emission measure, 2Ln 2 , is distributed uniformly over the temperature interval 0.74T a ≤ T ≤ T a . The average temperatureT falls exactly in the middle of this interval. Note that our approximation is not critical, because the DEM(T ) of a coronal observation is determined primarily by the distribution of different strands rather than by the variation along each strand.

Transition Region Model
The situation is very different in the transition region, where the temperature and density vary dramatically over a short distance along the magnetic field. Emission from the transition region is a critically important component in many observations, such as spatially unresolved observations of stars or measurements of the full Sun spectral irradiance. Even high-resolution observations on the solar disk tend to have lines of sight that pass through both coronal and transition region plasmas. We do not attempt to model the detailed spatial structure of the transition region, but instead deal directly with the differential emission measure. We have developed two separate approaches. The one we now discuss is the easiest to implement and has been used for all of the examples shown in the paper. The second approach, presented in the Appendix, has the advantage of being physically more revealing. It treats the limiting cases of strong evaporation, strong condensation, and static equilibrium, and provides simple analytical expressions for DEM(T ) in each case. The two approaches produce similar results. Neither is valid when nonthermal particles are important.
We begin with the steady state version of the energy equation, equation (4), in the absence of local heating: which should be approximately correct in the transition region. We next assume that the heat flux term can be approximated by −κ 0 T 3/2 (∂T /∂s) 2 . This is strictly valid when the scale lengths of temperature and heat flux are the same. Rewriting the enthalpy term using the ideal gas law and constant mass flux, and assuming that pressure is the same in the corona and transition region, the energy equation becomes This is quadratic in ∂T /∂s and can be solved trivially. DEM(T ) then follows directly from equation (28).
We have computed the radiative loss rate from the transition region by integrating the product DEM(T ) × Λ(T ) over the transition region temperature interval and find that it is similar to the value obtained from R tr = c 1 R c . The only exception is during times of strong evaporation, when the c 1 assumption is unimportant.

Example 1
We have coded up EBTEL in the Interactive Data Language (IDL) and now examine several simulations that were run on a desktop computer. The first example considers an impulsive energy release in a static equilibrium strand of half length L = 7.5×10 9 cm. An average coronal temperature of 0.52 MK in the initial equilibrium is produced by a heating rate of 10 −6 erg cm −3 s −1 . We obtain the equilibrium by guessing at the values ofT ,n, andP using scaling law theory (Rosner, Tucker, & Vaiana 1978;Craig, McClymont, & Underwood 1978) and letting the strand evolve while holding the heating rate constant.T changes very little during the relaxation, whilen andP decrease by about a factor of two.
We impose a nanoflare energy release on top of the steady background heating. It has a triangular profile with a total duration of 500 s and a peak value of 1.5×10 −3 erg cm −3 s −1 , 1500 times stronger than the background. Nonthermal electron beams are excluded from all but the last our examples. The solid curves in Figure 1 show howT ,n, andP respond to the event. The generic behavior is well documented (Cargill 1994;Klimchuk 2006). Temperature and pressure rise abruptly as the nanoflare energy is converted into thermal energy at a roughly constant density. An intense heat flux drives strong evaporation and the strand begins to fill with plasma. The temperature then declines as the nanoflare shuts off, but evaporation continues, and the peak density is not reached until well after the nanoflare has ended. Radiation becomes progressively more important as the temperature falls and density rises. It eventually takes over from thermal conduction as the dominant cooling mechanism. The strand then enters a long phase of draining and condensation.
We have run an exactly corresponding simulation with our sophisticated 1D hydro code called the Adaptively Refined Godunov Solver (ARGOS). As described in Antiochos et al. (1999), the code uses an evolving numerical mesh to resolve steep gradients wherever they may occur. We use the same radiative loss function used in EBTEL. For the 1D simulation, we make the additional assumptions, not required of EBTEL, that the strand is semi-circular, lies in a vertical plane, and is heated in a spatially uniform manner. We earlier defined the boundary between the corona and transition region to be the location where the divergence of the heat flux changes sign. This is not practical in the 1D simulation due to the more complicated temperature structure associated with waves and even shocks that are excited by the impulsive energy release. We therefore compute coronal averages by averaging over the upper 80% of the strand. These averages are indicated by dashed lines in Figure 1. The small wiggles are due to the aforementioned waves.
There is good agreement between the EBTEL and 1D results. This is highly encouraging given that EBTEL requires approximately four orders of magnitude less computing time. This run took only about 10 s. The biggest differences are in the density and pressure, where the EBTEL values are about 20% too high for the first 2000 s. The way that the 1D averages are computed is a contributing factor, since density and pressure are highest in the lower part of the strand leg that is excluded from the averages. Another contributing factor is that EBTEL assumes that all plasma energy is thermal ( 3 2 P ). In fact, some of the energy is kinetic, so the pressure is artificially inflated. Figure 2 shows the differential emission measure distribution for the full strand averaged over the first 10 4 s of the simulation. Both the corona and transition region are included. One reason for averaging over time is to simulate the observation of a multi-stranded loop. If the strands are heated randomly, then the time average of a single strand is equivalent to an instantaneous snapshot of an unresolved bundle. As long as the strands get reheated after they cool, then at all times there exists one strand in the bundle for each small time interval from the full simulation. As in Figure 1, EBTEL is represented by the solid curve and the 1D model is represented by the dashed curve. The agreement is once again very encouraging, especially considering that the DEM(T ) spans more than 3 orders of magnitude.
The EBTEL values are too high by factors of 2-3 at the higher temperatures. This is partly because of the enhanced densities discussed above and partly because temperature decreases more gradually below 3 MK in the EBTEL simulation (see Figure 1). The slower cooling rate could be because our assumption of a constant c 1 = 4.0 is not quite correct. We are currently investigating this issue.
The DEM(T ) plotted here and defined in equation (28) indicates the amount of plasma, n 2 ∆s, that is present in temperature interval ∆T . It has units of cm −5 K −1 . Some authors instead use which indicates the amount of plasma present in the logarithmic temperature interval ∆ ln T and has units of cm −5 . The two definitions differ by a factor T : DEM ln (T ) = T ×DEM(T ). Figure 3 shows DEM ln (T ). Still other authors define the differential emission measure in terms of the base 10 logarithm: DEM log (T ) = ln(10)×T ×DEM(T ). Figure 4 separates out the contributions to DEM(T ) from the coronal (dashed) and transition region (dot-dashed) sections of the EBTEL simulation and the coronal section of the 1D simulation (dotted). The transition region contribution of course dominates at low temperatures, but it is also significant at hotter temperatures that are normally associated with the corona. The transition region and coronal contributions are equal at T = 1.0 MK, which is approximately 1/4 the temperature of the hottest significant emission measure. Although it is difficult to observe the transition region in isolation from the corona, since lines of sight that reach the transition region must pass through the corona, it is easy to observe the corona in isolation from the transition region simply by looking above the limb. The agreement between the coronal DEM(T ) curves from the EBTEL and 1D simulations is reasonably good, except below 0.5 MK, where the EBTEL values are unreliable (since the DEM is rather arbitrarily cut off at 0.74T a ; Section 2.2). We discuss the agreement with actual observations in Section 5.
The results presented above assume a classical heat flux at all times. We have repeated the simulations with a heat flux that is allowed to saturate according to equation (22). The results differ only early during the nanoflare energy release, when saturation limits the thermal conduction cooling and the average coronal temperature rises to maxima of 4.7 MK and 7.6 MK in the EBTEL and 1D simulations, respectively (peak apex temperatures are of course higher). The densities are very low at this time, however, so the time-averaged DEM(T ) is minimally affected. The curves are nearly indistinguishable from those in Figure  2, the only difference being that the high temperature tail is extended by about 0.02 in the logarithm. We note that very hot emission, though generally very faint, provides the best diagnostics of nanoflare properties . Care should be taken to include heat flux saturation when studying the hottest emission.

Example 2
We next consider a much more impulsive nanoflare. It has the same total energy as the first example (5.625 ×10 9 erg cm −2 ), but the duration is ten times shorter (50 s) and the amplitude is ten times larger. This scenario provides a much better comparison with the Cargill model, in which all of the energy is deposited instantaneously. Figure 5 shows the evolution ofT ,n, andP as given by EBTEL (solid), the 1D model (dashed), and the Cargill model (dotted). Only classical heat flux results are presented, since the Cargill model does not include saturation effects. EBTEL again reproduces the 1D results quite well, although the densities and pressures are about 50% too high over the first half hour.
The EBTEL and Cargill results differ in several important respects.T reaches a maximum of 8.1 MK in the EBTEL model, similar to 7.7 MK in the 1D model, whereas the Cargill model peaks at only 4.0 MK. The Cargill model predicts substantially higher densities and pressures throughout most of the simulation. The primary reason is the assumption that radiation is ignorable during the first phase of cooling (ending at 1700 s). Since radiation is the only mechanism by which energy can leave the system, the thermal energy density and therefore the pressure are constant. Another reason for the excess pressures in both the Cargill and EBTEL models is the neglect of kinetic energy. All of the plasma energy is assumed to be thermal. This is reasonable only when the Mach number is small. The Mach number is mostly less than 0.15 after 500 s in the 1D simulation, but there are locations in the strand where it approaches 3 shortly after the nanoflare ends.
A final difference in the Cargill model is the prediction of a catastrophic cooling late in the evolution, at approximately 8000 s. This is not present in either the EBTEL or 1D simulations and is a consequence of the fact that no background heating is possible in the Cargill model. The radiative loss function, Λ(T ), is such that a thermal instability causes the temperature decline to accelerate in the Cargill model until a pre-set limit is reached (usually 0.1 MK). In the EBTEL and 1D models, the temperature asymptotically approaches the static equilibrium value corresponding to the background heating rate. We note that the 0D model of Fisher & Hawley (1990) predicts a catastrophic cooling even in the presence of background heating, but this appears to be a spurious result, at least in some cases.

Example 3
As a third example, we consider a qualitatively different heating scenario. The strand begins in static equilibrium with a uniform heating rate of 2×10 −4 erg cm −3 s −1 . The heating rate is slowly reduced by a factor of 100 over a period of 50,000 s, as shown at the bottom of Figure 6. It is maintained at the reduced level for 5000 s, then suddenly increased to the original level over 100 s. It is maintained at that level for 3900 s, then suddenly decreased again over 100 s. It remains at the reduced level for the remainder of the simulation.
The top two panels of the figure show the evolution of temperature and density for EBTEL (solid) and the 1D model (dashed). The 0D solution tracks the 1D solution very well. Temperature is systematically high, but the detailed shapes of both the temperature and density profiles are faithfully reproduced. This shows that our assumption c 1 = 4 is reasonable for situations other than impulsive heating.

Example 4
The final two examples are modifications of example 1. The 500 s nanoflare is ten times more intense in example 4. Figure 7 shows the time-averaged DEM(T ) curves for the whole strand (solid), corona (dashed), and transition region (dot-dashed). Note that the coronal curve is strongly peaked near 3 MK, as observed in active regions, which we return to shortly.

Example 5
Example 5 differs from example 1 only in the form of the nanoflare energy release. Half of the nanoflare energy is assumed to go into direct plasma heating, and the other half is assumed to go into nonthermal electrons with a mean energy of 50 keV. Figure 8 shows the coronal DEM(T ) curve (solid) together with the corresponding curve from example 1 (dashed). The curves are nearly identical except that plasma hotter than 3 MK is missing from example 5. The reason for this difference is as follows. The amount of evaporated material is determined largely by the total energy that is released, regardless of its form. Temperature, on the other hand, depends strongly on the form of the energy release. With direct plasma heating, the coronal temperature rises until either the nanoflare ends or the downward heat flux balances the nanoflare heating rate. In contrast, nonthermal electrons have no direct effect on the coronal temperature. Note that direct heating of the transition region, discussed after equation (10), has a similar effect on the coronal evolution as does a nonthermal electron beam.

Additional Tests
We have tested EBTEL against two other 1D hydro codes and found good agreement in both cases. F. Reale kindly simulated example 1 using the Palermo-Harvard code (Peres et al. 1982), and K. Reeves kindly simulated a loop-top flare with a peak temperature of 29 MK using the NRLFTM code (Mariska et al. 1982). It is interesting that the plasma evolution is similar even though the Palermo and NRLFTM codes use different radiation loss functions than do EBTEL and ARGOS. This shows that the precise form of the loss function is not important whenever the heating is impulsive.

Discussion
As evidenced by these examples, our simple 0D model is an excellent proxy for more sophisticated and far more computationally intensive 1D hydro simulations. It improves substantially on the 0D models of Cargill (1994), Fisher & Hawley (1990), and Kopp & Poletto (1993). The Cargill model assumes that heating is instantaneous and that cooling occurs either by thermal conduction or by radiation, but not by both at the same time. The Fisher-Hawley model: (1) predicts abrupt evolutionary changes as the strand evolves between three distinct regimes; (2) does not account for the evaporation that continues well beyond the end of an impulsive heating event; and (3) cannot return to the pre-event state due to unphysical catastrophic cooling. The Kopp-Poletto model shares some similarities with EBTEL, but it treats the flows in a fundamentally different way. Like EBTEL, it equates the enthalpy carried by evaporative upflows with an excess heat flux, but the excess is determined relative to the pre-event state, rather than to the time-varying radiative losses from the transition region. Condensation downflows in the model are given by a density-dependent fraction of the free-fall velocity. In actuality, gravity plays no direct role in condensation, since the downflows are driven by pressure gradient deficits relative to hydrostatic equilibrium, in the same way that evaporative upflows are driven by pressure gradient excesses. Gravity sets the value of the hydrostatic gradient, but it is only the deficit or excess relative to this value that is important for the flows. Inclined strands experience essentially the same condensation and evaporation as do upright strands of the same length. Finally, EBTEL has advantages over all three of the other models in that it provides the DEM(T ) of the transition region and treats nonthermal electron beams and heat flux saturation.
One obvious application of EBTEL is to investigate the idea that the basic structural elements of the corona are very thin, spatially unresolved magnetic strands that are heated impulsively. Loops may be bundles of such strands, as reviewed in Klimchuk (2006), and the diffuse corona may be similarly structured. Differential emission measure distributions are one important test of this idea. Observed DEM(T ) curves from active regions and the quiet Sun tend to be peaked near 10 6.5 and 10 6.1 K, respectively, and to have a slope (temperature power law index) ≥ 0.5 coolward of the peak (Raymond & Doyle 1981;Dere & Mason 1993;Brosius et al. 1996). This is consistent with the coronal DEM(T ) curves of examples 1 and 4 (Figures 4 and 7). The full loop curves are discrepant, on the other hand, due to the strong contribution from the transition region. The cited observations were made on the disk and should in principle include the transition region component. However, it is possible that absorption from chromospheric material such as spicules significantly attenuates the intensities of transition region lines used to construct the DEM(T ) curves (e.g., Daw, DeLuca, & Golub 2005;DePontieu et al. 1999;Doschek & Feldman 1982;Schmahl & Orrall 1979). We are currently investigating the magnitude of this effect.
One of the great mysteries of coronal physics that has come to light in the last few years is the discovery that warm (∼ 1 MK) coronal loops are much denser than expected for quasi-static equilibrium and live for much longer than a cooling time. The loops are therefore neither steadily heated nor cooling as monolithic structures. It has been shown that the observed densities and timescales can be explained by bundles of nanoflare heated strands, as long as nanoflares do not all occur at the same time (see Klimchuk 2006 and references cited therein). Neighboring strands will therefore have different temperatures, and loops are predicted to have multi-thermal cross sections. In particular, emission should be produced at temperatures hotter than 3 MK. Hot loops are sometimes observed at the locations of warm loops, but not always. Example 5 suggests that nonthermal electron beams are a possible explanation for the lack of hot emission. As we have discussed, beams can produce excess densities through evaporation without the need for high temperatures. We have just begun to explore this possibility. For now, we note that the coronal DEM(T ) curve of example 5 ( Figure 8) bears a close resemblance to the observed curves reported by Schmelz et al. (2001) for a loop seen above the limb.
In conclusion, EBTEL is a powerful new tool that can be applied to a variety of problems in which large numbers of evolving strands must be computed. For example, it is now feasible to construct multiple models of nanoflare-heated active regions or entire stars and therefore to examine a wide array of nanoflare parameters (magnitude, lifetime, occurrence rate, dependence on field strength and strand length, etc.). By determining which parameters best reproduce the observations, we can place important constraints on the heating and thereby gain insight into the physical mechanism (e.g., Mandrini, Démoulin, & Klimchuk 2000;Schrijver et al. 2004;Warren & Winebarger 2006). EBTEL is currently being used to study the emission characteristics of coronal arcades (Patsourakos & Klimchuk 2007b), to explain the light curves of solar flares (Raftery, Gallagher, & Milligan 2007), and to model coronal loops as self-organized critical systems (López Fuentes & Klimchuk 2005; We can combine these three limiting cases into a single expression that applies at all times: (A13) This expression reduces to the desired forms in the relevant limits. We have confirmed that the temperature dependencies in equations (A3), (A7), and (A12), are present in the differential emission measure distributions from 1D simulations.      , and transition region (dot-dashed) differential emission measure distributions for an EBTEL simulation of a nanoflare that is ten times larger than that of example 1.