© ESO, 2017

1. Introduction

Coronal loops are commonly subject to transverse oscillations (e.g. Aschwanden et al. 1999; Nakariakov et al. 1999; Verwichte et al. 2004; White & Verwichte 2012). The basic properties of observed transverse oscillations can be modelled by assuming homogeneous flux tube model (Nakariakov & Verwichte 2005). This approximation works reasonably well for predicting basic oscillation characteristics such as fundamental mode period; however for more complex oscillatory loop behaviour, the effects of the density variation within the loop should be considered.

The effect of longitudinal density stratification on the period and spatial structure of the fundamental kink mode has been studied for straight flux-tube coronal loop models (Andries et al. 2005a,b; Dymova & Ruderman 2005, 2006). For line-tied flux tubes, the effect of the density structuring on the fundamental kink mode properties is most significant when it is located near the loop apex and almost negligible in the case of the density structuring near the loop footpoints (Díaz et al. 2006).

The transverse loop oscillations have been observed to occur in two polarisations; horizontal (Aschwanden et al. 2002; Schrijver et al. 2002) and vertical (Wang & Solanki 2004; Mrozek 2011; White et al. 2013). The horizontally polarised oscillations are typically excited by a nearby flare (e.g. Aschwanden et al. 1999; Verwichte et al. 2004), while the source of excitation of the vertical oscillations remains unclear. They are usually explained by placing an external driver below the loop, by the reconnection process or by footpoint motion (Brady & Arber 2005; Gruszecki & Murawski 2008; Selwa et al. 2010; White et al. 2013). Subsonic siphon flows have also been proposed as a possible excitation mechanism (Ofman et al. 2012).

Recent observations show that many coronal loops are subject to thermal instability undergoing catastrophic cooling (Antolin & Rouppe van der Voort 2012). This leads to the formation of coronal rain which consists of numerous cool and dense plasma condensations formed near the top of the unstable loop. These condensations then fall towards the solar surface under the influence of gravity, often with sub-ballistic speeds due to effects of plasma pressure and magnetic field (Mackay & Galsgaard 2001; Oliver et al. 2014; Kohutova & Verwichte 2017). The thermally unstable loops are therefore likely to contain an overdense region in the part of the loop where the condensation occurs which is typically near the loop apex. Small amplitude transverse oscillations are often present in rainy coronal loops (Antolin & Verwichte 2011; Kohutova & Verwichte 2016; Verwichte & Kohutova 2017), with the observed loops in question typically tilted with respect to the photospheric normal.

A one-dimensional (1D) mechanical model for the excitation of small amplitude transverse oscillations by the coronal rain condensations via their concentrated mass was proposed by Verwichte et al. (2017). These authors also showed that small amplitude transverse oscillations observed by Kohutova & Verwichte (2016) could be excited by the coronal rain. In the special case of a coronal loop lying in a plane perpendicular to solar surface, this mechanism would lead to the excitation of oscillation polarised in vertical direction only.

Fig. 1 Left: initial density configuration of the coronal loop with cool condensation region for μ = 0.7. The white lines show the B-field direction. Right: Alfvén speed profile along the centre of the coronal loop.

The first observational evidence of excitation and evolution of vertical coronal loop oscillations caused by catastrophic cooling and coronal rain formation was shown in Verwichte & Kohutova (2017). The observations presented therein also showed a change of the oscillation period over the duration of the loop evolution caused by the drainage of the coronal loop mass by the coronal rain.

We carried out 2.5D magnetohydrodynamic (MHD) simulations of a coronal loop with cool plasma condensation region at the loop apex. This is found to excite sustained, small amplitude, vertically polarised transverse loop oscillations. We analysed the oscillation parameters as a function of condensation region mass. We further calculated the expected change in the period of the fundamental mode due to changing density profile along the loop. This change is a result of the condensations falling down towards the loop footpoints under the influence of gravity.

2. Numerical model

We solved the nonlinear MHD equations using Lare2d (Arber et al. 2001) assuming perfectly ionised, fully conductive plasma and using the ideal equation of state. We included the effect of gravity and shock viscosity. The role of shock viscosity is to introduce dissipation at strong gradients only, preventing false oscillations behind the shock that would arise from applying second order accurate finite difference scheme when the pressure difference exceeds the jump condition. Thermal conduction and radiative transport terms were not included in the energy equation. The equations were solved on a square, uniform 1024 × 1024 grid with the extent − 100 Mm Mm in horizontal direction and 0 Mm Mm in the vertical direction. Grid convergence study has been carried out in order to check the convergence of the numerical results. We set up the problem so as to represent a relatively long coronal loop embedded in a magnetic arcade which contains a high-density, low-temperature condensation region at its apex (Fig. 1). The gravity was assumed to be uniform and along negative y-direction. Both the coronal loop and the ambient plasma are gravitationally stratified.

The equilibrium magnetic field is given by current-free magnetic arcade model determined by the potential ${A}\mathrm{=}\mathrm{-}{\mathit{B}}_{\mathrm{0}}{\mathit{H}}_{\mathrm{B}}\cos \mathrm{(}\mathit{x}\mathit{/}{\mathit{H}}_{\mathrm{B}}\mathrm{)}{\mathrm{e}}^{\mathrm{(}\mathrm{-}\mathit{y}\mathit{/}{\mathit{H}}_{\mathrm{B}}\mathrm{)}}\hat{{y}}$ (Priest 1982) such that the magnetic field components are given by: $\mathrm{\left(}{\mathit{B}}_{\mathit{x}}\mathit{,}{\mathit{B}}_{\mathit{y}}\mathit{,}{\mathit{B}}_{\mathit{z}}\mathrm{\right)}\mathrm{=}{\mathit{B}}_{\mathrm{0}}\mathrm{\left[}\cos \mathrm{\right(}\mathit{x}\mathit{/}{\mathit{H}}_{\mathrm{B}}\mathrm{)}\mathit{,}\mathrm{-}\sin \mathrm{(}\mathit{x}\mathit{/}{\mathit{H}}_{\mathrm{B}}\mathrm{\left)}\mathit{,}\mathrm{0}\mathrm{\right]}{\mathrm{e}}^{\mathrm{-}\mathit{y}\mathit{/}{\mathit{H}}_{\mathrm{B}}}\mathit{,}$(1)where H_B is the magnetic scale height given by H_B = W/π with W = 200 Mm being the horizontal extend of the arcade and B₀ = 70 G is the magnetic field at y = 0. This results in the magnetic field of ~20 G at coronal heights, which is representative of real coronal conditions.

We assumed the background temperature to be constant in the x-direction and in the y-direction we create a smoothed step function temperature profile representative of an atmosphere consisting of a cool chromosphere, transition region layer and hot corona (Cargill et al. 1997): $\mathit{T}\mathrm{\left(}\mathit{y}\mathrm{\right)}\mathrm{=}\frac{\mathrm{1}}{\mathrm{2}}\mathrm{(}{\mathit{T}}_{\mathrm{cor}}\mathrm{+}{\mathit{T}}_{\mathrm{ph}}\mathrm{)}\mathrm{+}\frac{\mathrm{1}}{\mathrm{2}}\mathrm{(}{\mathit{T}}_{\mathrm{cor}}\mathrm{-}{\mathit{T}}_{\mathrm{ph}}\mathrm{)}\tanh (\frac{\mathit{y}\mathrm{-}{\mathit{y}}_{\mathrm{t}}}{\mathrm{\Delta }\mathit{y}})\mathit{,}$(2)with photospheric temperature T_ph = 6 × 10³ K, coronal temperature T_cor = 10⁶ K, y_t = 4 Mm and Δy = 1 Mm. The temperature variation controls the pressure scale height Λ(y): $\mathrm{\Lambda }\mathrm{\left(}\mathit{y}\mathrm{\right)}\mathrm{=}\frac{{\mathit{k}}_{\mathit{b}}\mathit{T}\mathrm{\left(}\mathit{y}\mathrm{\right)}}{\mathit{m}{\mathit{g}}_{\mathrm{\odot }}}\mathit{,}$(3)where g_⊙ = 274 m s^-2 is the average solar surface gravity and m is the mean particle mass. The density profile for the non-isothermal stratified atmosphere is then determined by numerically solving for a hydrostatic pressure balance: $\mathit{p}\mathrm{\left(}\mathit{y}\mathrm{\right)}\mathrm{=}{\mathit{p}}_{\mathrm{0}}\exp (\mathrm{-}{\mathrm{\int }}_{\mathrm{0}}^{\mathit{y}}\frac{\mathrm{d}{\mathit{y}}^{\mathrm{\prime }}}{\mathrm{\Lambda }\mathrm{\left(}{\mathit{y}}^{\mathrm{\prime }}\mathrm{\right)}})\mathit{,} \mathit{\rho }\mathrm{\left(}\mathit{y}\mathrm{\right)}\mathrm{=}\frac{\mathit{mp}\mathrm{\left(}\mathit{y}\mathrm{\right)}}{{\mathit{k}}_{\mathit{b}}\mathit{T}\mathrm{\left(}\mathit{y}\mathrm{\right)}}\mathrm{·}$(4)The coronal loop was modelled as a density enhancement along a magnetic field line defined by: ${\mathit{y}}_{\mathrm{L}}\mathrm{\left(}\mathit{x}\mathrm{\right)}\mathrm{=}\frac{\mathrm{1}}{\mathrm{2}}{\mathit{H}}_{\mathrm{B}}[\log (\frac{\cos \mathrm{(}\mathit{x}\mathit{/}{\mathit{H}}_{\mathrm{B}}\mathrm{)}}{\cos \mathrm{(}\mathit{h}\mathit{/}{\mathit{H}}_{\mathrm{B}}\mathrm{)}})\mathrm{+}\log (\frac{\cos \mathrm{(}\mathit{x}\mathit{/}{\mathit{H}}_{\mathrm{B}}\mathrm{)}}{\cos \mathrm{\left(}\mathrm{\right(}\mathit{h}\mathrm{-}\mathit{a}\mathrm{)}\mathit{/}{\mathit{H}}_{\mathrm{B}}\mathrm{)}})]\mathit{,}$(5)where h = 90 Mm and a = 3 Mm is the loop scale width. The density variation between the loop and the background medium in is given by the symmetric Epstein profile (Nakariakov & Roberts 1995): $\mathit{\rho }\mathrm{\left(}\mathit{x,y}\mathrm{\right)}\mathrm{=}{\mathit{\rho }}_{\mathrm{e}}\mathrm{\left(}\mathit{x,y}\mathrm{\right)}\mathrm{+}\mathrm{\left(}{\mathit{\rho }}_{\mathrm{i}}\mathrm{\right(}\mathit{x,y}\mathrm{)}\mathrm{-}{\mathit{\rho }}_{\mathrm{e}}\mathrm{(}\mathit{x,y}\mathrm{\left)}\mathrm{\right)}{\mathrm{sech}}^{\mathrm{2}}(\frac{\mathrm{(}\mathit{y}\mathrm{-}{\mathit{y}}_{\mathrm{L}}\mathrm{(}\mathit{x}\mathrm{\left)}\mathrm{\right)}}{\mathit{a}})\mathit{,}$(6)where ρ_e and ρ_i are external and internal densities respectively. We assumed a constant density contrast χ = ρ_i/ρ_e = 10 along the whole loop. The density stratification was calculated using base density ρ₀ = 5 × 10^-8 kg m^-3, resulting in the densities in the upper half of the loop of the order of 10^-11 kg m^-3, that is, representative of coronal values.

The plasma condensation region was represented as density enhancement superimposed on top of the equilibrium density profile, accompanied by the corresponding low temperature region necessary to maintain the plasma pressure balance and to prevent rapid initial expansion. The two-dimensional (2D) Gaussian density enhancement is positioned at x₀ = 0 Mm and y₀ = y_L(x₀) and has widths σ_x = 7 Mm in horizontal direction and σ_y = 1 Mm in vertical direction and height ρ_blob = r_bcρ_bg(x₀,y₀), with r_bc being the density contrast between the peak condensation density and the density ρ_bg of the hot background plasma at the loop axis at the same position in height. The bottom domain boundary was fixed to create a line-tied loop and the boundary conditions along the remaining three boundaries are symmetric (i.e. gradients set to zero).

Without the presence of the condensation region, the whole system is in near-equilibrium state (it should be noted that we are not explicitly forcing plasma pressure balance inside and outside of the coronal loop, this is however negligible, since the plasma-β is small and the magnetic pressure dominates, preventing loop expansion). The finite mass of the condensation region, however, displaces the axis of the coronal loop downwards which triggers vertically polarised transverse loop oscillations.

3. Vertical loop oscillations

The excitation of the vertical oscillations is studied for six different values of peak density of the condensation region, with the density contrast ranging from r_bc = 10 (ρ = 9.2 × 10^-11 kg m^-3) to r_bc = 200 (ρ = 1.8 × 10^-9 kg m^-3). This corresponds to fraction of rain mass relative to total loop mass including rain ranging from μ = 0.19 to μ = 0.83.

Initially, the loop top is displaced downwards by the mass of the condensation region, resulting in onset of a vertically polarised transverse loop oscillation. For low condensation region densities (μ = 0.19 to μ = 0.54) the cool plasma is confined to the loop top by combination of the pressure from the underlying plasma and the ponderomotive force resulting from the loop oscillations directed towards the antinode of the oscillation, that is, the loop apex for the fundamental harmonic (Verwichte et al. 2017). The amplitude and period of the loop oscillation stays approximately constant for the duration of the simulation. For the high condensation region densities (μ = 0.62 to μ = 0.83) the cool plasma falls towards the loop footpoints along both sides of the loop similar to a real coronal rain shower scenario. This results in long period, large amplitude oscillations transitioning into persistent, smaller amplitude, shorter period oscillations once the condensations have fallen into the lower loop legs, where their mass has negligible effect on the loop oscillation period. The downfalling plasma blobs elongate as they fall due to the differential effective gravity component acting along the finite length of the blob. The blobs develop elongated wings, leaving an evacuated region behind. The development of this v-shape is likely caused by the fact that the motion of the blob is subsonic in the background plasma but supersonic for its own temperature, leading to shock behaviour. The blobs are then gradually slowed down and ultimately rebound multiple times in the lower loop leg due to the combined effect of the underlying plasma pressure and magnetic tension force counteracting the effective gravity (Mackay & Galsgaard 2001). This type of oscillatory motion of plasma condensations in longitudinal direction was studied in detail in Kohutova & Verwichte (2017).

Fig. 2 Top: ripples due to Rayleigh-Taylor instability formed at the interface between the condensation region and ambient coronal plasma near at the top of the loop at t = 5544 s. Bottom: subsequent fragmentation of tails of downfalling plasma condensations at t = 6720 s.

The density difference between the condensation region and the background plasma leads to onset of the Rayleigh-Taylor instability at the loop apex. Ripples of cool plasma are formed in the condensation region and later lead to fragmentation of the tails of the downfalling plasma blobs (Fig. 2). A similar phenomenon can be seen in coronal rain observations, which typically show a lot of structuring in the longitudinal direction, resulting in a clumpy appearance of the downfalling coronal rain material (e.g. Antolin & Rouppe van der Voort 2012; Antolin et al. 2015). Small fraction of the cool condensation plasma remains confined at the loop apex. This is caused by the action of the ponderomotive force pushing the plasma in the vicinity of the loop top towards the apex point where the effective gravity is zero. As the simulation setup is symmetric, the plasma remains trapped here for the duration of the simulation.

In order to analyse the evolution of the oscillation in detail, at each timestep we take a cut along the domain centre perpendicular to the loop axis to create density time-distance plots. The loop displacement time series is obtained by fitting a Gaussian to the density profile at each timestep. For the two peak profile seen at the higher condensation densities due to evacuation of the loop centre this approach leads to tracking one of the loop edges. Damped sine function of the form ξ(t) = ξ₀exp( − t/τ)sin(ωt + φ) is fitted to the loop top displacement timeseries at the beginning and end of the simulation to obtain oscillation periods and amplitudes shown in Table 1. The periods and amplitudes of the oscillation during the initial phase increase with increasing mass of the condensation region (Fig. 3). After most of the condensation mass has fallen from the loop top towards the loop footpoints, the oscillation period decreases to less than half of the original value. This is a result of change in the longitudinal density distribution along the loop, even though the total mass of the plasma in the loop remains conserved. A small degree of periodic modulation of the period can be observed during the latter half of the simulation as the plasma blobs oscillate in lower loop legs. The oscillation displacement profile of the fundamental harmonic also evolves with time, starting with a concentrated peak at the loop top and gradually broadening into a sine-like profile as the dense plasma moves towards loop footpoints (Fig. 4).

Fig. 3 Dependence of amplitudes (left) and periods (right) of the vertical oscillation on the rain mass fraction μ at the beginning (ξ01, P1) and end (ξ02, P2) of the simulation. Values for low mass condensation regions for which the plasma remains trapped at the loop apex are shown in red. The solid line shows amplitudes predicted by 1D mechanical model by Verwichte et al. (2017).

Using an expression from Verwichte et al. (2017) for the amplitude of a vertical oscillation excited by a dense plasma blob ξ₀ based on a 1D mechanical coronal rain model: ${\mathit{\xi }}_{\mathrm{0}}\mathrm{=}\mathrm{(}\mathrm{5.0}\mathrm{\pm }\mathrm{0.6}\mathrm{)}\mathrm{\times }{\mathrm{10}}^{-3}\sqrt{\mathrm{1}\mathrm{-}\frac{\mathrm{2}\mathit{\theta }}{\mathit{\pi }}}\mathit{\mu L,}$(7)where θ is the angle between the loop plane and the photospheric normal, μ is the rain mass fraction and L is the length of the loop we estimate expected oscillation amplitudes for all condensation masses. We find that amplitude values predicted by this model agree with the oscillation amplitudes determined from the simulation for lowest masses of the condensation region, while at higher masses the two diverge and the 1D mechanical model underestimates the oscillation amplitude (Fig. 3). Here we note that the simulated scenario does not comply with number of assumptions made in the 1D mechanical model, namely the assumption of small size of the plasma blobs and the constant loop length. The size of the condensation region used in the simulation is extended in the longitudinal direction, that is, along the loop axis and therefore cannot be approximated by a point-like mass. Similarly the coronal loop is significantly deformed in the simulation, thus violating the constant loop length assumption. The applicability of the 1D kinematic model to the cases with high condensation region masses studied here is therefore limited. It is however encouraging that the two models agree within their common range of validity at low condensation region masses. Furthermore, these results not only confirm that the presence of coronal rain is a possible excitation mechanism of transverse loop oscillations but also show that its effect is even greater than originally predicted.

Fig. 4 Loop axis displacement from its original position as a function of time and position along the loop for μ = 0.7.

Fig. 5 Time-distance plots at the loop apex for different values of condensation region masses. White solid lines show the centre of the loop profile determined by Gaussian fitting. Black dotted lines show best-fit damped sine function.

We further determine the dependence of the loop oscillation periods on the longitudinal density profile of the loop semi-analytically using a 1D model. The vertical oscillation of a loop with longitudinal density dependence can be modelled in the thin flux tube limit using the equation (Dymova & Ruderman 2005, 2006): $\frac{{\mathrm{d}}^{\mathrm{2}}\mathit{\xi }}{\mathrm{d}{\mathit{s}}^{\mathrm{2}}}\mathrm{+}\frac{{\mathit{\omega }}^{\mathrm{2}}}{{\mathit{C}}_{\mathrm{k}}^{\mathrm{2}}}\mathit{\xi }\mathrm{=}\mathrm{0}\mathit{,}$(8)where ξ(s) is the loop displacement, ω is the mode frequency, ${\mathit{C}}_{\mathrm{k}}^{\mathrm{2}}\mathrm{=}\mathrm{2}{\mathit{B}}^{\mathrm{2}}\mathit{/}{\mathit{\mu }}_{\mathrm{0}}\mathrm{\left(}{\mathit{\rho }}_{\mathrm{i}}\mathrm{\right(}\mathit{s}\mathrm{)}\mathrm{+}{\mathit{\rho }}_{\mathrm{e}}\mathrm{(}\mathit{s}\mathrm{\left)}\mathrm{\right)}$ and s is the coordinate along the loop. To model a line-tied coronal loop, ξ(s) must satisfy the boundary conditions ξ(0) = ξ(L) = 0. Equation (8) is derived for the case of the uniform magnetic field. In the case of magnetic field variation in the longitudinal direction ξ(s) in Eq. (8) modifies to $\mathit{\xi }\mathrm{\left(}\mathit{s}\mathrm{\right)}\mathit{/}\sqrt{\mathit{B}\mathrm{\left(}\mathit{s}\mathrm{\right)}}$ (Ruderman et al. 2008). This change of the normalisation however only affects the displacement profile of the oscillation; the oscillation period remains unchanged.

Fig. 6 Internal density profiles used for the calculation of the oscillation period in thin flux tube approximation for the μ = 0.7 case. Solid and dashed line show initial and final state of the loop respectively.

Equation (8) for the thin flux tube model is solved numerically, with the mode frequency subsequently determined using a shooting method. Internal and external plasma density and magnetic field strength profiles along the loop are determined using the same magnetic field and density configuration as described in Sect. 2 and by averaging along the transverse direction of the loop. We solve the equation for all condensation region masses and for each mass we use two different density profiles; one corresponding to the initial state of the loop with density enhancement located at the loop apex, and one with the density enhancement split into two components each located in one loop leg, representative of the final state of the loop, where the condensation mass has fallen down under the influence of gravity (Fig. 6). The resulting initial and final periods are shown in Table 1. For all condensation region masses, the thin flux tube model is found to overestimate the periods of the fundamental oscillation mode by about a factor of two compared to the periods determined from the simulation. The simulation setup used in this work is equivalent to modelling the loop as a curved magnetic slab rather than thin flux tube model used in the above calculation. We therefore attribute this difference to a different geometry used by the two models. The initial to final period ratios however agree with the simulation results for the higher condensation region masses (μ = 0.62 to μ = 0.83). We refrain from calculating the final state periods (i.e. with two density enhancements in the loop legs) for the low condensation region masses (μ = 0.19 to μ = 0.54) as in these cases there is no significant mass redistribution along the length of the loop during the simulation.

We further proceed to 1D calculation done for the case of a magnetic slab. The main caveat in this approach is that in the slab case, the kink speed is not well defined and the phase speed of a fast kink mode in a slab instead tends to external Alfvén speed in the long wavelength limit (Edwin & Roberts 1982). External Alfvén speed of the loop in our setup however does not contain the information about the high density region. In order to incorporate the density variation inside the loop into a slab model, C_k in Eq. (8) is replaced by ${\mathit{V}}_{\mathrm{Ae}}\mathrm{\left(}\mathit{s}\mathrm{\right)}\mathrm{=}\mathit{B}\mathrm{\left(}\mathit{s}\mathrm{\right)}\mathit{/}\sqrt{\mathit{\mu }{\mathit{\rho }}_{\mathrm{i}}\mathrm{\left(}\mathit{s}\mathrm{\right)}\mathit{/}\mathit{\chi }}$, where ρ_i is the density along the centre of the loop and χ = 10 as specified in Sect. 2. This is not strictly correct as it assumes constant density contrast between the internal and external density profile which is equivalent to assuming that the structure of the external medium self-adjusts such that it reflects the density structuring along the centre of the loop. This is further supported by the fact that both the oscillation periods and period ratios calculated using this approach diverge from the simulation values with increasing mass of the condensation region (Table 1).

4. Discussion and conclusions

We presented 2.5D MHD simulations of evolution of coronal loop with a cool and dense condensation region located at the loop apex. The mass of the condensation region was found to excite vertically polarised transverse loop oscillations. We analysed the dependence of the amplitudes and periods of the excited oscillations as a function of condensation region mass. The observed shift in the fundamental kink mode period of the loop following a coronal rain shower seen in the simulation agrees with the period shift described in Verwichte & Kohutova (2017), where it was explained as a consequence of the drainage of the total loop mass by the coronal rain into the chromosphere. In this work we have shown that period shifts can also be explained by redistribution of the mass along the length of the coronal loop as the cool dense plasma falls from the top of the loop.

Long term confinement of the plasma at the loop apex seen in the simulation is not likely for real loops undergoing condensation formation due to lack of symmetry being common in these cases. This lack of symmetry can be due to a number of reasons such as the asymmetric nature of footpoint heating mechanism or due to siphon flows caused by the pressure difference between the footpoints. It has indeed been observed that the latter can lead to shifting the condensation region away from the loop apex such that the coronal rain is seen to fall along one loop leg only (Kohutova & Verwichte 2016).

The development of Rayleigh-Taylor instability along the interface between the dense condensation plasma and coronal plasma at the bottom edge of the loop apex resulting from the density difference between the two was seen during early stages of the loop evolution. This resulted in a formation of ripples along the interface and subsequently lead to fragmentation of the tails of the downfalling condensations. The Rayleigh-Taylor instability has been observed to develop in prominence-corona transition region in cases where dense prominence plasma is suspended above a coronal cavity (Berger et al. 2010; Hillier et al. 2012). It is however unclear if such scenario is realistic in coronal loops undergoing condensation formation, given the realistic density estimates and taking into account typical magnitudes and orientation of the magnetic field (which at the loop top is parallel to the fluid interface). The fragmentation of the condensed plasma seen in our simulation resulting from the Rayleigh Taylor instability is however in line with observed clumpy structure of the coronal rain.

The approach used here imposes the existence of the condensation region as an initial condition and hence the density enhancement is introduced abruptly, which is equivalent to assuming that the catastrophic cooling leading to condensation formation occurs instantaneously. Here we point out that in the real scenario, the condensation formation process will be gradual and the coronal loop can be expected to readjust to some extent while the condensation region is being formed. This will limit the effect of the density enhancement on the oscillatory behaviour of the loop. We can estimate whether the vertical oscillations will still be triggered during the gradual formation process by comparing the characteristic radiative cooling time scale with the typical oscillation amplitude. The radiative cooling time scale τ_rad can be estimated using the RTV radiative loss function: $\begin{array}{ccc}{\mathit{n}}_{\mathrm{e}}^{\mathrm{2}}{\mathit{\chi }}_{\mathrm{R}}{\mathit{T}}^{\mathit{\alpha }}\mathrm{\simeq }\frac{{\mathit{k}}_{\mathrm{b}}}{\mathrm{(}\mathit{\gamma }\mathrm{-}\mathrm{1}\mathrm{)}}{\mathit{n}}_{\mathrm{e}}\frac{\mathit{T}}{{\mathit{\tau }}_{\mathrm{rad}}}\mathit{,}& & \end{array}$(9)where n_e is the plasma electron density, γ is the ratio of specific heats, χ_R and α are radiative loss function coefficients as given in Rosner et al. (1978). Using the typical coronal values T = 10⁶ K and ρ = 10^-11 kg m^-3 as employed in the simulation, we obtain τ_rad = 300 s, i.e. similar but less than the typical oscillation periods of 400−900 s. This is the higher bound on the τ_rad since as the condensation forms, the density increases from coronal to chromospheric values and the characteristic radiative cooling time scale will further decrease. In this particular scenario it can be concluded that the condensation formation is abrupt enough for the loop oscillations to be triggered. In addition, we also use the temperature and density values determined from the observations in Verwichte & Kohutova (2017) and obtain τ_rad ~ 500 s as compared to the typical oscillation period of ~150 s observed therein, i.e. the same order of magnitude but larger. Here it should be pointed out that the determination of the plasma electron density from the observations depends on the assumption of line-of-sight depth of the loop and is therefore subject to large uncertainties. The oscillations and catastrophic cooling can therefore be expected to occur on similar time scales, exact balance of which is subject to specific conditions in the loop.

A decrease in the oscillation period similar to the one studied in this work was found in cooling loops (Morton & Erdélyi 2009; Ruderman 2011; Magyar et al. 2015). Both have the same fundamental cause, the change in the longitudinal density profile of the loop, namely the redistribution of the loop plasma from the loop top to the loop footpoints. However, the processes leading to the mass redistribution are different in the two cases; in our case the cause are downfalling condensations, modelled as a result of catastrophic cooling process with short timescale τ_cool. The change in period due to downfalling condensations is therefore abrupt and occurs on a timescales of ~100 s, that is, the time it takes for the condensations to move significantly under the influence of gravity. In the latter scenario it is the gradual flow of mass from the loop top to the loop points caused by temperature decrease leading to scale-height decrease as a result of gradual cooling process with long τ_cool. The period change due to loop cooling obviously depends on the cooling timescale, but is usually more gradual and in the above studies occurs on timescales of ~1000 s.

A natural extension of the work presented here is to model the coronal rain self-consistently, that is, model the whole process of thermal instability onset and coronal rain formation starting with footpoint-heated stratified loop and look for signatures of vertical oscillations. The excitation of the vertical oscillations will then depend on the fraction of the total loop mass that condenses to form coronal rain and on the spatial distribution of the condensed plasma as well as on the temporal evolution of the condensation process. The approach used in this work however enables us to carry out a parameter study. Imposing the existence of the condensation region as an initial condition gives us control over its properties, which makes it easier for us to understand the conditions under which the excitation of the vertical oscillations is significant and how the condensation mass affects the periods and amplitudes of the excited oscillations. The self-consistent scenario including the coronal rain formation will be addressed in a follow-up paper.

Acknowledgments

P.K. acknowledges the support of the UK STFC Ph.D. studentship. E.V. acknowledges financial support from the UK STFC on the Warwick STFC Consolidated Grant ST/L000733/I.

References

All Tables

All Figures

	Fig. 1 Left: initial density configuration of the coronal loop with cool condensation region for μ = 0.7. The white lines show the B-field direction. Right: Alfvén speed profile along the centre of the coronal loop.
In the text

	Fig. 2 Top: ripples due to Rayleigh-Taylor instability formed at the interface between the condensation region and ambient coronal plasma near at the top of the loop at t = 5544 s. Bottom: subsequent fragmentation of tails of downfalling plasma condensations at t = 6720 s.
In the text

	Fig. 3 Dependence of amplitudes (left) and periods (right) of the vertical oscillation on the rain mass fraction μ at the beginning (ξ01, P1) and end (ξ02, P2) of the simulation. Values for low mass condensation regions for which the plasma remains trapped at the loop apex are shown in red. The solid line shows amplitudes predicted by 1D mechanical model by Verwichte et al. (2017).
In the text

	Fig. 4 Loop axis displacement from its original position as a function of time and position along the loop for μ = 0.7.
In the text

	Fig. 5 Time-distance plots at the loop apex for different values of condensation region masses. White solid lines show the centre of the loop profile determined by Gaussian fitting. Black dotted lines show best-fit damped sine function.
In the text

	Fig. 6 Internal density profiles used for the calculation of the oscillation period in thin flux tube approximation for the μ = 0.7 case. Solid and dashed line show initial and final state of the loop respectively.
In the text