Anomalous Emission from Li- and Na-like Ions in the Corona Heated via Alfvén Waves

The solar ultraviolet intensities of spectral lines originating from Li- and Na-like ions have been observed to surpass the expectations derived from plasmas with coronal approximation. The violation of the coronal approximation can be partially attributed to nonequilibrium ionization (NEI) due to dynamic processes occurring in the vicinity of the transition region. To investigate the impact of these dynamics in the Alfvén wave-heated coronal loop, a set of equations governing NEI for multiple ion species was solved numerically in conjunction with 1.5-dimensional magnetohydrodynamic equations. Following the injection of Alfvén waves from the photosphere, the system undergoes a time evolution characterized by phases of evaporation, condensation, and quasi-steady states. During the evaporation phase, the ionization fractions of Li- and Na-like ions were observed to increase when compared to the fractions in ionization equilibrium, which led to an enhancement in the intensity of up to 1.6. This over-fractionation of Li- and Na-like ions was found to be induced by the evaporation process. While collisions between shocks and the transition region temporarily led to deviations from ionization equilibrium, on average over time, these deviations were negligible. Conversely, under-fractions of the ionization fraction led to a reduction in intensity down to 0.9 during the condensation phase and the quasi-steady state. Given the dependency of the over/under-fractionation on mass circulations between the chromosphere and the corona, these observations will serve as valuable benchmarks to validate not only Alfvén wave models but also other existing mechanisms on coronal heating.


INTRODUCTION
The investigation of extreme ultraviolet (EUV) emission lines stands as a foundational approach for comprehending high-temperature, low-density plasmas present in the solar corona.This comprehension serves as vital information essential to understanding coronal heating mechanisms.A substantial portion of deductions drawn from EUV line analysis relies on ionization equilibrium calculations.These calculations have been effectively applied to establish models for the solar corona and transition region, estimate chemical abundances, and propose diagnostic methodologies for assessing coronal density and temperature (see a review by Del Zanna & Mason (2018) and references therein).
Among the generally successful EUV observations, a noteworthy exception lies in the observations of Li-and Na-like ions, which exhibit abnormal behavior compared to the other ions, as well documented in Del Zanna et al. (2002).In essence, the anomalous behavior refers to a consistent underestimation or overestimation, by substantial factors, of the theoretical intensities of lines emitted by these ions, following a Differential Emission Measure (DEM) analysis using lines from various isoelectronic sequences.This inconsistency has been noted by multiple researchers (Burton et al. 1971;Dupree 1972), and efforts to address it by incorporating density dependence in dielectric recombination and ionization from meta-stable levels initially appeared to resolve the anomalies (Vernazza & Raymond 1979;Raymond & Doyle 1981).However, improvements in atomic databases and photometric accuracy revealed that these anomalies persist (Judge et al. 1995).Despite the persistence of anomalies, notable improvements have been achieved in reducing the inconsistency by considering the effects of high densities, solar radiation, and charge transfer on ion formation (Dufresne et al. 2023).
The anomalous behavior exhibited by Li-and Na-like ions has prompted certain researchers to explore so-lutions beyond the assumptions commonly applied to numerous other spectral lines.One of the main assumptions is that all the ions are in a state of ionization equilibrium.However, when the plasma encounters a rapid increase in temperature, the assumption of instantaneous response in ion populations becomes invalid.This circumstance could lead to the potential existence of ions with relatively low charge at much higher temperatures than typically expected, consequently resulting in anomalous emission.Studies have investigated line formation during sudden temperature fluctuations, examining factors such as velocity redistribution, steady flow (Noci et al. 1989;Spadaro et al. 1990;Pietarila & Judge 2004), and abrupt heating or cooling (Hansteen 1993;Bradshaw & Mason 2003;Bradshaw et al. 2004).Efforts with multi-dimensional simulations have been undertaken to reconcile the disparities observed between theoretical predictions and actual observations (Olluri et al. 2013(Olluri et al. , 2015;;Martínez-Sykora et al. 2016).
However, performing multi-dimensional simulations remains a challenging task due to the thin transition region, where the majority of Li-and Na-like ions form.
The anomalous emission from Li-and Na-like ions, at least in part, appears to be linked to heating mechanisms and subsequent mass movements.Investigating this behavior holds potential importance in addressing the coronal heating problem.Reproducing this anomalous emission is important as a benchmark test for existing coronal heating models.This paper aims to explore the ability to reproduce the anomalous line formation process in a coronal model heated by Alfvén waves.To accomplish this, we employed a 1.5-dimensional magnetohydrodynamics (MHD) model (Moriyasu et al. 2004;Moriyasu & Shibata 2004) to predict the emissions from ions in non-equilibrium ionization (NEI).The hot corona is maintained via shocks produced through nonlinear mode conversions from Alfvén waves (Hollweg et al. 1982;Kudoh & Shibata 1999), and its behaviors in quasi-steady states (Antolin & Shibata 2010) and during long-term evolution (Washinoue & Suzuki 2019) have been extensively investigated.By solving a set of equations for NEI, we can quantitatively explore the behaviors of Li-and Na-like ions in the coronal loop with complex dynamics.Because the transition region, where most Li-and Na-like ions form, is a very thin layer, even small mass flows associated with evaporation, condensation, and shock propagation can occur in time scales shorter than the ionization and recombination time scale.This phenomenon may contribute to deviations in ionization fractions from ionization equilibrium, although detailed numerical simulations are re-quired to quantitatively estimate the degree of departure in ionization fraction and the associated emissions.Please note that our focus lies on the effects of bulk flows on NEI, although, in realistic situation, we can not decouple the impacts of alterations in charge states under high densities (Vernazza & Raymond 1979) and/or non-Maxwellian electrons (Dudík et al. 2014(Dudík et al. , 2017)).Including these effects will be the future studies.
This paper is organized as follows: Section 2 presents models and assumptions, Section 3 details our simulations and analysis with discussions, and Section 4 provides the conclusions of the study.

MODELS AND ASSUMPTIONS
In this study, we investigated the impact of NEI on emergent intensity by employing a coronal loop model heated by Alfvén waves.We concurrently solved the dynamic heating process governed by MHD equations and the evolution of the ion fraction.Subsequently, we reconstructed UV radiations based on the obtained ion fractions, electron density, and temperature.
We employed a 1.5-dimensional numerical model for the coronal loop, building upon the pioneering works of Moriyasu et al. (2004) and Moriyasu & Shibata (2004).This model naturally reproduces warm loops as a consequence of Alfvén wave injection from the photosphere.The hot corona is achieved through MHD shocks generated by nonlinear mode conversion from Alfvén waves (Hollweg et al. 1982;Kudoh & Shibata 1999).The model's properties have been extensively investigated, including parameter dependencies (Antolin & Shibata 2010;Matsumoto & Shibata 2010), and differences from the nanoflare model (Antolin et al. 2008).
The fundamental equations solved are as follows: The equation of mass conservation: the equation of momentum conservation along the field line: the equation of angular momentum conservation: the induction equation: and the equation of total energy conservation: In these equations, ρ and P g represent mass density and gas pressure, respectively; v s is the velocity along the field line; v φ is the toroidal velocity; B φ is the toroidal magnetic field strength normalized by √ 4π.E is the total energy density given by where γ is the ratio of specific heats and is taken to be 5/3.The variables g s and A represent gravitational acceleration along the field line and the cross-section of the flux tube, respectively.These variables depend on the shape of the flux tube, and we adopted the same shape as used in Moriyasu et al. (2004).The radius of the flux tube, r ≡ A/π, is described by the following set of equations: where z denotes the length along the axis of the flux tube and θ is defined as The parameters, k, z d , and H 0 , are set to 3.48, 6H 0 , and 200 km, respectively.The field line below the chromosphere is significantly inclined from the vertical direction, leading to p-mode leakage.Equation (5) includes source terms accounting for radiation, thermal conduction, gravity, and torque.For radiation, we employed a composite model that considers both optically thick and thin radiative loss mechanisms (Shoda & Takasao 2021).The radiative loss functions from optically thin plasma were estimated via CHIANTI (Del Zanna et al. 2021), assuming photospheric element abundances.This choice was made to ensure consistency, as these abundances are also utilized in the EUV synthesis, as elaborated upon in subsequent paragraphs.We did not account for the feedback of NEI in the radiative loss function, although it has been suggested that the impact of NEI on the radiative loss function is projected to be at most 5% in dense atmospheres after solar flares (MacNeice et al. 1984) or at most a factor of 2 to 4 in coronal loops (Mariska et al. 1982;Hansteen 1993).
We selected the Spitzer-type thermal conduction.We did not use any flux limiter (e.g.Klimchuk 2006) because the heat flux in our simulations consistently smaller than the saturated flux.To reduce numerical diffusion near the lower transition region, we implemented a broadening technique for temperatures below 0.15 MK, denoted as T c (Lionello et al. 2009).This technique involves adjustments to both the thermal conduction coefficient and radiative cooling.By implementing this technique, the width of the transition region below T < T c will be broadened by a factor of ∼ (T c /T ) 5/2 .Consequently, we achieved a reduction in the relative temperature difference between adjacent grid points, ∆ ln T , to less than 5% for the current grid size, thereby ensuring the allocation of more than 20 grid points for the transition region.We will discuss the limitations of these modifications in subsequent sections.
The amplitude of the torque at the footpoints is denoted as L trq .This torque was enforced only at the footpoints, and its amplitude was determined such that the root mean square of the toroidal velocity amplitude reached ∼ 1 km s −1 .This velocity amplitude is consistent with the horizontal velocities in the granular cells (van Ballegooijen et al. 1998;Matsumoto & Kitai 2010;Chitta et al. 2012;Oba et al. 2017).
We adopted an approximated equation of state in our model to include the effect of a change in molecular weight (Matsumoto & Suzuki 2014).Although this slightly modified the atmospheric structure below the chromosphere (T < 10 4 K) compared to constant molecular weight, the effects on the transition region and the corona should be subtle.
We investigated the impact of NEI for the elements C, N, O, Ne, Mg, Si, S, and Fe, as these elements harbor the most abundant ion species in the transition region and the corona.In this investigation, we present results specifically for C, N, O, Si, and S, which demonstrate notable deviations from ionization equilibrium.In this regard, we solved a series of NEI equations for these elements, expressed as follows: where N i denotes the number density of a specific element in the ith ionization stage while N e indicates the electron number density.The coefficients on the righthand side, S i and α i , indicate temperature-dependent collisional ionization and recombination rate coefficients obtained from CHIANTI atomic database 10.0.2.These rates are obtained in low density limit, although the density dependence removes some discrepancies between theory and observations on Li-and Na-like ions (Vernazza & Raymond 1979;Raymond & Doyle 1981;Dufresne et al. 2023).Note that we did not include the feedback effects on MHD variables through the radiative cooling function.
We developed an original MHD code capable of conducting precise and stable simulations, even in the low beta region surrounding the transition region.To calculate numerical flux, we employed the HLL-approximated Riemann solver (Einfeldt et al. 1991).Furthermore, conservative variables were reconstructed in each cell using a 3rd order weighted essentially non-oscillatory (WENO) scheme and subsequently integrated in time through the 3rd order Arbitrary Derivative Riemann Problem (ADER) scheme (Balsara et al. 2009).Considering that the time scale of thermal conduction is generally much shorter than that of dynamics, we adopted an operator split method, implicitly solving thermal conduction using the super time-stepping method (Meyer et al. 2012).Because the shortest time scale in equations ( 13) is often smaller than dynamical time scales, we implicitly solved the right-hand side of equations ( 13).
We applied the same boundary and initial conditions as outlined in Antolin & Shibata (2010).In brief, the initial conditions maintained hydrostatic equilibrium up to a height of 2 Mm, above which an artificially denser atmosphere was assumed to avoid severe CFL condition.The numerical domain spanned 103 Mm, including 2 Mm of subphotospheric layers at both ends to mitigate artificial oscillations arising from the boundaries.Grid spacing started at 10 km for the initial 16 Mm from both boundaries and gradually increased to 100 km in the central region.We initially computed the entire evolution of coronal loops using this coarse grid spacing.We have confirmed that this coarse run converged well with 10 km grid size (see appendix).Subsequently, we concentrated on three specific 20-minute time intervals and conducted additional calculations with finer grid sizes, reaching down to 1.25 km near the transition regions.This approach enabled us to reduce computational costs while effectively resolving the narrow transition region.We also have conducted a convergence test for this fine grid simulation (see appendix).

RESULTS AND DISCUSSIONS
In this study, we performed 1.5-dimensional MHD simulations of a coronal loop model subjected to Alfvén wave heating.Simultaneously, we solved the NEI equations for selected elements.Subsequently, by tracking the dynamic evolution of temperature, electron density, and ionization fractions, we calculated the emergent intensity employing the CHIANTI database.Our results unveil that specific Li-and Na-like ions display higher intensities during ionizing process when contrasted with the predictions based on the coronal approximation: the plasma is optically thin, the plasma have Maxwellian distribution function, the plasma is in ionization equilibrium, ionization and recombination rate coefficients do not change with density (zero density limit in this study).

Properties of coronal model
The system underwent both the evaporation (t < 200 min) and condensation (200 < t < 500 min) phases before attaining a quasi-steady state (t > 500 min) (Fig. 1).The coronal mass in Fig. 1(b) was normalized by the cross-sectional area at the base, denoted as A 0 , and defined as ρA/A 0 ds, where integration was performed over the region where the temperature exceeded 10 5 K.The particle flux during the evaporation and conden- sation phase was at most 5×10 14 and -1×10 14 cm −2 s −1 , respectively.In the quasi-steady state, the model reproduced a warm loop with an apex temperature of approximately 1.1 MK, featuring sharp transition regions between the chromosphere and the corona (the minimum temperature scale height of ∼ 70 km on average).The average electron density at the apex was 6.2 × 10 8 cm −3 , while the average length of the coronal loop was 79 Mm.These characteristics aligned with those of warm loops that have been extensively studied since the work of Aschwanden et al. (1999).At T=10 5 K, the velocity fluctuation ( v 2 LOS ∼ v 2 φ /2) measured approximately 25 km s −1 , consistent with the typical non-thermal line width observed in the transition region of the quiet sun (e.g.Peter 2001).There are no noticeable trends in the vertical flow (i.e.| v s | < 1 km s −1 ) within the transition region, although lowertemperature EUV lines typically exhibit a redshift (e.g.Peter & Judge 1999).
Two significant dynamics pertaining to the transition region are noteworthy: the collision of shocks and evap-oration/condensation processes.Firstly, the transition region exhibited vertical motion in response to the interaction with MHD shocks emanating from the photosphere (see Fig. 2).The typical altitude of the transition region measured approximately 9 Mm, with temporal variations of up to 3 Mm.These motions have previously been interpreted as spicule motion (Hollweg et al. 1982;Kudoh & Shibata 1999).This spicular dynamics did not change significantly through the whole calculation after the formation of the hot corona.Secondly, the shock-induced heating coincided with evaporation and condensation processes (see Fig. 1b), resulting in the circulation of materials between the corona and the chromosphere.These phenomena contributed significantly to the dynamic ionization and recombination processes occurring within the region, rendering them pivotal factors for consideration in the examination of emission originating from the transition region.Subsequent subsections of this study will delve into these aspects in greater detail.

Ionization fractions in the evaporation phase
Significant departures from ionization equilibrium were observed in the middle of corona as well as near the transition region.Figure 3 presented a snapshot of the temperature, density, and ionization fractions of C ions during the evaporation phase.The spatial distribution of ionization fractions exhibited two notable features.First, in the middle corona, the spatial profile of ionization fractions appeared considerably smoother than that in the equilibrium state (Fig. 3b).This phenomenon can be attributed to the longer ionization time scale, which prevents the ionization fraction from promptly responding to temperature fluctuations induced by shocks and acoustic waves.The corresponding temperature and density fluctuations can be found in the figure 3a.Previous simulations with time-dependent heating rates have demonstrated that NEI effects can arise due to temperature fluctuations, even in the absence of mass motion (Mariska et al. 1982).Observations have also indicated that the variability time scales are often constrained by ionization processes, regardless of the underlying atmospheric dynamics (Golub et al. 1989).Second, the peak of the ionization fraction for Li-and Na-like ions (e.g., C IV in Fig. 3c) was frequently displaced to higher altitudes.Consequently, the distributions of C III and C V were also shifted to higher altitudes.
The aforementioned deviations in the transition region were also evident in the probability distribution function (PDF) of the ionization fraction in temperature space (Fig. 4).Due to the effects of NEI, the ionization fraction did not solely depend on the local temperature, leading to a certain distribution with finite width at a specific temperature.We discretized the temperature space into bins of ∆ log T = 0.05 and computed PDFs, denoted as F i , for each temperature bin as follows: where i serves as the index for ion species; P (n i ≤ x) represents the probability that a particular ionization fraction, n i , is smaller than x at temperature T .To estimate PDF in the evaporation phase, we used the data in ionization phase (135 ≤ t ≤ 145 min).The solid lines in Fig. 4 depicted the ranges of 1-sigma levels for PDF (F i ∈ [0.17, 0.83]), whereas the dashed lines illustrate the ionization fraction at the equilibrium state.The overfraction and shift toward higher temperature was found for C IV while other C ions almost followed the ionization fraction in equilibrium.These characteristics differ from those documented in previous studies (Olluri et al. (2013(Olluri et al. ( , 2015)); Martínez-Sykora et al. ( 2016)), where substantial deviations in the ionization fraction were observed, even for ions among non Li-and Na-like ions.This discrepancy can be mainly attributed to the lower density in the transition region in their simulations, while the impact of numerical diffusion cannot be disregarded.The outcomes obtained from the coarse grid (∼10 km width), demonstrating finer resolution than the earlier multi-dimensional simulations, revealed marked deviations from ionization equilibrium due to numerical diffusion, notably noticeable among non Li-and Na-like ions.
The reason for the enhancement of C IV fraction can be primarily attributed to the evaporation from the chromosphere to the corona.In figure 5, we depicted evolution of a certain plasma parcel co-moving with fluid to trace its temperature, electron density, position, and ionization fractions.From t = 137.5 min, the ionization fraction of C IV increased as the plasma temperature increased to its formation temperature.As time went on, the plasma parcel evaporated to coronal temperature from t = 142.5 min in 40 sec, and then, the ionization fraction started to decrease.Because the evaporation time scale of ∼ 40 sec was comparable to the recombination time scale, the ionization fraction stayed larger than that in equilibrium during the evaporation process.The similar over-fractionation in Li-and Na-like ions can be found in the steady flow solution (Noci et al. 1989;Spadaro et al. 1990;Pietarila & Judge 2004).The over-fraction of C IV is also shown in the evaporation phase in nanoflare-heated loops (Hansteen 1993) and impulsive phase in flaring loop (Bradshaw et al. 2004).Throughout the evaporation process, this plasma parcel encountered shocks at t ∼ 140 min and 140.8 min.Upon interaction with the shocks, the C IV ionization fraction increased spontaneously.Despite these collisions occurring within a 10-second timeframe, which is shorter than the ionization time scales, the C IV ionization fraction closely tracked the equilibrium fraction.Considering that C III and C V displayed over-and under-fraction, respectively, the C IV ionization fraction was likely maintained by an enhancement of ionization and a reduction in recombination during the shock passage.Nano-flare-generated waves entering to the transition region from the corona could also modify the ionization fraction (Hansteen 1993), although we did not observe this phenomena probably because our model lacks nanoflares.

Ionization fractions in the condensation and quasi-steady phase
During the condensation phase, the C IV ions exhibited, on average, an under-fraction.Figure 6 displayed the pertinent physical characteristics of the condensing plasma parcel.Due to the finite duration required for the recombination process from C V to C IV, which exceeded the condensation timescale, the ionization fraction of C IV remained below that of the equilibrium state.This under-fraction is qualitatively consistent with the simulation by Hansteen (1993) during the condensation phase in a nanoflare-heated loop and by Noci et al. (1989); Spadaro et al. (1990) in the down flowing siphon flow.Contrary, the over-fraction is reported in the inflowing model of Pietarila & Judge (2004), which might be due to the difference in atmospheric model used, boundary conditions at the corona, or amplitude of mass flux.It was observed that the plasma parcel experienced a collision with a shock at t ∼ 311 min, during which the C IV ionization fraction once again closely followed the equilibrium fraction, as was noted during the evaporation phase.Despite the plasma undergoing repeated evaporation (e.g., at t ∼ 307.8 min) and condensation episodes during the whole condensation phase, the average ionization fraction of C IV was consistently smaller than that in the equilibrium state.
The shocks in themselves do not alter the average ionization fraction of C IV when averaging over time.In the quasi-steady state, where there is no tendency towards continuous evaporation or condensation, we can delineate the impacts of temperature and density fluctuations caused by shocks on NEI.Since there are no substantial deviations from ionization equilibrium in the quasi-steady state, we assert that the processes of evaporation/condensation are necessary to induce over/under fractions within the present model.

Synthetic UV intensity
The synthetic UV intensities obtained using the coronal approximation were generally smaller than those derived with NEI by a maximum of 40% during the evaporation phase (Table 1).The UV intensity computation under the coronal approximation is governed by the following formula: Here, Ab(X) represents the abundance of element X, assuming photospheric abundance.The function C(T, λ, N e ) corresponds to the contribution functions for each spectral line, calculated using the CHIANTI software.For NEI plasma intensities, we modified the integrand in eq. ( 15) by multiplying it with the ionization fraction ratio, N i;NEI /N i;EI .Among several spectral lines, we specifically highlighted lines with a ratio I EI /I NEI less than 0.8 during the evaporation phase and an intensity greater than 10 erg cm −2 s −1 sr −1 .While the effect of NEI significantly enhanced the intensity through the over-fraction during the evaporation phase, our model still exhibited discrepancies between R (defined as I EI /I NEI ) and I th /I obs , as observed in quiet sun (Judge et al. 1995;Dufresne et al. 2023).
The observational data utilized to determine the intensity ratio is obtained from a complex combination of several instruments documented in Judge et al. (1995) and Dufresne et al. (2023).In short, these datasets were obtained from the quiet region.However, it is noteworthy that the data employed in Judge et al. (1995) are spatially unresolved, potentially allowing contributions from active regions to be present.The underfraction observed during the condensation and quasisteady phases resulted in the ratio R < 1.1 at most.The observed intensity ratios are still significantly smaller than those from our model, although we can not directly compare the observations in quiet regions and the model in active region loop.It is anticipated that the influence of non-equilibrium ionization (NEI) will be more pronounced in the regions characterized by lower density, particularly in the quiet region.

The effect of broadening technique of the transition region
In our study, we implemented a numerical broadening technique for the transition region to mitigate numerical diffusion effects near the lower transition region, which was identified as a factor reducing the impacts NEI.To investigate the impact of this technique, we conducted simulations with a broader transition region, specifically setting T c = 0.2 MK.The results of this simulation revealed a increase in the ratio of R, such as an increase from 0.6 to 0.8 for C IV ions 384.17 Å, and similar increases were observed for other Li-and Na-like ions.This increase in the intensity ratio can be attributed to the broadening of the lower transition region, which results in a longer dynamical time scale for traversing the layer compared to realistic conditions, thus potentially reducing the influence of NEI.Consequently, we consider the intensity ratio obtained under our parameter settings to represent an upper limit in this context.

Mass circulation
The magnitude of over-fraction evolved through the mass circulation process occurring between the corona and the chromosphere, as our model has unveiled distinct ionization fractions across various phases.Our model disclosed that evaporation occurred prior to t < 200 min, followed by the onset of condensation until approximately t ∼ 500 min.The duration of these mass cycles could be affected by the choice of abundance thorough the radiative cooling function.Because we used the photospheric abundance in this study, the cooling time scale could be shorter when using the coronal abundance instead.NEI effect on the radiative cooling function may slightly increase the cooling time scale of the coronal loop (Bradshaw & Mason 2003) if we include the feedback effects.Subsequently, the system maintained a quasi-steady state, at least until t = 800 min, though it is worth noting that the continuity of this quasi-steady state beyond that point is not guaranteed.Similar models with significantly extended computational duration have elucidated the cyclic evolution of coronal loops (Washinoue & Suzuki 2019), a phenomenon that could be attributed to the thermal instability inherent in coro-  nal loops (Kuin & Martens 1982).While providing a quantitative estimate for the amplitude of over-fraction remains challenging, it is reasonable to anticipate that the ratio R presented during the evaporation phase (as detailed in Table 1) would undergo increase when taking the condensation phase into consideration.

CONCLUSIONS
In this study we have conducted 1.5-dimensional MHD simulations for Alfvén-wave-heated coronal loops simultaneously solving a series of NEI equations for several ion species.After introducing the Alfvénic fluctuations at the foot point of the loop, the system experiences the evaporation and condensation phases before it reaches the quasi-steady state.During the evaporation phase, the over-fractionation of Li-and Na-like ions results in the intensity ratios, R, as low as 0.6.Conversely, during the condensation and quasi-steady phases, the underfractionation leads to R values as high as 1.1.These pronounced fluctuations in ionization fractions are primarily attributed to the processes of evaporation and condensation between the corona and the chromosphere.Interestingly, the collision with shocks do not significantly deviate the C IV ion fraction from the ionization equilibrium.
While our model has successfully demonstrated a reduction in the intensity ratios R during the evaporation phase, a noticeable gap between the model and observational data still persists.This may suggest that the heating mechanism in quiet regions exhibits more sporadic behavior, leading to a more dynamic mass circulation between the chromosphere and the corona.In such instances, deviations from equilibrium ionization could be pronounced, resulting in a reduced intensity ratio during the evaporation phase.Concurrently, this leads to an increased intensity ratio during the condensation phase.Therefore, it is not straightforward to assert that stronger mass circulation unequivocally leads to a smaller intensity ratio when averaged over the entire period.Bridging this gap also necessitates further investigations into atomic physics aspects, such as density effects, photoionization, and charge transfer, which collectively have the potential to align the theoretically predicted UV intensities of transition region spectral lines with observations (Dufresne et al. 2023).
The observed over-fractionation of Li-and Na-like ions holds significant scientific interest, as it may serve as an indicator of mass motions closely linked to coronal heating mechanisms and mass loss processes.This phenomenon gains relevance in the context of the wealth  of supporting evidence for impulsive heating events that drive the cycle of evaporation and condensation in the solar atmosphere (Klimchuk 2006).Consequently, the extent of over-fractionation could potentially offer valuable constraints on the amplitude and frequency of such impulsive heating events (Bradshaw & Klimchuk 2011).Furthermore, it's worth noting that the solar wind are known to expand the UV-observed solar atmosphere via the effects of NEI (Neupert 1964;Mariska et al. 1978;Dupree et al. 1979).Consequently, the degree of overfraction could be considered as a metric for mass loss processes in the Sun.Importantly, anomalies in the ionization states of Li-and Na-like ions have also been identified in stellar atmospheres (Del Zanna et al. 2002), suggesting that the observation of these ions may provide valuable insights into the heating and mass-loss mechanisms in other stars.
Multi-dimensional simulations that incorporate NEI effects are crucial for conducting a quantitative comparison between theoretical models and complex observational data (Olluri et al. 2013(Olluri et al. , 2015;;Martínez-Sykora et al. 2016).However, such simulations remain a formidable challenge, even with the current computational resources at our disposal.This challenge arises from the necessity to accurately resolve the narrow transition region when solving the advection equations for Li-and Na-like ions.The width of the transition region where those ions form, determined by the Field length or κT /|ρL| (Begelman & McKee 1990), was sometimes going down to ∼ 20 km at T ∼ 10 5 K in our simulation.Notably, this width is broadened by a factor of (T c /10 5 K) 5/2 ∼ 2.8 with our broadening technique, effectively reducing the width to approximately 7 km in the realistic situation.To achieve a high-resolution representation of this thin layer, typically requiring at least 10 grid points to mitigate numerical diffusion, we find it necessary to employ the adaptive mesh refinement schemes, even in one-dimensional simulations (Bradshaw & Mason 2003).The TRAC method proposed in Johnston & Bradshaw (2019) is a robust technique for accurately determining coronal density regardless of grid size.However, it may not improve the situation when solving the advection equations in ionization fractions because the transition region is resolved with only a few grids in this method.Insufficient resolution to capture the transition region leads to pronounced numerical diffusion, giving rise to unphysical deviations from ionization equilibrium.
The extent of over-fractionation is likely contingent on the mass flux or the rate of mass exchange between the corona and the chromosphere.However, our current model does not comprehensively explore this particular parameter.To conduct a thorough investigation of these parameters, it would be advantageous to employ timesteady solutions (Noci et al. 1989;Spadaro et al. 1990;Pietarila & Judge 2004;Gilly & Cranmer 2020).By deriving the over-and under-fractionation behaviors of Liand Na-like ions as functions of mass flux, we could potentially establish valuable constraints on coronal heating mechanisms and mass loss rates.These constraints would be particularly informative when compared with UV observations in future observations from Solar Or-biter/SPICE (Fludra et al. 2013) or Solar-C/EUVST (Shimizu et al. 2020), which cover large wavelength range.
We extend our heartfelt gratitude to the reviewer for their insightful comments and constructive feedback, which greatly enhanced the quality of this manuscript.This work was supported by JSPS KAKENHI Grant Number JP23K03456.This study was carried by using the computational resource of the Center for Integrated Data Science, Institute for Space-Earth Environmental Research, Nagoya University through the joint research program.CHIANTI is a collaborative project involving George Mason University, the University of Michigan (USA), University of Cambridge (UK) and NASA Goddard Space Flight Center (USA).

Software: CHIANTI (Del Zanna et al. 2021) APPENDIX
We performed two distinct convergence tests to explore the resolution dependency of our findings.The first test pertains to the coarse simulation, which served as the initial conditions for the fine-grid simulations.In this coarse simulation, we employed a grid size of 5 km near the transition region, whereas the original grid size was 10 km. Figure 7 illustrates the entire evolution of the maximum temperature (solid black line) and the density at the apex (dashed black line), alongside the results obtained with a 10 km grid size (depicted by red lines).Based on this comparison, we infer that the maximum temperature and the coronal density exhibit minimal sensitivity to resolutions finer than 10 km.
The second test pertains to the fine-scale simulations employed for computing the ionization fractions.We varied the resolution from 10 to 1.25 km and depicted the ionization fraction at t=145 min (∼ 20 min after the starting time of fine-scale simulation) in Figure 8.The colors represent the ionization fraction of Carbon ions, following the same manner as depicted in Figure 3. Results obtained with grid sizes of 10, 5, 2.5, and 1.25 km are denoted by solid, dashed, dotted, and dash-dotted lines, respectively.A comparison between the results obtained with grid sizes of 2.5 km and 1.25 km indicates that the ionization fraction has nearly converged with a grid size of 1.25 km.While the ionization fraction in the coronal region (C V, C VI, C VII) exhibits convergence with a grid size of 10 km, a resolution of 1.25 km is necessary to achieve convergence for transition region ions such as C III and C IV.

Figure 1 .
Figure 1.Temporal evolution of (a) the maximum temperature and (b) the coronal mass normalized by the cross sectional area at the photosphere.

Figure 2 .
Figure 2. Time distance diagram of (a) temperature and (b) velocity along the field line.The black solid line in panel (b) indicates the contour at T = 10 5 K that represent the transition region.

Figure 3 .
Figure 3. (a) coronal temperature (solid line) and density (dashed line) as a function of s at t=135.8 min.(b) A snapshot of ionization fractions of C ions.(c) Zoomed-in view of (b) focused on the transition region.For panels (b) and (c), solid lines depict the NEI fractions, while dashed lines represent the ionization fraction in equilibrium state.

Figure 4 .
Figure 4. Ionization fractions of C ions as a function of temperature.Solid lines depict 1 sigma intervals of PDF of the NEI fractions.Dashed lines represent the fractions in ionization equilibrium.

Figure 5 .
Figure 5. (a) Temperature (dashed line), electron density (dotted line), position (solid line), and (b) ionization fractions of a traced particle in the evaporation phase.The solid and the dotted lines in panel (b) indicate the ionization fraction in non-equilibrium and equilibrium, respectively.

Figure 6 .
Figure 6.The same format as in the figure 5 in condensation phase.
−2 s −1 sr −1 ]; Re, Rc -ratio of IEI to INEI during the evaporation and condensation phase, respectively; I th /I obs -ratio of theoretical predictions of intensity to observed intensity ( a Judge et al. (1995) and b Dufresne et al. (2023)).

Figure 7 .Figure 8 .
Figure7.Resolution dependence on the evolution of the maximum temperature (solid) and coronal density at apex (dashed).The red and black lines represent the results from coarse simulation runs with 10 and 5 km, respectively.

Table 1 .
Comparison of intensities from ionization fractions in equilibrium and non-equilibrium during the evaporation phase.Note.Ion -emitting ion; Seq -ion isoelectronic sequence; λ -wavelength; IEI and INEI -predicted intensity in EI and NEI during the evaporation phase [ergs cm