Modeling of AR 12760 with GX Simulator and Evidence for the Extended Transition Region in Peripheral Active Region Loops

TL;DR

This study models a solar active region using GX Simulator to understand heating mechanisms and finds that the transition region extends farther along longer loops than previously assumed, affecting EUV emission modeling.

Contribution

It introduces a new modeling approach that accounts for the extended transition region in longer loops, improving the accuracy of EUV emission predictions in active regions.

Findings

Best fit heating model depends on magnetic field and loop length

Models underestimate emission in the legs of longer loops

Extended transition region significantly influences EUV emission

Abstract

In order to understand solar atmospheric heating it is important to test heating models against spatially resolved data from solar active regions. Here we model a small active region, AR~12760 observed on 2020 April 28, with the GX Simulator package by fitting the EUV intensities in wavebands observed by Solar Dynamics Observatory's Atmospheric Imaging Assembly. We assume the temporally and spatially averaged heating rate along a loop has a power-law dependence on loop length, $L$ and average magnetic field strength along the loop, $B_{avg}$. We find that the best fit heating model for the 211~\AA\ band is $<\!\!Q\!\!>\approx 7\times 10^{-3} (B_{avg}/{100 \mbox{G}})^{1.5}(L/{10^9\mbox{cm}})^{-1}$ erg cm$^{-3}$ s$^{-1}$ but that there is a range of parameters that give qualitatively reasonable fits, which we conclude is due to a correlation between $B_{avg}$ and $L$. In addition, we find that the models of the bands including cooler emission (131 and 171~\AA) greatly underestimate the extent of the emission in the legs of the longer loops at the peripheries of the active region that are the strongest contributors of the emission in those bands. We conclude that this is because the modeling assumes that all transition region emission is confined to the loop foot points, but in reality the upper transition region of longer loops extends significantly farther into the loop. It is important to consider this aspect of the transition region in future efforts to model EUV emission.

Tables

Table 1. Table 1: 211 and 355 Å best-fit comparisons

$a$	$b$	$q_0$ 211 Å	$q_0$ 335 Å	211 Å quality	355 Å quality
		(${10}^{-3}$ erg cm-3 s-1)
1.5	1	7.1	4.8	Best fit	within 10% of min $\sum {\delta }^2$
0.5	1.5	3.9	6.1	within 10% of min $\sum {\delta }^2$	Best fit

Equations

$$\langle Q\rangle=q_{0}\bigg(\frac{B_{\mbox{\scriptsize avg}}}{B_{0}}\bigg)^{a}\bigg(\frac{L}{L_{0}}\bigg)^{-b}\mbox{erg~cm${}^{-3}$~s${}^{-1}$}$$

$$B_{\text{mbox}avg}\propto L^{-\gamma }$$

$$b=b_0-\gamma a.$$

$$⟨Q⟩\propto L^{-b_0}.$$

Full text

Sun: Sun:corona

††journal: The Astrophysical Journal††software: SolarSoft

show]therese.a.kucera@nasa.gov

show]gnita@njit.edu

show]james.a.klimchuk@nasa.gov

show]gfleishm@njit.edu

I Introduction

Active regions are the chief source of solar emission in X-rays and much of the extreme ultraviolet (EUV), which are, in turn, critical for understanding how the Sun affects Earth’s upper atmosphere. Understanding the source of this emission and how to predict it is also closely connected to the long standing question of how the corona and transition region are heated. There are many different aspects of this problem still to be addressed, but one key test of our understanding is the spatially resolved modeling of the emission from a specific active region in three dimensions.

Active region modelers often take the approach of applying and constraining scaling laws describing heating to a set of model loops. Here, “loop” refers to a magnetic flux tube, not necessarily a bright feature in an image. There is a diffuse component of the corona comprised of observationally indistinguishable loops. There are many theoretically based relationships connecting the heating to quantities like loop geometry, magnetic field, and currents (Mandrini et al., 2000; Cranmer and Winebarger, 2019) to which the results can be compared (Gudiksen and Nordlund, 2002, 2005; Schrijver et al., 2004; Mok et al., 2005, 2008, 2016; Warren and Winebarger, 2006; Lundquist et al., 2008a, b; Ugarte-Urra et al., 2017) for the purpose of bolstering or discarding particular heating models. Of special interest to us is the use of scaling laws based explicitly on episodic “nanoflare” heating models, employed, for instance, by Barnes et al. (2019); Mondal et al. (2025).

Another technique is to heat the corona using self-consistent MHD simulations in which coronal heating is produced by numerical resistivity and viscosity (Gudiksen and Nordlund, 2002, 2005; Bourdin et al., 2013; Rempel, 2017; Warnecke and Peter, 2019). This has the advantage of creating a self-consistent model, sometimes called an ab initio model, which does not rely on assumptions about coronal heating. However, the drawback of this approach is that the resistivity and viscosity in the models are orders of magnitude larger than solar values, and the mechanism of heating may be fundamentally different from the true solar mechanism, e.g., relatively passive Ohmic dissipation versus explosive magnetic reconnection.

In contrast to the works above that have their main (although not exclusive) focus on heating parameters, Plowman (2021, 2023) made fairly simple assumptions about the distribution of emission in loops so as to use the EUV emission to constrain the magnetic field.

In this paper we focus on testing the uses of the GX Simulator model to simulate the emission of an active region observed in EUV in multiple wave bands by the Solar Dynamics Observatory’s Atmospheric Imaging Assembly (SDO/AIA) in order to determine the dependence of emissions on loop length, magnetic field strength and average heating rate. GX Simulator was originally developed for the purpose of modeling flare loops in microwaves and X-rays (Nita et al., 2015), but has since been expanded to model active regions in radio and EUV (Nita et al., 2018; Fleishman et al., 2021b; Nita et al., 2023; Fleishman et al., 2025a; Kuznetsov et al., 2025). GX Simulator uses scaling-law based heating combined with magnetic field extrapolation calculated using a vector magnetogram (see Section II). One advantage to GX Simulator over related approaches is that it computes the emission from every voxel in the volume modeled rather than along particular loops extrapolated upward from a grid of footpoints on the photosphere. Models of the latter type either have large emission-free gaps between loops in the corona or non-physical loop overlap at low altitudes.

In Section II we describe the process used to model AR 12760 with GX Simulator, including the magnetic field, plasma properties, and fitting of the model to the EUV intensity distribution of the active region. Section III describes the fit results, with special focus on the 211 and 171 Å bands. In Section IV we discuss the results including their indication of a relationship between the magnetic field and loop length, the importance of correctly modeling the upper transition region, and ways in which the modeling might be improved. In Section V we summarize our results.

II Modeling Process

II.1 Overview

GX Simulator is described extensively by Nita et al. (2018, 2023) and is available as the gx_simulator package in Solar Soft (Freeland and Handy, 1998). Here we discuss the points most relevant to this paper, starting with a very general overview and adding details in Sections II.2-II.5.

Our goal is to calculate the dependence of the temporally and spatially averaged volumetric heating along a loop, $\langle Q\rangle$, as a function of the mean magnetic field, $B_{\mbox{\scriptsize avg}}$, and loop length, $L$, on the assumption that it can be described by:

[TABLE]

where $B_{0}=100$ G, $L_{0}=10^{9}$ cm and $q_{0}$ is the typical heating rate.

The basic process is set out in Figure 1, and details are presented in the sections below. For a given active region, a 3D magnetic field model is developed via a nonlinear force-free field (NLFFF) extrapolation. The values for $L$ and $B_{\mbox{\scriptsize avg}}$ for each voxel are determined from the field line closest to the center of the voxel. Then, for a given $a$, $b$, and $q_{0}$ it is possible to calculate $\langle Q\rangle$ using Equation 1. The code utilizes a lookup table, to determine the differential emission measures (DEMs) of the corona and transition region for the loop as a function of $\langle Q\rangle$ and $L$. The DEM values are convolved with the instrument response function for a given EUV band to determine the amount of emission in that band for a particular voxel. The emission in each voxel along the line of sight is then integrated to produce a 2D model image. This image is compared to the EUV data to determine the best fit as a function of $a$, $b$, and $q_{0}$.

Files needed to reproduce the modeling, including the input SDO/AIA data, fitting masks, magnetic field models, and GX simulator result files are posted at doi:10.5281/zenodo.16896516.

II.2 EUV Data

We selected for study Active Region (AR) 12760, a small compact active region that crossed the meridian on 2020 April 28. It was chosen because of its compact structure, location near Sun center, lack of sunspots, and low level of activity. It was also the only significant active region on the Sun at the time, making it a good candidate for comparison with full-sun spectral irradiance measurements that, although not a part of this paper, pertains to a longer-term goal of a larger project.

In order to constrain the heating parameters we use data from the 94, 131, 171, 193, 211, and 335 Å bands of Solar Dynamics Observatory’s Atmospheric Imaging Assembly (SDO/AIA; Pesnell et al. (2012); Lemen et al. (2012)). Because we wish to simulate the general, time averaged behavior of the active region rather than its emissions at any one moment, we use an average of the AIA data at a 5 minute cadence from 9:15-13:15 UT (see Figure 2). This time period is slightly before the region crossed the meridian, but was selected to avoid small instances of emerging flux that occurred later in the day and could have affected the magnetic model. Each image was corrected for solar rotation before averaging, aligning it to the region’s location at 11:15 UT.

At any given time, the observed active region corona is composed of a diffuse component and a component of bright loops (on the order of 1 Mm, or a few arcseconds in width in the AIA images). There is more total emission in the diffuse component than the loop component, as loops are typically only modest (10-40%) enhancements over the background (Klimchuk and DeForest, 2020). These bright, thin loops appear and disappear at different locations with lifetimes of several thousand seconds. Over a 4-hour average they tend to combine to form larger structures (with widths on the order of $10\arcsec$), as can be seen in Figure 3, which compares individual images taken at 11:15 UT with the average over 4 hours. A successful GX Simulator model would reproduce both the diffuse component and this medium-scale structure. Note that there is an abundance of theoretical evidence that both components are comprised of ultra-thin, spatially unresolved loops (e.g., Klimchuk et al., 2023; Johnston et al., 2025).

II.3 Magnetic Field Model

To determine the coronal magnetic field we use the Automatic Model Production Pipeline (AMPP, Nita et al., 2023) initiated with the SDO/Helioseismic and Magnetic Imager (HMI; Scherrer et al. (2012)) vector magnetic field map from 11:12 UT.

The magnetic field is calculated using a nonlinear force free field extrapolation as described in Fleishman et al. (2017) and based on the methods of Wiegelmann (2004). Following tests performed by Fleishman et al. (2017) we did not apply any preprocessing to the photospheric boundary condition. The top and side boundary conditions are taken from the initial potential-field extrapolation that is initiated by the $B_{z}$ bottom boundary condition. To minimize the effect of these boundary conditions on the NLFFF solution, we follow Wiegelmann (2004) and employ buffer zones at the top and side boundaries that comprise 10% of the corresponding size of the data cube as detailed in Nita et al. (2023). The C++ source code, along with set of compiled libraries for Windows, Linux, and Mac OS platforms, their calling IDL wrappers, and the GUI application, AMPP, are included in the GX Simulator SSW distribution package. These source codes and platform specific compiled libraries are automatically updated from their independently maintained GitHub development repository111Magnetic-Field_Library. A Python utility based on the same library is also available from GitHub (pyAMPP, Yu et al., 2025). In addition, this NLFFF package may be directly downloaded from a Zenodo digital repository (Stupishin, 2023)222Magnetic Field Library: NLFFF and magnetic lines. The magnetic field strength in this active region is fairly low, making it somewhat difficult to model in detail, as will be discussed in Section III.

We developed two models, each of which is a compromise between resolution and field of view (FOV). Both models were centered at 2020 April 28 11:12 UT at $x=-75\arcsec$, $y=-50\arcsec$, and the point of view derives from those coordinates. The AIA data were adjusted to match that time. Thus, the view of the active region we are modeling is quite close to the actual view of AIA and not a simple top-down point of view. The models consisted of cubes with 180 cells in the $x$ and $y$ directions and 200 cells in the $z$ direction. Magnetic Model 1 (Figure 4a) has a cell size of 1400 km so as to better model the core field lines but some of the longer loops extending from the active region periphery are “open” and thus not included in our modeling of the emission (see Section II.4 below). Magnetic Model 2 (Figure 4b), conversely, covers an area twice as large while increasing the cell size to 2800 km to remain within the memory constraints of the computer being used. Model 2 includes the long loops, but produces a more simplified model of the core of the active region. Sample field lines from both models are shown in Figure 4 in which field lines that close within the model box are shown in green, while ones that extend outside the box are in yellow.

II.4 Loop Plasma

The active region volume is divided into voxels and a value for $B_{\mbox{\scriptsize avg}}$ and $L$ is assigned to each voxel using the NLFFF model. Each voxel is characterized using the field line passing closest to its center, providing a value for $B_{\mbox{\scriptsize avg}}$ and $L$. The plasma properties in the loop are determined using a look-up table as a function of $L$ and $\langle Q\rangle$, determined using the Enthalpy-Based Thermal Evolution of Loops (EBTEL, Klimchuk et al. (2008); Cargill et al. (2012a, b)) model with a specific distribution of heating events. EBTEL is a “0D” model that computes the evolution of the coronal-averaged temperature, pressure, and density in a loop without spatial resolution. The instantaneous coronal DEM is assumed to be distributed uniformly over a temperature interval centered on the average temperature and having a width equal to twice the difference between the apex and average temperature. Because the time-averaged coronal DEM is determined by the behavior of many unresolved loops over time, it is reflective primarily of longer term loop evolution and not the instantaneous temperature distribution. Thus, we use the time average of long-duration EBTEL simulations that encompass many nanoflare energy releases. The time-averaged DEM is spread uniformly along the entire coronal portion of the loop. The EBTEL runs currently used to calculate GX Simulator input tables assume no loop expansion.

For the transition region EBTEL computes a separate lookup table with the detailed DEM distribution at each time in the simulation based on energy balance and averages them to produce DEM curves for the $\langle Q\rangle$ and $L$ grid. EBTEL takes the temperature boundary between the corona and transition region to be the value where thermal conduction switches from a cooling term above to a heating term below. This temperature is roughly 60% of the apex temperature. The lower and middle transition region by temperature are confined to low altitudes, but 1D loop modeling indicates that the upper transition region extends an appreciable distance from the loop footpoint. This plays an important role in the model images, as we discuss in Section IV.2.

In order to produce a lookup table of coronal and transition-region DEMs as a function of $\langle Q\rangle$ and $L$, a large number of EBTEL runs are done. The $60(Q)\times 49(L)$ grid provides logarithmically spaced values in the range $L/2=1$ to $6.3\times 10^{4}$ Mm. For each $L/2$ value there is a different range of $\langle Q\rangle$  with $\langle Q\rangle$ ranging from 0 to 794 erg cm*-3* s*-1* for $L/2=1$ and from $1.5\times 10^{-8}$ to 0.0126 for $L/2=6.3\times 10^{4}$ Mm. As discussed above, for each set of parameters a long run is averaged in order to represent a loop as a group of smaller strands heated randomly by an ensemble distribution of heating events. The currently available EBTEL lookup tables in GX Simulator are described in Nita et al. (2023). For this work we use ebtel_scale=1_alpha=-1.sav, which assumes a random distribution of nanoflares with an average time between nanoflares equal to the loop cooling time (scale = 1) and a power law distribution in the delay between events with index $\alpha=-1$. The energy of the events is proportional to the time between them so that the energy also has a power-law distribution with the same index.

GX Simulator calculates the transition region and corona emission independently and then sums them. For each voxel the DEMs are determined from the EBTEL table using a bilinear interpolation. GX Simulator places all the transition-region emission in a single layer at the bottom of the model cube. We assume a loop with uniform cross-sectional area. This can affect the emission and the ratio between coronal and transition-region emission in particular (Cargill et al., 2022). The transition region is relatively fainter in expanding loops, but how adding expansion would affect our model is a complex topic, as discussed in Section IV.3.

Currently the way GX Simulator can model magnetic field lines that leave the model box is to either assume they contain plasma in hydrostatic equilibrium based on an assumed temperature and base pressure or to exclude them entirely. We have elected to leave out emission from these “open” field lines (open in the sense that they leave the box, not necessarily open to the solar wind). In addition, the emission was also set to zero for voxels centered on a small number of very long field lines not within the boundaries of the $\langle Q\rangle$,$L$ lookup table and in cases in which the fitting algorithm was attempting to use unrealistically low values of $q_{0}$ in an effort to find a good fit. The boundaries of the EBTEL table do not affect the results greatly in the parameter ranges with good fits to the data, although we discuss a possible issue related to open field lines in Sections III.2 and IV.2.

To calculate the emission in each voxel the DEM was convolved with the temperature response function currently available in the routine provided by the AIA team, aia_get_response.pro, which uses the Feldman et al. (1992) coronal abundances.

II.5 Fitting of the Model

The fit is done by comparing on a pixel by pixel basis the SDO/AIA image with the intensity image calculated by GX Simulator. We fit the $240\times 240\arcsec$ region show in Figure 2. Most of the calculated model images had a resolution of 2″ (resulting in a $120\times 120$ pixel field of view), and were convolved with a 1.2″ point spread function, but, as will be discussed in Section III.1, we also did fits at 8″ resolution. The AIA data were rescaled to match, with the intensities rescaled to account for the change in the pixel area.

GX Simulator does any needed alignment of the AIA maps with the magnetic field model (see Section II.3). Unlike the radio modeling done with GX Simulator (e.g., Fleishman et al., 2021b), the relative pointings of the data from the two SDO instruments is well understood, so no cross-correlation based alignment was done.

When calculating the model image we constrain the transition-region emission to areas where $|B_{z}|\geq 200$ G in the magnetogram, as in Nita et al. (2018). The resulting transition-region mask is shown in Figure 5. As stated in the pervious section, all transition-region emission is restricted to the lowest layer of the model data cube. With this additional mask, we also confine the transition-region emission to only those areas where the photospheric field is strong and mostly vertical (Klimchuk, 1987; Nita et al., 2018). The mask accounts approximately for the rapid expansion of the field in the upper photosphere and chromosphere, which occurs because of the transition from very large to very small plasma beta. This transition occurs over a thin but finite layer. In contrast, the magnetic field extrapolation assumes an infinitely sharp boundary, which leads to unrealistically high magnetic flux in areas with lower magnetic field, and thus strong transition-region emission everywhere, contrary to observations. Thus we suppress the model transition-region emission in these low $B_{z}$ regions. The physical motivation for the mask is discussed further by Nita et al. (2018) and in Appendix A.

To emphasize the bright regions of the map, while still incorporating constrains from the the dark areas surrounding the active region, we assigned a constant uncertainty to all pixels in a given wave band, so that our model-to-data comparison metric is the sum of the squared residuals between the model and data, $\sum\delta^{2}=\sum(I_{\mbox{\scriptsize obs}}-I_{\mbox{\scriptsize mod}})^{2}$, where $I_{\mbox{\scriptsize obs}}$ and $I_{\mbox{\scriptsize mod}}$ are the observed and modeled intensities respectively. This means that the brighter areas of the image are weighted more heavily, but the emission is still constrained to low values outside the active region. This metric of success allows us to compare the quality of the fit within a given wave band, but it does not allow us to provide quantitative goodness-of-fit contours in our parameter space, or compare the goodness of fits between wave bands.

For each ($a,b$) pair, we employed the Coronal Heating Modeling Pipeline utility (CHMP, Nita et al., 2023), which is part of the GX Simulator package, to compute maps for a range of values of $q_{0}$ converging to a minimum $\sum\delta^{2}$ value, thus determining the best fit for that set of ($a,b$) values. For areas of the ($a,b$) grid with the lower $\sum\delta^{2}$ values the $q_{0}$ vs. $\sum\delta^{2}$ curves were smooth with clear minima. We considered a range of values for $a$ from $-1$ to 3 and $b$ from $-1$ to 2.5, in increments of 0.5.

We did two general types of data-model comparisons. Most included both the modeled corona and transition-region emission. As mentioned above, the modeled transition-region emission is limited to the bottom layer of the model cube in regions with photospheric $|B_{z}|\geq 200$ G as shown in Figure 5. All points in the array were used in the data-model comparison except for a one pixel wide band around the edge of the array to allow for shifts of the mask when interpolating to align with the AIA data. We also did fits modeling coronal emission only. GX Simulator calculates the corona and transition region separately (as discussed in Section II.4) and the transition-region emission can be turned off entirely. Thus in the corona-only fit the model transition-region emission was set to zero and the $\sum\delta^{2}$ calculation excluded the points designated as transition region (Figure 5) as well as the edge pixels.

Below we focus on two AIA bands, 211 Å, which shows chiefly coronal emission, and 171 Å, which shows a combination of coronal and transition-region emission.

III Fit Results

III.1 Fit to the 211 Å band

We first present the fits to the 211 Å band, which is dominated by Fe XIV with a peak emission temperature $\log T\mbox{(K)}\approx 6.3$ (Boerner et al., 2012). In Figure 6 we show a grid of thumbnails showing the best $q_{0}$ solutions for each of the ($a,b$) pairs considered for Model 1, including emission contributed by both the model corona and transition region. The best fit corresponds to $a=1$, $b=1$, $q_{0}=4.4\times 10^{-3}$ erg cm*-3* s*-1*. This model image has a title and border in green. Maps with $\sum\delta^{2}$ within 10% and 20% of the minimum have titles and borders in orange and red, respectively. The result shows a linear area with similar fit quality. We discuss this in more detail in Section IV.1, which includes a plot of $\sum\delta^{2}$.

However, the magnetic field in AR 12760 is weak, making it hard to model accurately because of observational uncertainties of the magnetic field inference. To try to compensate for that, we also compare the data and model maps at a lower resolution of 8″ to reduce the effect of differences on a pixel-by-pixel basis (see Figure 18 in the Appendix for the resulting thumbnail grid). With the comparison at this resolution we find the best fit at $a=1.5$, $b=1.0$, $q_{0}=7.1\times 10^{-3}$ erg cm*-3* s*-1*. We note that this best fit ($a,b$) pair corresponds to a pair that is within 10% of the minimum $\sum\delta^{2}$ for the higher-resolution fit (Figure 6). In general we would say that our uncertainty for $a$ and $b$ is around 0.5 because that is the step size in our ($a,b$) grid. However, the situation is rather more complicated due to the complex ($a$,$b$,$q_{0}$) parameter space, as will be discussed in Section IV.1.

The best fits for the 2″ and 8″ comparisons are shown in Figure 7 and Figure 8 respectively, along with the AIA data at the same resolution, and the residual, $\delta=I_{obs}-I_{mod}$. The plots also show the model coronal and transition-region contributions to the emission as separate maps. The coronal and transition-region emissions are calculated separately and summed by GX Simulator, and seperate transition-region and coronal maps can be synthesized in addition to the total emission. The model maps look reasonable on a general level, with most of the emission coming from the active region core. However, the detailed correspondence of the bright features is not as good as we would like. A significant part of this can be attributed to inaccuracies in the magnetic field extrapolation (Fleishman et al., 2019). Even if the hydrodynamics along field lines is modeled correctly, the spatial location of the emission will be off. We mitigate this in part by reducing the resolution to 8″, although differences are still apparent. The transition-region emission is a significant contributor to the emission near loop footpoints, but it is not the main source of the 211 Å emission.

The transition-region emission in the model seems especially bright and localized compared to the observations. We address the localization later. The brightness is a problem that has historically afflicted many models, including detailed 1D hydrodynamics models, i.e., “loop” models. As mentioned in Section II.4, the discrepancy is reduced in expanding loops, but the excess brightness remains an unsolved problem. For this reason, we also considered a fit to the corona only. The resulting best-fit model is shown in Figure 9. For this fit we turned off the transition-region emission and used only pixels not designated as including transition region (Figure 5) to calculate $\sum\delta^{2}$. In this case we found the best fit to be $a=2$, $b=0.5$, $q_{0}=5.8\times 10^{-3}$ erg cm*-3* s*-1*. The intensity scale is the same as in Figure 7. The model looks similar to the best-fit model for the combined corona and transition-region fit, although is a bit brighter than the corona alone in that fit. This ($a,b$) pair also corresponds to a pair that is within 10% of the minimum $\sum\delta^{2}$ for the corona and transition-region fit shown in Figure 6.

III.2 Fit to the 171 Å band

The fit results for the 171 Å band, which has emission produced principally by Fe IX with a peak at $\log T$(K)$\approx 5.8$, are rather different from those of the 211 Å band. Figure 10 shows the ($a,b$) grid of fits using Magnetic Model 1 with both the model corona and transition region included. The best fit is $a=0.5$, $b=0.0$, $q_{0}=2.0\times 10^{-4}$ erg cm*-3* s*-1*. These are substantially different values for $a$ and $b$ than those found for 211 Å. The best-fit 211 Å model for these $a$ and $b$ values (Figure 6) looks like a terrible fit to 211 Å, and, conversely, the best ($a,b$) values for 211 Å do not appear to lead to good fits for 171 Å, although there seem to be significant problems with all the 171 Å fits.

The difference appears to be connected to longer loops on the periphery of the active region. Figure 11 shows the details of the best-fit model, again with separate images showing the corona and transition-region components. As would be expected, the localized model transition-region emission is more evident in the 171 Å band than in 211 Å. However, in the “best fit” this seems to be exaggerated; almost all the emission is coming from the designated transition-region foot points with very little emission in the active region core. Also, although the 171 Å image is dominated by the bright areas in the long, peripheral loops (circled in red in Figure 11a), these features are entirely missing from the model images. This is not just true of the best fit, but of all the fits shown in Figure 10. Feature A appears to be the base of a fan-loop, a large scale structure that we think closes in the corona. Feature B seems to form a more (although not entirely) complete arc, although its exact size and altitude are not clear.

One possible explanation for the total lack of emission in the peripheral loops is that the relevant loops did not close within the modeled cube, and thus were not included in the emission (see Section II.3). This was the reason we developed Magnetic Model 2. Figure 12 shows magnetic field lines from both Magnetic Model 1 and Model 2 superimposed over the 171 Å reference image. It can be seen that in Magnetic Model 1 some of the longer field lines extending from the periphery of the active region do not close down again inside of the box (these are shown in yellow), and so are being excluded from our model emission. In contrast, in Magnetic Model 2, the lines from these areas do close down inside the NLFFF model cube, so we would expect them to be included in the emission.

However, this proposed explanation does not seem to be correct. In Figure 13 we show the best fit result for model maps calculated with Magnetic Model 2. The best fit is $a=1.0$, $b=0.0$, $q_{0}=4.7\times 10^{-4}$ rg cm*-3* s*-1*. The full ($a,b$) grid is shown in Appendix B, in Figure 20. The longer field lines in the periphery of the active region are still not present in the emission, although the footpoints of the loops in the transition region are very bright.

Thus, the use of the lower-resolution but larger-volume magnetic field model results in a less accurate distribution of emission in the core of the active region, but it does not result in the appearance of the longer loops, so the problem does not seem to be that the loops are being left out of the magnetic field model. As we will discuss in Section IV.2, we think the problem is in how the model transition-region emission is distributed along the loop.

III.3 Fits to other AIA bands

The other AIA bands we considered, 94, 131, 193, and 335 Å, show a temperature dependent pattern similar to that seen in 171 and 211 Å(see Figure 14). For the 211 and 335 Å bands, which are dominated by emission at $\log T\mbox{(K)}\gtrsim 6.3$ (Boerner et al., 2012), the fits look similar. Although they differ from the original AIA image in detail, they show similar amount of emission in the active region core, which dominates the emission in both data and model. The thumbnail plots for the comparisons of the model to the 335 Å data binned to 8″ are shown in Figure 21, and more discussion of the 211 Å and 335 Å results is in Section IV.1.

The 131 and 171 Å bands, on the other hand, which in the absence of flaring have peak emissions at $\log T\mbox{(K)}\lesssim 5.8$, both show bright peripheral loops that are missing from the models except at the loop bases. The 94 Å band, in this case dominated by its response peak at $\log T\mbox{(K)}\approx 6.05$, and 193 Å, with a response peaking at $\log T\mbox{(K)}\approx 6.1$, are in between, showing the peripheral loops that are, again, missing from the model fits, but the peripheral loops are not as bright as for the bands dominated by cooler emission.

IV Discussion

IV.1 Heating Parameters

There are many measurements of $a$, $b$ and related parameters in the literature. There is a summary of measured and theoretical values (the latter from Mandrini et al., 2000) presented by Ishigami et al. (2024), who include a summary plot of values that correspond to our $a$ and $1-b$. They vary quite a bit and our values of $b$ in particular are lower than those measured by others. It is important to note that these measured values - obtained by comparing models with observations as we do here - are for active regions in some studies and the global corona in other studies. Differences are not surprising because the dominant heating mechanism may be different in active regions and the quiet Sun.

Our most reliable fits are the low-resolution (8″) fits to the 211 and 335 Å data, because the lower resolution smooths out some of the inaccuracies in the magnetic field model. Figure 15 shows plots of our model-to-data comparison metric, $\sum\delta^{2}$, as a function of $a$ and $b$, where each fit has been optimized for $q_{0}$ for that fitted $a,b$ pair (see Figures 18 and 21). The best fit for 211 Å, shown by the asterisk, is $a=1.5$, $b=1$, and $q_{0}=7\times 10^{-3}$ erg cm*-3* s*-1*. If we consider the $q_{0}$ values for fits within 20% of the minimum $\sum\delta^{2}$ value then we get $q_{0}$ in the range $3\times 10^{-3}$ and $11\times 10^{-3}$ erg cm*-3* s*-1*. For 335 Å it is $a=0.5$, $b=1.5$, and $q_{0}=4\times 10^{-3}$ erg cm*-3* s*-1*. The fits within 20% of the minimum include $q_{0}$ values from $3\times 10^{-3}$ to $7\times 10^{-3}$ erg cm*-3* s*-1*.

A recurring feature in these plots of model-to-data comparison metrics is the clear diagonal feature ranging from low $a$, high $b$ values on the top-left to high $a$, low $b$ values towards the bottom-right. This feature is consistent with other studies with EUV (Ugarte-Urra et al., 2019) and microwave data (Nita et al., 2023; Fleishman et al., 2025b; Kuznetsov et al., 2025). For the 8″ resolution 211 Å fit we found that this feature roughly corresponds to $b\approx 1.6-0.6a$ and for the 335 Å fits we get $b\approx 1.7-0.6a$. That our best fits for each wave band do not fall exactly on these lines is indicative of the coarseness of our ($a,b$) grid and uncertainties in our fitting of diagonal features. More specifically, if we assume a 0.25 uncertainty (half a ($a,b$) grid step) on both $a$ and $b$ we have $b$ intercepts of $1.61\pm 0.12$ and $1.70\pm 0.13$ and slopes of $-0.56\pm 0.07$ and $-0.59\pm 0.08$ for 211 Å and 335 Å, respectively, so the results for the two bands are the same within uncertainties.

We note that although the best $(a,b)$ values are different for the two wave bands, they both fall close to these best-fit lines. Table 1 compares the best model-to-data comparisons of the two wave bands. We see that for both wave bands the best ($a,b$) values fall within 10% of the minimum $\sum\delta^{2}$ for ($a,b$) values of the other wave band. This gives an idea, albeit somewhat qualitative, of our confidence in our fit values. The model-to-data comparisons for the other bands are sufficiently poor (see Sections III.2 and III.3) that we do not include them in this analysis.

This linear feature in the model-to-data comparison metric parameter space indicates that $L$ and $B_{\mbox{\scriptsize avg}}$ are correlated. This is, of course, no surprise. It is well-known that the coronal magnetic field is on average stronger in the short loops of the active region core than in the long loops of the periphery. We can show that this feature is reasonably consistent with simple relationship between $B_{\mbox{\scriptsize avg}}$ and $L$.

As stated in Equation 1, we are assuming $\langle Q\rangle\propto B_{\mbox{\scriptsize avg}}^{a}L^{-b}$. So if

[TABLE]

(Klimchuk and Porter, 1995) then we have $\langle Q\rangle\propto L^{-\gamma a-b}$. Model fits should be similarly good when $\gamma a+b=$ constant; call it $b_{0}$. Good fits should therefore fall along a straight line in a plot of goodness of fit in the $(a,b)$ plane:

[TABLE]

Thus, the slope of the best fit line in the goodness-of-fit parameter space gives us the relationship between $B$ and $L$ (Equation 2), and the relationship for $\langle Q\rangle$ and $L$ is dependent on the y-intercept of the line:

[TABLE]

Our results shown in Figure 15 suggest that $b_{0}\approx 1.6$ or 1.7 and $\gamma\approx 0.6$. As stated above, other analyses also show that similarly good fits fall along a straight line. Fleishman et al. (2025b) using GX Simulator to model optically thick gyroresonant microwave emission observed by the Siberian Solar Radio Telescope (SSRT) and Nobeyama Radioheliograph (NoRH) found $b=2.4-1.0a$ at 5.7 GHz and $b=2.2-0.7a$ at 17 GHz and a best solution at ($a,b$) = (0.6, 1.8). Kuznetsov et al. (2025) analyzed eight active regions in 6, 8, and 10 GHz with the Siberian Radioheliograph and found that $b_{0}$ was mostly between 1.6 and 3.2, with the distribution peaking at 2.5, and $\gamma$ ranged from 0.2 to 1.2, with a most probable value of 0.55. It is also interesting that Porter and Klimchuk (1995) concluded that $\langle Q\rangle\propto L^{-2}$ from a study of soft X-ray loops observed by Yohkoh, which can be compared to the $b_{0}$ values derived from GX Simulator.

Mandrini et al. (2000) examined magnetic field extrapolation models and found $\gamma=0.88$ using Equation 2 directly. Fleishman et al. (2021a) also found a $B_{\mbox{\scriptsize avg}}$-$L$ dependence based on a magnetic field model, although their fitting model was more sophisticated than a simple power law, i.e. $\log^{2}B_{\mbox{\scriptsize avg}}=48.61-1.84\log L$.

With this in mind, we fit the $B_{\mbox{\scriptsize avg}}$ vs $L$ dependence in our Magnetic Model 1, as shown in Figure 16, using a least-squares fit in log-log space. We derived $\gamma=0.8$. The plot also shows a line with a slope corresponding to $\gamma=0.6$ as suggested by the GX Simulator modeling. The value derived directly from the magnetic model is slightly higher than the value derived from $\sum\delta^{2}$ plots, but it is clear from the distribution of points in Fig. 16 that the distribution of the field lines in the model is more complex that a simple power-law, with noticeable substructure. We expect Equation 2 to be a general trend rather than an equation applying precisely to all parts of the active region. These results are, however, consistent with the idea that the diagonal feature in GX Simulator model-to-data comparison metrics is due to a dependence of $B_{\mbox{\scriptsize avg}}$ on $L$.

We should note that Fleishman et al. (2025b) and Kuznetsov et al. (2025) are looking at different active regions, which in addition to a difference in wave bands and a different model-to-data comparison metric could lead to differences in the properties of the zones of best fit. Kuznetsov et al. (2025) in particular find a range of values over the many active regions they modeled.

However, the broad agreement of these results is very encouraging. This validates the power-law dependencies assumed in Equations 1 and 2. Had similarly good fits not fallen along a straight line, one or more of the dependencies would be invalidated.

IV.2 Transition-Region Emission

An important issue in our modeling effort is that the emission from long features at the periphery of the active region is missing from the models. For this active region these are dominant features in the cooler bands, 131, and 171 Å, and significant contributors to other bands with peak contributions below about 1.6 MK (94 and 193 Å). We think this is due to an unrealistic distribution of the transition-region emission in these long loops. In this work and other modeling studies, the transition region is defined to be the locations in the loops where thermal conduction serves as a heating term rather than a cooling term (Vesecky et al., 1979). As discussed in Klimchuk et al. (2008), this is a precise, physically motivated definition that avoids ambiguities that arise from vague observational definitions based on temperature or altitude. As physically defined, the transition region extends from the top of the chromosphere to about 10% of a loop half-length ($0.05L$) (see Klimchuk et al., 2008; Cargill et al., 2022). An example of this is shown in Figure 17, which depicts the temperature, electron density, and normalized 171 Å intensity as a function of position of a 1D hydrodynamic model of a loop in equilibrium with steady uniform heating. The model was computed with the ARGOS code (Antiochos et al., 1999) and is broadly representative of static equilibrium loops with a wide range of lengths and heating rates. The vertical line shows the top of the transition region at $0.054L$ from the loop footpoints. This same spatial fraction applies approximately to other loops, but the temperature of the top of the transition region varies, always about 60% of the apex temperature.

For relatively short loops it is sufficient to assume that the transition-region emission would be confined to a limited range approximated by a single voxel along the loop in our model cube. In our higher-resolution Model 1 each voxel is 1.4 Mm on a side. That would be 5% of a loop with $L=28$ Mm, within a factor of two or so of loop length for the core of an active region. However, based on Magnetic Model 2, the loops at the periphery of the active region likely have $L=200-700$ Mm, implying transition regions on the order of 10-50″ in the AIA images, much larger than the pixels immediately around the loop base and roughly consistent with the large cool loops we see in the data. Thus, in order to get a better fit to EUV images showing extended emission in loop legs, it will be necessary to account for the finite spatial extent of the transition region. Note that it is only the upper transition region, the upper third as determined by temperature, where this is important. The temperature scale lengths in the lower and even middle transition region are extremely short. The amount of transition-region emission predicted by EBTEL is reasonable (although see caveats about general transition-region modeling below), but because all of it is assigned to the footpoints by GX Simulator, the footpoints are too bright and the extended legs do not appear in model images.

This explanation would explain Feature A, the fan-loop feature in Figure 11. It is not as clear if it accounts for Feature B in that plot, in which the emission appears to extend farther than $0.05L$ along a magnetic tube. We think it unlikely Feature A extends into the solar wind because we would expect the emission along open magnetic field lines to be fainter, as they are in coronal holes. It is possible the actual feature is on closed magnetic field lines that extend outside the volume we are modeling. However, that would still leave the question of why we would see emission on such a structure but not along the model field lines that are in Model 2.

Our group is working on multi-strand MHD simulations that will make it possible to determine DEMs for both corona and transition region as a function of location along the loop, which, if implemented in GX Simulator, should take care of this problem (Johnston et al., 2025; Sow Mondal et al., 2025).

As discussed in Fleishman et al. (2025b), one possible way to obtain more information about the transition region is to combine analysis of microwave emission, which is sensitive to only the corona, with EUV emission, which includes both transition region and coronal components. This will help to separate out the emission from the two sources. Another approach would be to consider active regions from multiple points of view using AIA combined with either STEREO/EUVI (Solar TErrestrial RElations Observatory/Extreme UltraViolet Imager; Howard et al., 2008) or Solar Orbiter/EUI (Extreme Ultraviolet Imager; Rochus et al., 2020) near quadrature, so that lower-altitude transition-region emission could be separated out. However, for both of these approaches it would be important to consider the extended transition-region emission from the long loops discussed in this paper.

IV.3 Loop Expansion

We note an inconsistency in our models. Loops expand with height in the magnetic field extrapolation of GX Simulator, whereas the EBTEL simulations assume a constant loop cross section. Often, multiple voxels in a cut across the magnetic field map magnetically to a single transition-region voxel. Each of the coronal voxels is assigned a coronal DEM that corresponds to a companion transition-region DEM from the EBTEL simulation. The transition-region DEM appears one time in the active region model, but the coronal DEM appears multiple times. Consequently, these geometrical considerations suggest that the volume integrated coronal emission in the model is overestimated relative to the transition-region emission by a factor approximately equal to the area expansion factor of the loop.

However, this geometric effect is at least partially offset by the different thermodynamic properties of expanding and non-expanding loops. Cargill et al. (2022, Fig. 2) found that the ratio of coronal to transition-region emission per unit area is approximately independent of the expansion factor for a fixed heating rate. Therefore, for a given amount of transition-region emission, the volume integrated coronal emission should scale roughly with the expansion factor. This means that the thermodynamic effects missing from the EBTEL models (because they assume a constant cross section) fortuitously offset the geometric effects not accounted for in GX Simulator. In this sense, the current modeling is approximately correct. These effects must, of course, be investigated and verified more rigorously in future work. We note that EBTEL has the capability of simulating expanding loops, and this feature may need to be implemented in next-generation models. This is a significant undertaking as it adds a third dimension to the model lookup table.

We note also that Cargill et al. (2022) find that expansion can affect the ratio of the length of the transition region to $L$, so incorporating loop expansion may also affect the treatment of the location of the transition region in the model (see discussion in Section IV.2 above.)

IV.4 Further Considerations

In addition to the treatment of the transition region and loop expansion there are other areas in which our modeling could be improved. As stated above, the details of the core magnetic field lines do not reflect the data well. This is likely because this active region had a rather weak magnetic field. It was selected because we wished to have a very quiet active region so as to avoid solar flares in the time-average images. Also, we wished to avoid sunspots because it is highly unlikely that sunspots obey the same heating scaling law that applies throughout most of the active region. For one thing, loops rooted in sunspot umbrae have much different emission properties from other active region loops; they are cooler and fainter (Tiwari et al., 2017). For another, many coronal heating models are based on the small-scale stressing of the field by photospheric convection, and convection is suppressed in sunspots. We therefore preferred to avoid this complication in our initial study, but of course sunspots must be included in an eventual complete model. A more typical active region, with stronger fields, may be easier to model magnetically. It would optimally make sense to have an iterative fit in which both the plasma and magnetic field parameters could be varied (see Fleishman et al. (2021a) for a case where fits to radio data were done using multiple magnetic models), but that is beyond the scope of this paper and our readily available computing power.

Another area to consider is the assumptions going into the heating calculations in the lookup tables. Fleishman et al. (2025b) and Kuznetsov et al. (2025) considered the various lookup tables currently available in GX Simulator and, in analyses utilizing microwaves, concluded that the best fits for the diffuse coronal component come from assuming heating parameters based on higher frequency nanoflares than we used for this study or even constant heating. So, ultimately, heating with various nanoflare distributions should be considered for the EUV as well. Furthermore, the heating frequency may vary across the active region (Mondal et al., 2025), which should be incorporated into future models. Power-law forms for the nanoflare energy distribution are commonly assumed in the modeling community (e.g., Pauluhn and Solanki, 2007; Hood et al., 2016; Knizhnik et al., 2018) and are motivated by observed distributions (see review by Reale, 2014). However, the observed distributions are typically based on microflares and full-sized flares, and it is not obvious that the distributions can be extrapolated to smaller energies. Furthermore, the measured events are volume-integrated, whereas the nanoflares inputs to the EBTEL models are energy per unit cross-sectional area. Our multi-strand MHD simulations indicate that the energy per unit area actually obeys a log-normal distribution (Sow Mondal et al., 2025). The MHD simulations also indicate that the delay and energy of the nanoflares are uncorrelated (Sow Mondal et al., 2026), contrary to the assumption in the current GX Simulator lookup tables. These properties must be included when new lookup tables are computed.

Although we concluded that it would not have helped in our case, it would be good to develop a way to estimate the loop lengths of field lines that exit the model box. At the moment the options are to either ignore the plasma in such lines or assume it is hydrostatic, neither of which is consistent with the modeling of the loops that close in the box. The calculation and manipulation of large volume NLFFF models needed to include more loops within the modeled box uses a lot of computer resources, so a way of estimating $L$ and $B_{\mbox{\scriptsize avg}}$ for large loops exiting smaller boxes would be helpful.

V Summary and Conclusions

We used GX Simulator to model the EUV flux from a small, quiescent active region averaged over a 4 hour period in order to determine the dependence of the volumetric heating, $\langle Q\rangle$, on the average magnetic field strength along the loop, $B_{\mbox{\scriptsize avg}}$, and loop length, $L$, assuming $\langle Q\rangle\propto B_{\mbox{\scriptsize avg}}^{a}L^{-b}$. We found that for the AIA bands dominated by emission with $\log T\mbox{(K)}\gtrsim 6.3$ (211 and 335 Å) the larger spatial scale intensity distribution was reasonably well fit, although the detailed structure was different due to difficulties in the magnetic modeling of an active region with relatively weak field. Using a relatively low-resolution version of the model (8″) we found the best fit for 211 Å has $\langle Q\rangle\approx 7\times 10^{-3}(B_{\mbox{\scriptsize avg}}/{100\mbox{G}})^{1.5}(L/{10^{9}\mbox{cm}})^{-1}$ erg cm*-3* s*-1*, and for 335 Å the best fit is $\langle Q\rangle\approx 4\times 10^{-3}(B_{\mbox{\scriptsize avg}}/{100\mbox{G}})^{0.5}(L/{10^{9}\mbox{cm}})^{-1.5}$ erg cm*-3* s*-1*. However, for all fits there is a range of parameter values that give qualitatively reasonable fits lying roughly in a line in ($a,b$) space, indicative of a dependence between $L$ and $B_{\mbox{\scriptsize avg}}$. This line was approximately defined by $b\approx 1.6-0.6a$. This relationship is consistent with a $B_{\mbox{\scriptsize avg}}$ dependence on $L$.

We also found that for EUV lines showing a cooler (but still optically thin) plasma the models did not replicate the emission from loops along the edge of the active region, which dominates the emission in the 171 and 131 Å bands and also contributes significantly to the 94 and 193 Å bands. We conclude that this emission is due to relatively high-altitude transition-region emission in these long loops. EBTEL defines the transition region to be the section of the loop in which thermal conduction serves to heat the plasma rather than cool it (as occurs in the coronal section), and as such the transition region generally fills a length about 5% of the full loop length from the footpoint. GX Simulator and the EBTEL model informing it assume that all transition-region emission comes from one voxel at the footpoint, which is not sufficient to model the spatial distribution of the emission in these AIA bands. Thus we conclude that it is important in modeling active regions to consider the spatial extent of the upper transition region as a function of temperature. The cooler part of the transition region is confined to low altitudes close to the footpoints.

T.A.K. and J.A.K. were supported by the GSFC Heliophysics Internal Scientist Funding Model (competitive work package program). This work is also part of the effort of the Committee on Space Research International Space Weather Action Team (COSPAR-ISWAT) S2-06: Origins of the Spectral Irradiance and its Intermediate Timescale Variability. GX Simulator is a package available as part of Solar Soft (Freeland and Handy, 1998) and also on github at https://github.com/Gelu-Nita/GX_SIMULATOR. G.D.F. and G.M.N. were partly supported by the NSF grants RISE-2324724 and AST-2206424, NASA grant 80NSSC23K0090 to the New Jersey Institute of Technology and by CSTR/NJIT. The authors gratefully acknowledge the use of data from SDO/AIA and HMI. Thanks to the referee for thoughtful suggestions that have significantly enhanced the paper. SDO (AIA, HMI)

Appendix A Transition-Region Mask

We present here a short physical justification for the transition-region mask (see Section II.5), which is expanded upon in more detail by Nita et al. (2018) in their Section 4.1 and Appendix A.

For a given loop length and heating rate, the brightness density of the transition-region emission at the solar surface decreases as the loop leg becomes more inclined to vertical, as discussed after Equation 5 of Nita et al. (2018). One way to understand this is as follows. The radiative energy losses from the transition region per unit cross-sectional area depend on the loop length and heating rate and are independent of the base inclination. By definition, the cross section is perpendicular to the loop axis. It is different from the footprint of the loop on the solar surface, which is larger than the cross section by an amount that increases with the tilt from vertical. Because the fixed radiative losses are spread out over a larger area, the brightness density of the footprint decreases with increasing tilt.

Photospheric magnetograms of active regions reveal regions of strong vertical magnetic field surrounded by regions of much weaker vertical field (e.g., Klimchuk, 1987). This gives rise to a magnetic canopy, where the strong vertical fields bend over sharply to become nearly horizontal in the chromosphere and transition region above the weak fields. From the discussion above, the transition region is observed to be much fainter above weak fields than strong fields.

Canopies do not arise in extrapolated magnetic field models because the transition from large to small plasma beta happens discontinuously at the lower boundary of the model. The field is much less inclined above regions of photospheric weak field than on the real Sun and therefore the predicted transition-region emission is much brighter. To correct this, we impose a transition-region mask. We set the transition-region emission to zero wherever the photospheric field is weaker than a critical value, which we take to be 200 G.

Appendix B AIA Model Grids

This appendix contains more examples of calculated models that give the lowest $\sum\delta^{2}$ values as a function of $q_{0}$ for each of a grid of ($a,b$) values, Figures 18-21.