Oscillation and Evolution of Coronal Loops in a Dynamical Solar Corona

Observations have revealed two regimes of kink oscillations of coronal loops. Large amplitude oscillations excited by impulsive energy releases such as coronal mass ejections are characterized by their strong damping by resonant absorption. Lower amplitude oscillations may also be excited and sustained by the ubiquitous motions present in the corona and so are characterized as being decayless. We perform numerical simulations to study the oscillation and evolution of coronal loops in a dynamical environment. We investigate the observational signatures of kink oscillations and the Kelvin-Helmholtz instability in terms of high-resolution seismological and spatial data analysis techniques. We find that low amplitude kink oscillations are capable of generating significant changes in the loop profile which can affect estimates of the transverse loop inhomogeneity based on seismological and forward modeling methods. The disparity between methods may be indicative of non-linear evolution of coronal loops. The influence on forward modeling estimates could also account for previous observational evidence favoring loops having wider inhomogeneous layers.


INTRODUCTION
Kink oscillations of coronal loops are periodic displacements of the loop axis. They were first observed using the Transition Region And Coronal Explorer (TRACE; Handy et al., 1999) in an active region perturbed by a solar flare (Aschwanden et al., 1999;Nakariakov et al., 1999). Hundreds of observations of standing kink modes have now been studied (e.g., catalogs by Goddard et al., 2016;Nechaeva et al., 2019), aided by the full disk coverage and enhanced temporal resolution of the Atmospheric Imaging Assembly (AIA) on board the Solar Dynamics Observatory (SDO; Lemen et al., 2012).
Large amplitude kink oscillations are observed to be strongly damped which is attributed to the process of resonant absorption (e.g., De Moortel et al., 2016, and references therein). This is a robust damping mechanism since it only requires that the transition between a higher density coronal loop and the lower density background plasma occurs over a finite spatial scale. In analytical and numerical models this transition is usually described in terms of a boundary layer surrounding the core of the loop. The inhomogeneous layer width (commonly normalized by the loop radius to produce the parameter ǫ) represents a transverse spatial scale which is also important for other physical processes such as phase mixing (Heyvaerts and Priest, 1983) and the Kelvin-Helmholtz instability. It is therefore a key physical parameter but difficult to observe directly, which motivates its seismological inference through the damping of kink oscillations.
The length L of a coronal loop is typically much larger than its radius R. Under this thin tube approximation, the phase speed of the kink mode is the kink speed C k , which is a density-weighted average of the internal and external Alfvén speeds. Under the further thin boundary approximation corresponding to the limit ǫ → 0, the period of oscillation for the fundamental standing kink mode is then P = 2L/C k . However, parametric studies (Van Doorsselaere et al., 2004;Soler et al., 2014;Pascoe et al., 2019a) find that the period of oscillation does also depend on the width of the boundary layer when it is sufficiently large.
The boundary layer of a coronal loop provides a continuous variation of the local Alfvén speed, and the resonance condition is satisfied where the kink mode frequency matches the local Alfvén frequency. A transfer of wave energy from collective kink motions to localized azimuthal perturbations (m = 1 Alfvén waves) will appear as a damping of the periodic loop displacement. The azimuthal motions will typically appear as unresolved Doppler velocity perturbations due to line-of-sight (LOS) integration of multiple waves and structures (e.g., De Moortel and Pant et al., 2019). Although this transfer is a linear and ideal process, subsequent phase-mixing of the Alfvén waves in the inhomogeneous layer can generate small spatial scales which enhance dissipative processes such as viscosity and resistivity (e.g., Pagano and De Moortel, 2017;Pagano et al., 2018).
Initial descriptions of the resonant absorption of kink oscillations took the form of exponential damping profiles (Goossens et al., 2002(Goossens et al., , 2012Ruderman and Roberts, 2002;Arregui et al., 2007) based on analytical derivations for the asymptotic state of the system under the thin tube thin boundary (TTTB) approximation. Numerical simulations  found poor agreement with exponential damping behavior at early times, where a Gaussian damping profile appeared more suitable. The analytical description of Hood et al. (2013) accounted for this Gaussian damping regime as an approximation for the initial behavior of the kink mode, with exponential damping behavior at a later stage. The switch from Gaussian to exponential damping behavior was shown to depend on the loop density contrast ratio ζ , which was confirmed by the parametric study by Pascoe et al. (2013). It was also used to propose a general damping profile (GDP), combining both the Gaussian and exponential damping regimes, as a convenient approximation for seismology. The Gaussian damping behavior was initially studied for propagating kink waves but is also applicable to standing waves (e.g., Ruderman and Terradas, 2013;Magyar and Van Doorsselaere, 2016a;Pagano et al., 2018). In particular, Pascoe et al. (2019a) performed 300 numerical simulations to produce a lookup table (LUT) as a convenient means of estimating the damping profile beyond the applicability of the thin boundary approximation. Resonant absorption therefore accounts not only for the strong damping of kink oscillations but also the observational evidence for their non-exponential damping behavior Ireland and De Moortel, 2002;Goddard et al., 2016;Morton and Mooroogen, 2016;Pascoe et al., 2016c), and allows the seismological inference of the transverse density profile of the loop (Pascoe et al., 2016b(Pascoe et al., , 2017a. The Kelvin-Helmholtz instability (KHI; e.g., Browning and Priest, 1984;Ofman et al., 1994) has been demonstrated by numerous numerical studies to cause redistribution of loop plasma during kink oscillations (e.g., Terradas et al., 2008a;Soler et al., 2010;Antolin et al., 2014;Karampelas and Van Doorsselaere, 2018). Since the rate of the instability depends on the shear flows it is also affected by the width of the inhomogeneous layer (e.g., Magyar and Van Doorsselaere, 2016a;Goddard et al., 2018). The "decayless" regime of kink oscillations refers to low amplitude oscillations which appear to be undamped (Anfinogentov et al., 2013(Anfinogentov et al., , 2015Nisticò et al., 2013) or even growing (Wang et al., 2012), in contrast with the large amplitude and strongly damped oscillations discussed above. Only coronal loops with a discontinuous boundary would not exhibit damping due to resonant absorption. However, the strong shear velocities that would be present in such a scenario would also lead to the rapid development of KHI and the establishment of an effective boundary layer by mixing of internal plasma. Forward modeling of EUV emission from coronal loops (e.g., Van Doorsselaere et al., 2016) suggests observational effects associated with KHI might contribute to the appearance of weak damping in some circumstances (Antolin et al., , 2017, but generally the presence of damping by resonant absorption must be compensated by the driving mechanism, such as persistent footpoint oscillations (e.g., Afanasyev et al., 2019Afanasyev et al., , 2020Karampelas et al., 2019) or flows Karampelas and Doorsselaere, 2020). Decayless kink oscillations might be connected with the observation of propagating transverse velocity perturbations (e.g., Tomczyk et al., 2007;Tomczyk and McIntosh, 2009;Verth et al., 2010;Pascoe et al., 2015) since both appear to be ubiquitous and require a persistent source of wave energy, though spatial analysis of the decayless kink modes confirms they are indeed standing modes (Duckenfield et al., 2018) as does the period of oscillation being proportional to the loop length (Anfinogentov et al., 2015).
In this paper we consider the effect of the Kelvin-Helmholtz instability on the distortion of the loop density profile and its appearance in optically thin EUV images. Our simulation of a coronal loop experiencing multiple perturbations which excite and maintain standing kink oscillations approximating the decayless regime is presented in section 2. In section 3, we consider the observational signature of our simulation and compare the results with a statistical study of the density profiles of coronal loops. Conclusions are presented in section 4.

NUMERICAL SIMULATION OF NON-LINEAR EVOLUTION
In this section, we describe our numerical simulation of a coronal loop which experiences multiple low-amplitude impulsive perturbations designed to excite the fundamental standing kink mode. Our use of multiple perturbations as opposed to a continuous driver facilitates analysis of the evolution of the loop and its oscillation in time by enabling us to use established techniques for the damping behavior of impulsively-generated decaying oscillations. The coronal loop is modeled as a straight flux tube embedded in a uniform magnetic field with no gravitational stratification. The only non-uniformity is thus provided by the coronal loop. We use a linear density profile in the inhomogeneous layer for comparison with previous work, in particular the parametric study by Pascoe et al. (2019a) which describes the effect of a thick inhomogeneous layer. We note that the numerical simulations in that paper studied the linear regime of kink oscillations and used cylindrical coordinates, with the assumption of m = 1 symmetry for the azimuthal dependence corresponding to the kink oscillation and Alfvén waves generated by resonant absorption (e.g., Pascoe et al., 2010). In this paper, we study non-linear evolution of the loop, which includes the breaking of m = 1 azimuthal symmetry by the development of KHI (e.g., the large deformations in Karampelas and Van Doorsselaere, 2018). We use LARE3D (Arber et al., 2001) which solves the non-linear MHD equations in Cartesian geometry.
The loop has a length of 100 Mm, radius of 1 Mm, and the uniform magnetic field strength corresponds to an external Alfvén speed of 1 Mm/s. The loop has a density contrast ratio of ζ = 3 for which the TTTB estimate for the period of oscillation of the fundamental kink mode is P k ≈ 283 s, allowing several cycles of the oscillation to be observed between each perturbation (1, 000 s apart). However, each subsequent perturbation will be applied at a time which does not correspond to a zero in the oscillation of the current oscillation. The numerical domain extends to [−4, 4] Mm in the x and y directions perpendicular to the loop axis with a resolution of 400 × 400 × 200 grid points. Numerical tests indicate this is sufficient for the problem, and the need for significantly higher resolution in the transverse directions. Boundary conditions are periodic in the x-and ydirections, and line-tied in z. We consider a low-beta plasma (β = 0.001) and the loop is initially set in pressure balance, corresponding to a cooler loop. We do not consider nonadiabatic effects (e.g., Aschwanden and Terradas, 2008;Magyar et al., 2015;Ruderman et al., 2019).
An impulsive fundamental standing kink mode perturbation is applied to the x-component of the plasma velocity every 1, 000 s. Our aim is to demonstrate the evolution of the loop and oscillation parameters over several iterations of the perturbation. For this purpose we use a inhomogeneous layer width ǫ = 0.3 which is sufficiently large to provide significant (readily measurable) damping for the oscillation while not too large as to inhibit the development of KHI. Similarly, the amplitude of the velocity perturbations is approximately 1% of the Alfvén speed to allow KHI to develop at a moderate rate. Our density profile parameters ζ = 3 and ǫ = 0.3 are reasonable for a coronal loop, for example the loop analyzed on 2012-03-09 by Pascoe et al. Our driver is intended to efficiently generate the fundamental standing kink mode in the loop, but does not correspond to the exact eigenfunction which would generally be unknown as the loop evolves during the simulation. We use the spatial form v x = Asech ν (r/R) cos kz with A = 0.015 Mm/s, k = π/L corresponding to the fundamental longitudinal harmonic, and ν = 3 providing an appropriate spatial scale in the radial (r = x 2 + y 2 ) direction. We note that the form of our driver and loop also excludes higher harmonics which can be readily generated by a broadband footpoint driver (e.g., Afanasyev et al., 2020). A persistent footpoint driver is likely more realistic but makes it difficult to accurately measure the period of oscillation due to the presence of the driver spectrum in addition to the natural frequencies. The presence of longitudinal stratification due to gravity, expansion, or curvature leads to coupling of longitudinal harmonics (e.g., Pascoe et al., 2009;Pascoe and Nakariakov, 2016) and changes to the ratio of their periods (e.g., Andries et al., 2005;McEwan et al., 2006;Safari et al., 2007;Verth and Erdélyi, 2008;Arregui et al., 2013) but does not significantly affect the damping rate due to resonant absorption Dymova and Ruderman, 2006). The period ratio for harmonics has been observed in many studies of large-amplitude oscillations (e.g., Verwichte et al., 2004;De Moortel and Brady, 2007;Van Doorsselaere et al., 2007;Pascoe et al., 2016aPascoe et al., , 2017aDuckenfield et al., 2019) and recently in decayless oscillations (Duckenfield et al., 2018). Figure 1 shows the evolution of the transverse loop density profile at the middle of the coronal loop, i.e., the antinode for the fundamental standing kink mode. The panels show the density profile in normalized units (external density ρ e = 1). The cross symbols in Figure 1 denote the center of mass of the loop. The dashed line represents the contour for which the local Alfvén speed matches the kink speed. Figure 2 demonstrates the evolution of the coronal loop and its oscillation during our simulation. The evolution of the loop profile is characterized by fitting the transverse density profile with a model based on a circular cross-section and linear density profile in the inhomogeneous layer, which is satisfied exactly for the initial condition but otherwise an approximation. They therefore correspond to a value which is averaged in the azimuthal direction, whereas it is known from previous work (e.g., Barbulescu et al., 2019;Hillier and Arregui, 2019) and our Figure 1 that KHI causes density perturbations to develop which exhibit the same m = 1 symmetry as the kink mode. For example, the loop becomes extended in the direction of its displacement (e.g., Karampelas and Van Doorsselaere, 2018). We note that this method for estimating loop parameters is only possible with simulation data, whereas the observational signatures will be considered in section 3.
The top left panel of Figure 2 shows the position of the coronal loop based on our fitting of the density profile (solid line). The dashed line is the position of the center of mass of the loop which is in good agreement. We apply a perturbation at t = 0 s and every 1, 000 s thereafter, indicated by the vertical dotted lines. The subsequent increase in amplitude of the oscillation depends on whether the perturbation is aligned with or opposed to the current direction of motion (which we choose arbitrarily to demonstrate both cases). We note that the abrupt changes in phase is a consequence of our driving mechanism, which is chosen for convenience of our analysis, and is not observed in decayless oscillations. The clear signature  Frontiers in Astronomy and Space Sciences | www.frontiersin.org of damping and then an increase in amplitude also does not represent decayless oscillations with an approximately constant amplitude, although we have deliberately chosen a sufficiently long interval between perturbations to allow seismological analysis and this effect would be reduced by more frequent and lower amplitude perturbations.
The maximum amplitude of 0.3 Mm is consistent with the decayless regime of kink oscillations which is lower than 1 Mm and an average of approximately 0.2 Mm (Anfinogentov et al., 2015). It is evident that the damping rate at the end of the simulation is greater than that at the start of the simulation (further discussion below). The other panels in Figure 2 show the evolution of the loop parameters. There is no significant change in the loop radius R or density contrast ratio ζ . The main effect of KHI is the increase in the width of the inhomogeneous layer, from its initial value of 0.3 to approximately 1 at the end of our simulation. The onset of this increase occurs after several cycles of the oscillation due to the time required for KHI to develop. The widening inhomogeneous layer and unchanged loop radius are consistent with the study of non-linear evolution of loops by Goddard et al. (2018), although that paper considered the changes in the transverse intensity profile rather than the density profile here. Goddard et al. (2018) also reported a decrease in the peak intensity of the loop, whereas our analysis here indicates an unchanged density contrast ratio. There is no inconsistency in these two results since the decrease in intensity can also be associated with the mixing of plasmas of different temperatures by KHI, and numerical simulations such as those in Goddard et al. (2018) typically impose temperature profiles which decrease inside the loop to maintain pressure balance. Figure 3 shows the oscillation of the coronal loop after the first applied perturbation and before the second. During this time there are minimal effects due to the non-linear evolution of the loop (see Figure 2) and so the damping rate is consistent with the black envelopes representing the linear solution calculated with the lookup table (LUT) of Pascoe et al. (2019a). The LUT envelope contains the correction for the boundary layer not being thin, whereas each of the colored envelopes are based on the thin boundary approximation and so overestimate the damping of the kink oscillation, particularly the exponential damping profile (red) which does not account for the Gaussian damping behavior during the first couple of cycles of the oscillation. For the general damping profile (GDP) the switch from a Gaussian profile (blue) to an exponential profile occurs after 2 cycles for a density contrast ratio of ζ = 3. The damping rate during subsequent cycles is increased due to KHI, as studied by Magyar and Van Doorsselaere (2016a). We also note that the fitted period of oscillation is P ≈ 1.01P k , consistent with the small increase due to the presence of a non-thin boundary layer also found in numerical simulations by Pascoe et al. (2019a), and contrary to the decrease in period reported by Soler et al. (2013) and Soler et al. (2014). Figure 4 shows the results of wavelet analysis (Torrence and Compo, 1998) of the oscillation of the loop axis (left panel). The contours for the wavelet power suggests a small decrease in the period of oscillation in time. This effect is also found by fitting the oscillation as a series of damped perturbations with a time-varying period of oscillation and inhomogeneous layer width (right panel) based on the LUT damping profile. The fitted period of oscillation starts as slightly larger than the TTTB approximation (dashed line in left panel) due to the thickness of the inhomogeneous layer as discussed above. The period of oscillation then decreases during the simulation, due to the non-linear evolution, to a value slightly less than the TTTB approximation. Overall, the TTTB approximation remains a good approximation for the period of oscillation.
The seismological estimate for ǫ is based on the damping due to resonant absorption calculated for different density profile parameters (ζ , ǫ) by Pascoe et al. (2019a) and summarized in the LUT 1 . This accounts for the effect of a thick inhomogeneous layer and the presence of the non-exponential damping regime. The loop model used in this paper is the same as that used to generate the LUT, and so the results may be directly compared with the relevant physical difference being the nonlinear evolution due to finite amplitude kink oscillations. We fit our oscillation with a sinusoidal oscillation comprised of five segments (corresponding to each of the applied perturbations) but a period of oscillation and inhomogeneous layer width which vary (as quadratic polynomials) continuously across all segments. The solid lines in Figure 4 represent the fitted dependence of P (t) and ǫ (t) required to account for the period and damping in our simulation. Within each individual segment the damping profile is calculated from the LUT using ζ = 3 and the local average for P (t) and ǫ (t).
The increase in ǫ required to account for the increase in damping rate in time is much larger (approximately twice) than the estimate based on fitting the transverse density profile. On the other hand, it is evident from Figure 1 that the smaller estimate is more reasonable since the loop has not become fully inhomogeneous by the end of the simulation as the seismological analysis with ǫ → 2 suggests. We can therefore conclude that the non-linear evolution increases the damping rate of kink oscillations above the level suggested by the increasing inhomogeneous layer width. This is not unexpected since our estimate of the inhomogeneous layer width only describes the large scale broadening effect of KHI due to the mixing of loop plasma and external plasma, while KHI is also known to generate small scale structures ("TWIKH rolls"). Figure 1 shows that these TWIKH rolls also modify the location where the condition for resonant absorption to occur is satisfied (dashed lines). In particular, near the locations of maximum mode coupling (x = 0, y = ±R), the TWIKH rolls cause the resonant contour to double back on itself, presenting additional locations at which resonant absorption can occur. Accordingly, we can caution that the seismologically-inferred value of ǫ is likely to be an overestimate in the case of loop oscillations with strong development of KHI. Previous studies have also estimated ǫ independently of seismological analysis using the transverse intensity profile of the loop. In the following section we consider how this estimate is also affected by the development of KHI.
The seismological model discussed above used a constant density contrast ratio and varying inhomogeneous width based on our understanding of the evolution of the loop due to KHI such as the results in Figure 2. However, approaching the oscillation from a purely seismological point of view, we could consider an alternative analysis of the increase in damping rate as due to a varying density contrast ratio with ǫ remaining constant. We find that this would require a significant increase to ζ > 15 and so could readily be refuted by a lack of brightening in EUV intensity images and due to a mechanism for large increases in loop mass (or decrease in surrounding plasma) to be theoretically accounted for. Such a large change in ζ would also significantly change the period of oscillation unless somehow balanced by a corresponding change in the magnetic field strength.

OBSERVATIONAL SIGNATURE OF LOOP EVOLUTION
In this section, we consider the evolution of our simulated coronal loop in terms of the established observational method to estimate the inhomogeneous layer width based on the shape of transverse intensity profile. The appearance of a coronal loop in optically thin EUV intensity images can be approximated using the square of the density integrated along the line-ofsight. This method was previously used by Aschwanden et al. (2003) and Aschwanden and Nightingale (2005) who modeled TRACE 171 Å data using a density profile with a sinusoidal inhomogeneous layer. It was extended by Aschwanden et al. (2007) to also consider rectangular (step) function and Gaussian density profiles, although these were found to typically be indistinguishable from the density profile with a sinusoidal inhomogeneous layer for TRACE data when including the smoothing effect of the instrumental point-spread function.
More recently, a similar method was applied by Pascoe et al. (2017b) to estimate the inhomogeneous layer width of a coronal loop observed using SDO/AIA 171 Å. The density profiles used in that study included a step function, Gaussian profile, and a inhomogeneous layer with a linear profile. Bayesian analysis found very strong evidence for the linear transition over the step profile (and very strong evidence for the step profile over the Gaussian), indicative of a small but finite inhomogeneous layer in the loop, and consistent with the separate seismological studies of the same loop by Pascoe et al. (2016bPascoe et al. ( , 2017a. A technique using both the seismological and intensity modeling methods simultaneously to estimate the inhomogeneous layer width of an oscillating coronal loop was produced by Pascoe et al. (2018). Goddard et al. (2017) performed a statistical study of 233 coronal loops comparing the step, Gaussian, and linear inhomogeneous layer density profiles. For 58% of cases there was strong evidence for the existence of an inhomogeneous layer, and there were no cases of strong evidence against there being an inhomogeneous layer. The inhomogeneous layer widths observed covered the full range of 0 to 2 and had no correlation with the loop radius. Figure 5 shows the value of inhomogeneous layer width inferred from our numerical simulation using the same method applied to observational data (Goddard et al., 2017;Pascoe et al., 2017b). This is calculated from the numerical data with the EUV intensity approximated as the square of the density integrated along the LOS, which we vary between 0 and 90 degrees relative to the plane of the loop oscillation. The observationally inferred value ǫ obs increases with time due to the non-linear evolution of the loop by development of KHI, as in Goddard et al. (2018), and with the LOS angle since the density perturbations associated with this evolution reflect the m = 1 symmetry of the kink mode driving it (see Figure 1).
Large LOS angles show the strongest variation in ǫ obs but would also exhibit the smallest transverse displacements. In such cases, the presence of kink oscillations could therefore be inferred indirectly by the increase of the inhomogeneous layer width in time, and directly using Doppler shift observations. For smaller LOS angles the increase in the inhomogeneous layer width is weaker but the transverse perturbations are larger, corresponding to the scenario for which decayless oscillations are detected by EUV imaging such as SDO/AIA. The estimate of ǫ presented in Figure 2 shows ǫ increasing to approximately 1 by the end of the simulation. This increase corresponds to an LOS angle of approximately 55 degrees in Figure 5. We can therefore expect that studies which consider the evolution of the loop coinciding with detectable oscillations (displacements) in EUV would typically underestimate the evolution which has actually taken place in terms of the evolution of the EUV intensity profile. This is in contrast to the results of the previous section which demonstrate how the seismologically-inferred value of ǫ would typically be an overestimate due to the generation of small scale density perturbations by KHI. Such a discrepancy between seismological and forward modeling estimates could therefore be indicative of non-linear evolution of the loop due to KHI, motivating further study to exclude other potential sources of discrepancy such as the influence of multiple loops along the LOS.
In Figure 6, we present a simple model to use the results of our numerical simulation to interpret an ensemble of observational data for a large number of loops such as the statistical study by Goddard et al. (2017). Those loops were chosen for a study of the transverse density profile with no consideration for oscillation. On the other hand, it has been established that decayless oscillations are ubiquitous in coronal loops (Anfinogentov et al., 2015). We can therefore consider that these loops will have undergone some range of non-linear evolution due to decayless oscillations, and that they also represent a range of LOS angles.
The top panel of Figure 6 shows a histogram of the observed value of the inhomogeneous layer width based on the study by Goddard et al. (2017). We note that in Figure 3 of that paper the values of ǫ correspond to the median value of the posterior probability distribution returned using Bayesian analysis with Markov chain Monte Carlo (MCMC) sampling (Anfinogentov et al., 2020). Here we also wish to consider the uncertainty associated with those values and so sum the full posterior probability distributions for all 233 loops. The study of each loop is based on 10 5 MCMC samples and the plotted histogram uses 100 bins hence the average histogram density being 233 × 10 5 /100 ∼ 200, 000.
Our model distribution for ǫ obs , shown as the dashed line in the top panel, is the combination of two probability distributions described below. We consider that for a loop with initial inhomogeneous layer ǫ 0 , the observed value ǫ obs may be any value greater than ǫ 0 (up to 2 by definition) depending on the amount of non-linear evolution and LOS angle as demonstrated in Figure 5. This is shown in the bottom panel for the example of ǫ = 0.3, with a uniform distribution assumed for simplicity. The middle panel represents the initial distribution of ǫ 0 . This will be determined by the mechanism which creates coronal loops and is presently unknown. We assume a half-normal distribution with the width of the distribution taken to be a free parameter which we fit by comparing our final model distribution to that from our observational data. The fitted width of approximately 0.5 implies loops are only required to have small initial inhomogeneous layer widths to account for the observed values, assuming the observational results are influenced by non-linear evolution similar to that shown in Figure 5.
The largest difference between the model and observational data is in the limit ǫ obs → 0. However, we note that the observational data is based on Bayesian analysis with a posterior probability density broadened to include values of 0 by observational error, whereas there is no consideration of errors in our model. The original statistical study by Goddard et al. (2017) found no cases for which the Bayesian evidence favored the step function density profile which corresponds to ǫ = 0. On the other hand, our assumption of a uniform probability distribution for ǫ obs under-represents the expected peak at ǫ 0 and especially penalizes the limit of ǫ 0 → 0 for which the range ǫ 0 ≤ ǫ obs ≤ 2 is greatest.
This simple model therefore demonstrates how observational results supporting larger values of inhomogeneous layer width are compatible with loops being generated with initially small values of ǫ under the assumption that the non-linear evolution of loops due to decayless oscillations acts only to increase the observed value of ǫ obs . This analysis also demonstrates how future observational studies might provide information about the mechanism by which coronal loops are generated if the probability distribution for ǫ 0 could be related to a theoretical model as opposed to our assumption of a half-normal distribution.

CONCLUSIONS
In this paper we have simulated the non-linear evolution of a coronal loop due to KHI and investigated the effect on the corresponding observational signatures. To facilitate analysis of the results our driver is impulsive which is typically associated with large amplitude decaying kink oscillations. On the other hand, the low amplitude of our driver corresponds to the decayless regime of kink oscillations and the repetitive perturbations may be considered as an approximation of some intermittent driving mechanism, for example, random buffeting of loop footpoints by photospheric motions. However, it is unlikely that decayless oscillations are actually generated in this way since it introduces abrupt changes in the phase of the oscillation which are typically not observed. Our resulting oscillation is decayless over the longer term, as kink mode energy lost by resonant absorption is replaced by subsequent additional perturbations. We find that these low amplitude kink oscillations approximating the decayless regime are capable of generating significant changes in the loop profile which can affect both the seismological and spatial estimates for the inhomogeneous layer width by different amounts. This can potentially be used to infer the presence of non-linear evolution of loop profiles through any disparity in estimates, although care would be required to ensure the direction and magnitude of the disparity is consistent with the physical model and observational conditions. Furthermore, we demonstrate that this effect could account for the observed distribution of inhomogeneous layer widths favoring larger values even if the mechanism by which loops are generated were to favor thin layers. This is similar to the argument by Magyar and Van Doorsselaere (2016b) for the nonexistence of multi-stranded coronal loops due to the mixing of plasma by KHI, but with the transverse plasma inhomogeneity being in the form of a boundary layer as opposed to a large number of small structures throughout the loop.
We can use the results of this paper to reconsider the case of a coronal loop which has been previously studied independently using seismological and spatial methods to determine to loop profile, and which exhibited both decayless and decaying kink oscillations. The seismological analysis of the loop, based on the large amplitude decaying oscillation, found an inhomogeneous layer width ǫ = 0.49 +0.23 −0.12 (Pascoe et al., 2017a) (the density contrast ratio was also estimated to be 3 as used in this paper). The spatial estimate based on forward modeling the transverse intensity profile produced a value of ǫ = 0.59 +0.14 −0.15 (Pascoe et al., 2017b). The loop was observed from the side, i.e., a LOS angle close to zero (see, e.g., Figure 1 of Pascoe et al., 2016b) and so a minimal increase in the width would be expected by studying the evolution of the intensity profile. The spatial study by Pascoe et al. (2017b) uses the loop intensity profile at the start of the large amplitude oscillation, however decayless oscillations were also detected in the same loop before, during and after the large amplitude oscillation; Nisticò et al. (2013) detected approximately 12 cycles of decayless oscillation before the onset of the large amplitude perturbation, while Pascoe et al. (2017a) estimated the amplitude of the decayless component as approximately 0.14 Mm (compared with a loop radius of approximately 3 Mm).
We can therefore consider two scenarios. Firstly, no KHI took place during the decayless regime, and so the spatial estimate of the inhomogeneous layer width is an accurate representation of the loop structure. Its consistency with the seismological estimate is due to the onset of any KHI arising from the large amplitude decaying oscillation being delayed by at least one cycle. Alternatively, KHI did take place during the decayless oscillations and in time to influence the results obtained for the large amplitude oscillation. The corresponding evolution of the loop profile means that the damping rate and inhomogeneous layer width would be overestimated by the seismological analysis, whereas the inhomogeneous layer width would be underestimated by the spatial analysis which has a very low LOS angle. The consistency of the two estimates for ǫ suggests the first scenario in which there was no KHI arising from the initial decayless stage of the oscillation. A recent study of the evolution of the temporal evolution of loops following large amplitude oscillations (Goddard and Nisticò, 2020) also found no clear evidence of KHI, but all loops analyzed did display changes after the onset of the oscillation. The extent to which this is due to observational issues, rather than genuine evolution of the loops, remains to be determined. That study was also based on DEM analysis which is known to produce overly broad spectra (e.g., Van Doorsselaere et al., 2018;Pascoe et al., 2019b) which may conceal the signatures of KHI.
The disparity between the abundant evidence for KHI in simulations compared to the lack of evidence in loop oscillations therefore remains to be accounted for. For the particular example considered above, the width of the loop being approximately three times larger than in our simulation may account for the apparent much slower evolution, with the non-linearity parameter being the displacement amplitude normalized by the loop radius (e.g., Van Doorsselaere et al., 2020). This reduces the 12 decayless cycles observed by Nisticò et al. (2013) to the equivalent of approximately 4 in our simulation, which is sufficiently early that the nonlinear evolution would not influence either the seismological or spatial results, and so they would remain consistent with each other and an accurate representation of the actual loop structure.
Our simulation demonstrates that the development of KHI is associated with an increasing damping rate and decreasing period of oscillation. However, the variation in the period of oscillation is small (approximately 3% during our entire simulation) and likely to be surpassed in observations by additional effects, such as changes in loop length or plasma properties, which can cause larger changes and either increases or decreases (e.g., Ireland and De Moortel, 2002;Nisticò et al., 2013;White et al., 2013;Pascoe et al., 2017c). Recent analysis of standing kink oscillations excited by two solar flares found the later oscillations damped faster than the earlier ones, consistent with KHI, while the period of oscillation increased by approximately 20% during the observation . The insensitivity of the period of oscillation to the transverse loop structure is consistent with the kink oscillation being a collective fast magnetoacoustic mode (e.g., Díaz et al., 2005;Pascoe et al., 2007). Previous studies have demonstrated that the damping of kink oscillations by resonant absorption is also not sensitive to loops being multi-stranded (Terradas et al., 2008b;Pascoe et al., 2011). However, the presence of such structuring is of interest since it can lead to a significant underestimation of the wave energy estimated using Doppler measurements Pant et al., 2019) or an overestimation of the broadness of DEM spectra . Due to LOS integration of optically thin EUV emission it is difficult to interpret fine structure in transverse intensity profiles, particularly when using a single channel as may be the case for coronal loops which are often identifiable only in the 171 Å channel of SDO/AIA. Observations using higher resolution instruments such as CRisp Imaging Spectro-Polarimeter (CRISP; Antolin et al., 2012;Scullion et al., 2014) and the High-Resolution Coronal Imager (Hi-C; Aschwanden and Peter, 2017;Williams et al., 2020) reveal that loop structuring indeed exists below the resolution of SDO/AIA. Our results demonstrate that the non-linear evolution of the loop profile by KHI produces a time-dependent (increasing) damping rate for kink oscillations that would not be expected for structuring due to loops being multi-stranded. This is due to KHI concentrating the changes to the loop profile in the inhomogeneous layer near the locations where resonant absorption occurs.
More generally, our results demonstrate the need for detailed comparison of non-linear numerical simulations and forward modeled observable signatures (preferably multichannel) to assist in the interpretation of physical effects such as KHI. The confirmation, or otherwise, of such effects allows additional physical conditions to be revealed, for example the inhibition of KHI due to the presence of magnetic twist (Terradas et al., 2018).

DATA AVAILABILITY STATEMENT
The raw data supporting the conclusions of this article will be made available by the authors, without undue reservation.

AUTHOR CONTRIBUTIONS
DP and TV contributed conception and design of the study. CG and DP performed the observational data analysis. DP wrote the first draft of the manuscript and received input from all co-authors.