SIMULATIONS OF THE KELVIN–HELMHOLTZ INSTABILITY DRIVEN BY CORONAL MASS EJECTIONS IN THE TURBULENT CORONA

Shear flows are ubiquitous in astrophysical problems, such as jet propagation in the interstellar medium (Ferrari et al. 1980; Begelman et al. 1984; Bodo et al. 1994), the dynamics of spiral arms in galaxies (Dwarkadas & Balbus 1996), cometary tails (Ershkovich et al. 1973; Brandt & Mendis 1979), and differential rotation in accretion disks (Balbus & Hawley 1998). It is also relevant in a variety of space physics problems, such as zonal flows in the atmospheres of rotating planets like Jupiter (Hasegawa 1985), the solar wind (Poedts et al. 1998), or the Earth's magnetopause (Parker 1958).

Shear flows often give rise to the well-known Kelvin–Helmholtz (KH) instability (Helmholtz 1868; Thomson 1871), which has been extensively studied by Chandrasekhar (1961). It is an ideal hydrodynamic instability, which converts the energy of large-scale velocity gradients into kinetic and/or magnetic energy at much smaller scales. The presence of a magnetic field component parallel to the shear flow has a stabilizing effect, and can even stall the instability if the parallel component of the Alfvén velocity becomes larger than one-half of the shear velocity jump (Lau & Liu 1980; Miura & Pritchett 1982). A similar instability condition was anticipated by Ershkovich et al. (1973) in connection with observational evidence of KH in comet tails. On the other hand, an external magnetic field pointing in any direction perpendicular to the shear flow has no effect on the linear regime of the instability, and it is simply advected by the flow.

The first observations of a KH pattern in the solar corona were reported by Foullon et al. (2011) for the 2010 November 3 event using data from the Atmospheric Imaging Assembly (AIA) on board the Solar Dynamics Observatory (SDO). Ofman & Thompson (2011) also reported observations of a KH pattern obtained by AIA/SDO for the 2010 April 8 event. AIA produces high spatial resolution (pixel size of 0.6 arcsec) and high temporal cadence (10–20 s) images of the Sun in several bandpasses covering white light, ultraviolet, and extreme ultraviolet. The observed pattern of the KH instability observed by Foullon et al. (2011) extends from about 70 Mm up to about 180 Mm above the solar surface (1 Mm = 10³ km). When a coronal mass ejection (CME) expands supersonically upwards from the solar surface, a bow shock is formed ahead of the CME and a strong shear flow develops across the contact discontinuity, separating the shocked ambient plasma from the ejected material. A similar configuration arises at the flanks of the Earth's magnetopause, where the KH instability has also been observed and studied (Fujimoto & Teresawa 1995; Fairfield et al. 2000; Nykyri & Otto 2001). More recently it was observed in connection to the magnetopause of other planets, such as Saturn (Masters et al. 2010) and Mercury (Sundberg et al. 2011). When the supersonic solar wind impinges on these magnetized planets, it first crosses a bow shock (and becomes subsonic in the reference frame of the planet) and then circumvents the planet, slipping through the outer part of a surface of tangential discontinuity known as the magnetopause, where a strong shear flow is generated.

The ambient corona is expected to be turbulent, as evidenced by measurements of nonthermal broadenings of highly ionized spectral lines. Most recent observations of nonthermal broadenings obtained by the Extreme-ultraviolet Imaging Spectrometer (EIS) on board Hinode correspond to nonthermal motions in the range of 20–60 km s⁻¹ (Doschek et al. 2014). The typical sizes of these nonthermal motions are sufficiently small to remain unresolved by EIS, whose pixel size is 2 arcsec, and therefore its only manifestation is an excess in the Doppler broadening of spectral lines (i.e., beyond the thermal Doppler broadening).

The goal of the present paper is to study the interaction between these two rather dissimilar features: the large-scale laminar pattern generated by the ongoing KH instability, and the small-scale nonthermal motions presumably corresponding to a well-developed turbulence. With this goal in mind, we set up three-dimensional simulations of the MHD equations to study the evolution of the KH instability in the presence of a turbulent ambient background. Nykyri & Foullon (2013) presented results from a large number of compressible 2.5D MHD simulations (without a turbulent background) for parameter values compatible with the observations of the 2010 November 3 event. This comparison is consistent with a magnetic field almost perpendicular to the flow plane, and therefore we make this assumption in our simulations. When a small-scale turbulent background is considered, the expected role of a large-scale flow is to produce the effect of an enhanced diffusivity, which can be characterized through an effective or turbulent viscosity. The effect of this extra diffusivity on an ongoing instability for the large-scale flow, as is currently the case for KH, is to reduce its growth rate or even to switch off the instability completely. We test and basically confirm this theoretical picture with a series of simulations of a KH-unstable shear flow superimposed on a small-scale turbulent background with different turbulence intensities. The AIA observations showing a KH pattern are described in Section 2.1 and the observed features of small-scale turbulence are summarized in Section 2.2. We introduce the MHD equations in Section 3 and describe the basic properties of the KH instability in Section 4. The characteristic features of the turbulent background generated in our simulations are discussed in Section 5 and our numerical results are shown in Section 6. The potential consequences of the results presented in this paper are discussed in Section 7, and our conclusions are listed in Section 8.

2.1. AIA Observations

2.2. Turbulent Broadening

The incompressible MHD equations for a fully ionized hydrogen plasma are the Navier–Stokes equation and the induction equation:

$$\begin{eqnarray}\displaystyle \frac{\partial {\boldsymbol{U}}}{\partial t} & = & -({\boldsymbol{U}}\cdot {\boldsymbol{\nabla }}){\boldsymbol{U}}+{v}_{{\rm{A}}}^{2}({\boldsymbol{\nabla }}{\boldsymbol{\times }}{\boldsymbol{B}}){\boldsymbol{\times }}{\boldsymbol{B}}\\ & & -{\boldsymbol{\nabla }}P+\nu {{\rm{\nabla }}}^{2}{\boldsymbol{U}}+{\boldsymbol{F}}\end{eqnarray} \tag{ 1 }$$

$$\begin{eqnarray}&&\displaystyle \frac{\partial {\boldsymbol{B}}}{\partial t}={\boldsymbol{\nabla }}{\boldsymbol{\times }}({\boldsymbol{U}}{\boldsymbol{\times }}{\boldsymbol{B}})+\eta {{\rm{\nabla }}}^{2}{\boldsymbol{B}}.\end{eqnarray} \tag{ 2 }$$

The velocity ${\boldsymbol{U}}$ is expressed in units of a characteristic speed U₀, the magnetic field ${\boldsymbol{B}}$ is in units of B₀, and we also assume a characteristic lengthscale L₀ and a spatially uniform particle density n₀. In general terms, the assumption of incompressibility is valid provided that the plasma velocity associated with the instabilities being considered (i.e., the fluctuating part of the velocity profile) remains significantly smaller than the speed of sound. Note that it is only the inhomogeneous part of the velocity field that should be much smaller than the speed of sound. This might be a good assumption for some KH events, while other KH events might require including compressible effects. Notwithstanding this, in the present paper we adopt incompressibility as a simplifying assumption. Because of quasi-neutrality, the electron and the proton particle densities are equal, i.e., ${n}_{e}={n}_{i}={n}_{0}$. The (dimensionless) Alfvén speed is then ${v}_{{\rm{A}}}={B}_{0}/\sqrt{4\pi {m}_{i}{n}_{0}}{U}_{0}$, while η and ν are respectively the dimensionless magnetic diffusivity and kinematic viscosity. Note that for simplicity we assume isotropic expressions for both dissipative effects, even though in the presence of magnetic fields a tensor representation would be a more appropriate model (Braginskii 1965). These equations are complemented by the solenoidal conditions for both vector fields, i.e.,

$$\begin{eqnarray}&&{\boldsymbol{\nabla }}\cdot {\boldsymbol{B}}=0={\boldsymbol{\nabla }}\cdot {\boldsymbol{U}}.\end{eqnarray} \tag{ 3 }$$

Let us assume that the plasma is subjected to an externally applied shear flow given by

$$\begin{eqnarray}&&{{\boldsymbol{U}}}_{0}={U}_{0y}(x)\hat{{\boldsymbol{y}}},\end{eqnarray} \tag{ 4 }$$

so that the total velocity field is now ${{\boldsymbol{U}}}_{0}+{\boldsymbol{u}}$, where

$$\begin{eqnarray}&&{U}_{0y}(x)={U}_{0}\left[\mathrm{tanh}\left(\displaystyle \frac{x-\frac{\pi }{2}}{{\rm{\Delta }}}\right)-\mathrm{tanh}\left(\displaystyle \frac{x-\frac{3\pi }{2}}{{\rm{\Delta }}}\right)-1\right].\end{eqnarray} \tag{ 5 }$$

The velocity profile given in Equation (5) simulates the encounter of largely uniform flows of intensities $+{U}_{0}\hat{{\boldsymbol{y}}}$ and $-{U}_{0}\hat{{\boldsymbol{y}}}$ through a parallel interface of thickness $2{\rm{\Delta }}$. The numerical setup is sketched in Figure 1, where the jump provided by the hyperbolic tangent is duplicated to satisfy periodic boundary conditions throughout the numerical box. Also, we assume the presence of an external and uniform magnetic field ${{\boldsymbol{B}}}_{0}$ tangential to the interface and almost perpendicular to the shear flow (see Figure 1), so that the total magnetic field is ${{\boldsymbol{B}}}_{0}+{\boldsymbol{b}}$. The assumption of a hyperbolic tangent velocity profile is often adopted (Drazin 1958; Chandrasekhar 1961; Miura 1992) as a way to model shear flows with a finite thickness. The velocity profile given in Equation (5) is an exact equilibrium of Equations (1) and (2) obtained through the application of the stationary external force ${{\boldsymbol{F}}}_{0}=-\nu {{\rm{\nabla }}}^{2}{U}_{0y}(x)\hat{{\boldsymbol{y}}}$ (see also Gómez et al. 2014), and therefore it is numerically implemented simply by the application of the volume force ${{\boldsymbol{F}}}_{0}$.

Zoom In Zoom Out Reset image size

In the KH event that took place at one of the flanks of the 2010 November 3 CME, the fluid is observed to move along the contact discontinuity, albeit at very different speeds on either side. We choose to describe the development of the KH instability from a reference frame moving along the interfase at the average between these two speeds. In this reference frame, the flow will display a hyperbolic tangent type of profile, for which the parameter U₀ (see Equation (5)) will be equal to one-half of the relative velocity.

A shear flow such as the one given by Equation (5) is subjected to the well-known KH instability, which is of a purely hydrodynamic nature, i.e., it occurs even in the absence of any magnetic field. Within the framework of MHD, the stability of a tangential velocity discontinuity (i.e., in the limit of ${\rm{\Delta }}=0$) was first studied by Chandrasekhar (1961). For the case of an external magnetic field aligned with the shear flow, the mode is stabilized by the magnetic field, unless the velocity jump exceeds twice the Alfvén speed. For the case at hand, we assume the parallel component of the external magnetic field to be sufficiently weak (i.e., ${v}_{{\rm{A}}}^{\parallel }\lt 1$), since otherwise the instability pattern would not have been observed in AIA images. A stability analysis of a sheared MHD flow of finite thickness (i.e., ${\rm{\Delta }}\ne 0$) in a compressible plasma has also been performed (Miura & Pritchett 1982), confirming the result of the purely hydrodynamic case. Compressibility has a stabilizing effect in the sense that the growth rate is reduced as the velocity jump approaches the speed of sound, and even stalls the instability when the Mach number becomes unity (Miura & Pritchett 1982). From Table 2 of Foullon et al. (2013), we derive a shear flow amplitude ${U}_{0}=340\;\mathrm{km}\;{{\rm{s}}}^{-1}$, which remains below the speed of sound on both sides of the CME-pause. For the sake of simplicity, we neglect the effect of compressibility, which would bring an extra parameter to the problem: the Mach number. Yet another effect that might become relevant for the evolution of the KH instability is the presence of a density contrast between both sides of the shear flow (Prialnik et al. 1986; Gonzalez & Gratton 1994; Wyper & Pontin 2013). However, for the particular event under consideration it is not expected to play a role, since the mass density on both sides of the CME-pause remains virtually the same (Foullon et al. 2013).

If we approximate the hyperbolic tangent profile given in Equation (5) by piecewise linear functions, the instability growth rate ${\gamma }_{\mathrm{KH}}$ arising from the linearized version of Equations (1) and (2) is (for details, see Drazin & Reid (1981))

$$\begin{eqnarray}&&{\left(\displaystyle \frac{{\gamma }_{\mathrm{KH}}{\rm{\Delta }}}{{U}_{0}}\right)}^{2}=\displaystyle \frac{1}{4}({e}^{-4{k}_{y}{\rm{\Delta }}}-{(2{k}_{y}{\rm{\Delta }}-1)}^{2}),\end{eqnarray} \tag{ 6 }$$

which attains its maximum at ${\lambda }_{\mathrm{max}}\approx 15.7\ {\rm{\Delta }}$ and ${\gamma }_{\mathrm{max}}\approx 0.2{U}_{0}/{\rm{\Delta }}$, as shown in Figure 2. More importantly, Figure 2 also shows that the KH instability only arises for large-scale modes, i.e., such that ${k}_{y}{\rm{\Delta }}\leqslant 0.64$, corresponding to $\lambda \geqslant 9.82{\rm{\Delta }}$.

Zoom In Zoom Out Reset image size

We perform numerical integrations of Equations (1) and (2) subject to the external force ${{\boldsymbol{F}}}_{0}=-\nu {{\rm{\nabla }}}^{2}{U}_{0y}(x)\hat{{\boldsymbol{y}}}$ (where ${U}_{0y}(x)$ is given in Equation (5)) on the cubic box of linear size $2\pi$ sketched in Figure 1, assuming periodic boundary conditions in all three directions. The number of gridpoints is 256³ and the dimensionless Alfvén speed was set at ${v}_{{\rm{A}}}^{\parallel }=0.2$ in all our simulations, indicating that the component of the external magnetic field parallel to the flow (i.e., ${B}_{0y}$, see Figure 1) is such that its associated Alfvén velocity component remains smaller than the maximum velocity U₀ of the shear profile, and it is therefore unable to quench the instability. This is indeed the case of the 2010 November 3 KH event. Nykyri & Foullon (2013) performed a series of 2.5D MHD simulations seeking to match the time development of the KH pattern observed by AIA/SDO. Their numerical quest is consistent with slightly different magnetic field intensities on either side of the shear layer within the range of 8–11 G, forming small angles with the $\hat{z}$-direction (between 1° and 10°, see Figure 1), which leads to values of ${v}_{{\rm{A}}}^{\parallel }$ in the range of ${v}_{{\rm{A}}}^{\parallel }\quad \approx$ 0.04–0.31.

In our simulations, we use a pseudospectral method to perform the spatial derivatives and a second order Runge–Kutta scheme for the time integration (see a detailed description of the code in Gómez et al. 2005). For the viscosity and resistivity coefficients we chose $\nu =\eta ={2.10}^{-3}$, which are small enough to produce energy dissipation only at very small scales, comparable to the Nyquist wavenumber. In particular, dissipative effects are certainly negligible for all wavenumbers becoming unstable due to KH (see Equation (6) and the text right below it). The values of all the dimensionless parameters adopted for our simulations are summarized in Table 1. In all simulations, the pressure in Equation (1) is obtained self-consistently by taking the divergence of the equation, using the incompressibility condition, and solving at each time step the resulting Poisson equation for the pressure.

Table 1.  Values of Dimensionless Parameters for the Simulations: N is the Linear Size, U₀ is the Velocity on Each Side of the Shear Layer, Δ is the Thickness of the Shear Layer, ${v}_{{\rm{A}}}^{\parallel }$ is the Parallel Component of the Alfvén Speed, η is the Magnetic Diffusivity, ν is the Kinematic Viscosity, l_turb is the the Length Scale of the Turbulence, and f_turb is the Strength of the Turbulent Forcing

N	U0	Δ	${v}_{{\rm{A}}}^{\parallel }$	η	ν	lturb	fturb
256	1	0.1	0.2	${2.10}^{-3}$	${2.10}^{-3}$	0.05	0, 2, 3, 4, 5, 10, 15

Download table as:  ASCIITypeset image

The evolution of the $\hat{{\boldsymbol{z}}}$-component of vorticity is shown in Figure 3 at three different times, displaying the characteristic pattern of the KH instability. The observed frame corresponds to the right half of the numerical box displayed in Figure 1, which covers the shear layer centered at ${x}_{0}=3\pi /2$, and has been rotated for better viewing. The observed pattern shows the presence of the largest Fourier mode in our numerical box, characterized by ${k}_{y}=1$, whose growth rate according to Equation (6) is ${\gamma }_{\mathrm{KH}}({k}_{y}=1)\simeq 0.87$. At the same time, the presence of harmonics is also apparent, judging by the smaller scale patterns showing up as the instability progresses. In fact, from Equation (6) (see also Figure 2) we can anticipate which ones would be the growing Fourier modes.

Zoom In Zoom Out Reset image size

To estimate the instability growth rate, we use the component ${u}_{x}({x}_{0},y,z)$ evaluated at ${x}_{0}=\pi /2,3\pi /2$ (i.e., in the central part of the shear flows) as a proxy (see Figure 4). A Fourier analysis performed on ${u}_{x}({x}_{0},y,z)$ for any fixed value z confirms that the exponentially growing modes belong to the interval ${k}_{y}=1,...,6;$ which is consistent with the theoretical prediction shown in Figure 2 for ${\rm{\Delta }}=0.1$. Since the KH pattern is a two-dimensional flow taking place at z = constant planes, we take the maximum velocity of the profile ${u}_{x}({x}_{0},y,z)$ at any given value of z, and then average in the $\hat{{\boldsymbol{z}}}$-direction, i.e.,

$$\begin{eqnarray}&&{U}_{x,\mathrm{max}}={\displaystyle \int }_{0}^{2\pi }\displaystyle \frac{{dz}}{2\pi }\mathrm{max}[{u}_{x}({x}_{0},y,z),0\leqslant y\lt 2\pi ].\end{eqnarray} \tag{ 7 }$$

Zoom In Zoom Out Reset image size

In Figure 5 we show the maximum value of the ${u}_{x}({x}_{0},y,z)$ profile (averaged with respect to the $\hat{{\boldsymbol{z}}}$-direction) for both ${x}_{0}=\pi /2$ and ${x}_{0}=3\pi /2$, although as expected the two curves are undistinguishable. The straight gray line corresponds to the theoretically predicted growth rate ${\gamma }_{\mathrm{KH}}\simeq 0.87$ for the Fourier mode ${k}_{y}=1$ (using Equation (6)), which is the one observed in the time sequence shown in Figure 3. The fact that our empirical determination of the growth rate so strongly resembles ${\gamma }_{\mathrm{KH}}({k}_{y}=1)$ even though (as mentioned) the observed pattern is more complex than a single Fourier mode arises as the combined result of the z-averaging and our choice of the maximum of the velocity profile, as defined in Equation (7). Note that even though the simulations include dissipative effects and the theoretical prediction does not, the coincidence between both curves during the linear regime of the instability is nonetheless remarkable. Considering that the attenuation effect of viscosity can be estimated by $\gamma \simeq {\gamma }_{\mathrm{KH}}-\nu {k}_{y}^{2}$, we can easily verify that the dissipative correction is absolutely negligible for the evolution of the KH instability, as expected.

Zoom In Zoom Out Reset image size

To generate a turbulent background in our simulations, we apply a stationary force to all modes within a thin spherical shell of radius ${k}_{\mathrm{turb}}=1/{l}_{\mathrm{turb}}$ in Fourier space, consisting of a superposition of harmonic modes with random phases. The nonlinear interactions between these Fourier modes that are being externally driven with a force of intensity f_turb will develop a stationary turbulent regime with its associated energy cascade involving all wavenumbers $k\geqslant {k}_{\mathrm{turb}}$. To make sure that it is a small-scale turbulence, we chose l_turb to be much smaller than the wavelength observed for the KH pattern, and even somewhat smaller than the thickness Δ of the shear layer (i.e., ${l}_{\mathrm{turb}}\lt {\rm{\Delta }}$).

The pattern of vorticity obtained when only the turbulent forcing is applied (i.e., a simulation with no KH driving) is shown in Figure 6. The observed pattern corresponds to a turbulent regime which is statistically stationary, homogeneous, and isotropic. Even though all spatial scales from l_turb down to the smallest scales available in the simulation participate in the dynamics and in the ensuing energy cascade, only those vortices of sizes comparable to l_turb can be identified, which is to be expected for a power-law power spectrum with a negative index such as Kolmogorov's. Therefore, these concentrations of vorticity can safely be associated with the energy-containing eddies of the turbulence. As mentioned in Section 1, the expected effect of this small-scale turbulence on a larger-scale flow is an effective or enhanced diffusivity. In the case at hand, its effect on the instability growth rate is expected to be

$$\begin{eqnarray}&&\gamma (k)={\gamma }_{\mathrm{KH}}(k)-{\nu }_{\mathrm{turb}}{k}^{2},\end{eqnarray} \tag{ 8 }$$

where ${\gamma }_{\mathrm{KH}}(k)$ is given in Equation (6) and ${\nu }_{\mathrm{turb}}$ is the aforementioned effective or turbulent viscosity. The effect of increasing turbulent viscosity on the instability growth rate is illustrated in Figure 7, showing that it is not only the growth rate that is reduced but also the range of unstable wavenumbers.

Zoom In Zoom Out Reset image size

Zoom In Zoom Out Reset image size

We performed simulations applying both the large-scale force ${{\boldsymbol{F}}}_{0}$ to drive the KH instability and the small-scale force of intensity f_turb to drive the turbulent regime. In Figure 8 we show the resulting distribution of the vorticity component ${\omega }_{z}(x,y)$, which can be compared with the one shown in Figure 3 for the KH instability on a laminar background and the one shown in Figure 6 for the purely turbulent run, with no KH pattern. We can qualitatively see that the role of turbulence is in fact an attenuation in the growth of the instability.

Zoom In Zoom Out Reset image size

One of the observable consequences of this turbulent regime is the nonthermal broadening of coronal spectral lines caused by the turbulent motion of the fluid elements emitting these (optically thin) spectral lines. Once this turbulence reaches a Kolmogorov-stationary regime, the rms value of the turbulent velocity u_turb is

$$\begin{eqnarray}&&{E}_{\mathrm{turb}}=\displaystyle \frac{{u}_{\mathrm{turb}}^{2}}{2}={\displaystyle \int }_{1/{l}_{\mathrm{turb}}}{dk}\ {\epsilon }^{2/3}\ {k}^{-5/3}\propto {(\epsilon {l}_{\mathrm{turb}})}^{2/3},\end{eqnarray} \tag{ 9 }$$

where E_turb is the (dimensionless) kinetic energy density of the turbulence and is its energy dissipation rate. Note that neither or E_turb are known a priori, since they arise as a result of the stationary regime attained by the turbulence. However, using heuristic arguments we can find how these quantities scale with the input parameters of this turbulence, namely, l_turb and f_turb. The fluid is energized by the work done per unit time by the external force of intensity f_turb at scale l_turb; energy then cascades down to smaller scales and it is dissipated by viscosity at the rate at dissipative scales. In a stationary regime, the power delivered by the external force should match the energy dissipation rate, i.e.,

$$\begin{eqnarray}&&\epsilon \propto {f}_{\mathrm{turb}}\ {u}_{\mathrm{turb}}.\end{eqnarray} \tag{ 10 }$$

Equations (9) and (10) can be combined to obtain both u_turb and in terms of f_turb and l_turb,

$$\begin{eqnarray}&&\epsilon \propto {({f}_{\mathrm{turb}}^{3}{l}_{\mathrm{turb}})}^{1/2},\end{eqnarray} \tag{ 11 }$$

$$\begin{eqnarray}&&{u}_{\mathrm{turb}}\propto {({f}_{\mathrm{turb}}{l}_{\mathrm{turb}})}^{1/2}.\end{eqnarray} \tag{ 12 }$$

On dimensional arguments, the turbulent viscosity introduced in Equation (8) has to be proportional to the turbulent velocity u_turb times the typical scale l_turb, i.e., ${\nu }_{\mathrm{turb}}\propto {u}_{\mathrm{turb}}\ {l}_{\mathrm{turb}}$, which, considering Equation (12)

$$\begin{eqnarray}&&{\nu }_{\mathrm{turb}}=C{({f}_{\mathrm{turb}}{l}_{\mathrm{turb}}^{3})}^{1/2}.\end{eqnarray} \tag{ 13 }$$

To quantify the role of turbulence in the evolution of the KH instability, we performed a sequence of simulations for which the only parameter being changed is the turbulent forcing f_turb. As the parameter f_turb is gradually increased, the corresponding turbulent velocity u_turb (observationally perceived as nonthermal broadening of spectral lines) is also increased, which in turn raises the turbulent viscosity ${\nu }_{\mathrm{turb}}$. As a result, the instability growth rate (see Equation (8)) is expected to be reduced. To estimate the instability growth rate from our simulations, we follow the same procedure described in Section 4, which amounts to following the temporal evolution of the profile ${u}_{x}(y)$ for the gridpoints centered at the shear layer. Note, however, that now the velocity at each grid point can be split into a part corresponding to the large-scale KH evolution plus another part corresponding to the turbulence.

Because of the geometrical setup of our simulations, the large-scale part of the flow at each z = constant plane is an exact replica of one another (KH is a two-dimensional flow) while the turbulent part is not, since it is a fully three-dimensional flow. The averaging procedure in the $\hat{{\boldsymbol{z}}}$-direction described in Equation (7) gets rid of the turbulent part of the flow, since the mean velocity of this turbulence is exactly zero. We can also compute the rms deviation of the velocity when averaging in the $\hat{{\boldsymbol{z}}}$-direction, which should exactly correspond to u_turb, since the KH part of the flow is identical for all z = constant planes. Therefore, this statistical strategy allows us to obtain the main features of both the large-scale (i.e., the KH instability) and small-scale (the turbulence) components of this complex flow.

Figure 9 shows the main result of the present study, which is the value of ${U}_{x,\mathrm{max}}$ (defined in Equation (7)) as a function of time in a lin–log plot, for runs corresponding to different turbulent intensities. The thick black lines correspond to ${U}_{x,\mathrm{max}}(t)$ for each simulation, the thin black lines indicate one standard deviation with respect to the average (i.e., ${U}_{x,\mathrm{max}}\pm {u}_{\mathrm{turb}}$), and the straight gray lines are the theoretical predictions for each case, as emerges from Equation (8). Note that the theoretical slopes (i.e., the gray lines in Figure 9) are not the best fits to each of the simulations, but the result arising from Equation (8), which contains only one free parameter for the whole set of simulations, namely the constant C. This constant is the only dimensionless parameter that remains undetermined by the dimensional analysis described above. We find that the value of C that best fits all our simulations is $C\approx 18.8$.

Zoom In Zoom Out Reset image size

In the previous section, we presented results from numerical simulations showing the role of a background turbulence in reducing the growth rate of an ongoing KH instability. These numerical results are intended to simulate the KH instability being developed at the interface between some CMEs and the ambient corona, which have been recently reported in the literature. There is also mounting observational evidence about the turbulent nature of the solar corona, mostly related with spatially unresolved motions leading to measurable nonthermal broadenings in coronal spectral lines.

To numerically model this turbulent background, we made a number of simplifying assumptions. For instance, we assume the turbulent regime to be spatially homogeneous and isotropic and also stationary. We maintain this turbulent state throughout the whole simulation by applying a stationary stirring force of intensity f_turb at a well-defined lengthscale l_turb. We deliberately chose this lengthscale to be much smaller than the wavelength of the KH-unstable mode, since the AIA images reporting the KH pattern do not show any observable evidence of a turbulent background. Also, the rotation period of the energy-containing vortices is of the order of ${\tau }_{\mathrm{turb}}\simeq {l}_{\mathrm{turb}}/{u}_{\mathrm{turb}}$, which remains shorter than the instability growth time for all the cases considered. The properties of this turbulent regime are therefore determined by only two input parameters: l_turb, which is kept fixed throughout the whole study, and f_turb, which is varied to give rise to cases with different turbulent velocities (u_turb) and effective viscosities (${\nu }_{\mathrm{turb}}$).

We can use Equations (12) and (13) to express the effective viscosity ${\nu }_{\mathrm{turb}}$ in terms of two measurable quantities such as u_turb and l_turb. A crude estimate of the dimensionless constant in Equation (12) leads to ${u}_{\mathrm{turb}}\approx 0.22\ {({f}_{\mathrm{turb}}{l}_{\mathrm{turb}})}^{1/2}$ and therefore

$$\begin{eqnarray}&&{\nu }_{\mathrm{turb}}\approx 85.4\ {u}_{\mathrm{turb}}\ {l}_{\mathrm{turb}}.\end{eqnarray} \tag{ 14 }$$

If we refer for instance to the KH event occurred on 2010 November 3 and reported by Foullon et al. (2011), they estimate a velocity jump at the interface of ${U}_{0}=340\;\mathrm{km}\;{{\rm{s}}}^{-1}$ and a wavelength for the KH pattern of $\lambda =2\pi {L}_{0}=18.5\;\mathrm{Mm}$ (corresponding to a length unit of ${L}_{0}=3\;\mathrm{Mm}$ and ${k}_{y}=1$ in our simulations). For ${k}_{y}=1$, the dispersion relation reduces to $\gamma \approx 0.87-{\nu }_{\mathrm{turb}}$, as shown in Equation (8). The instability growth rate estimated by Foullon et al. (2013) for this event is $\gamma \approx 0.033\;{{\rm{s}}}^{-1}$, which in our dimensionless units becomes $\gamma {L}_{0}/{U}_{0}=0.29=0.87-{\nu }_{\mathrm{turb}}$. From this expression we can estimate the value of ${\nu }_{\mathrm{turb}}$ required to explain the growth rate observed for this particular KH event. More interestingly, using Equation (14) we can obtain a level of turbulent velocity of ${u}_{\mathrm{turb}}\approx 47\;\mathrm{km}\;{{\rm{s}}}^{-1}$ (for the value of l_turb used in our simulations), which is well within the range reported by Doschek et al. (2014) from Hinode observations. It is important to recall that other effects besides turbulence might contribute to reduce the instability growth rate. Depending on the parameter values of the particular KH event being considered, the compressibility of the plasma or the strength of the magnetic field component along the shear flow might play a role.

Another consequence that we can derive from the present analysis is that, given the fact that the turbulence did not completely suppress the KH instability, we can in principle use Equations (8)–(14) to estimate an upper bound for l_turb for any observed value of u_turb. For the turbulent attenuation to be negligible (i.e., ${\nu }_{\mathrm{turb}}\ll 0.87$) and assuming a turbulent velocity of $60\;\mathrm{km}\;{{\rm{s}}}^{-1}$ (see Doschek et al. 2014), we obtain for l_turb an upper bound of 170 km. In general,

$$\begin{eqnarray}&&{l}_{\mathrm{turb}}\ll 170\;\mathrm{km}{\left(\displaystyle \frac{{u}_{\mathrm{turb}}}{60\mathrm{km}{{\rm{s}}}^{-1}}\right)}^{-1}.\end{eqnarray} \tag{ 15 }$$

In summary, in order for the invoked turbulent state to produce nonthermal broadening of spectral lines of the order of u_turb and at the same time not to affect the observed KH event in any appreciable manner, the typical size l_turb of its energy-containing eddies should satisfy Equation (15).

The study presented in this paper was motivated by two relatively recent observational findings on the nature of the solar corona. One of them is the apparent development of the KH instability as some CMEs expand in the ambient corona, as shown by AIA/SDO images (Foullon et al. 2011; Ofman & Thompson 2011; Foullon et al. 2013). The second one is that the coronal plasma seems to be in a turbulent state, as evidenced by the nonthermal broadening of coronal spectral lines measured from EIS/Hinode data (Doschek et al. 2008, 2014; Brooks & Warren 2011; Tian et al. 2012).

Our main goal has been to study the feasibility for these two apparently dissimilar features to coexist. Namely, the large-scale laminar pattern observed for the KH instability, and the small-scale spatially unresolved turbulent motions leading to the observed nonthermal broadenings. We therefore performed three-dimensional simulations of the MHD equations, to study the evolution of the KH instability in the presence of a turbulent ambient background for different intensities of this turbulence.

Theoretically, the effect of a small-scale turbulence on a large-scale flow would be to produce an enhanced diffusivity which can be modeled by an effective or turbulent viscosity. The impact of this small-scale turbulence on an ongoing large-scale instability such as KH would then be a reduction of its growth rate, as emerges from Equation (8). The degree of this reduction is controlled by the turbulent viscosity ${\nu }_{\mathrm{turb}}$ which we obtained from a dimensional analysis to be ${\nu }_{\mathrm{turb}}=C{({f}_{\mathrm{turb}}{l}_{\mathrm{turb}}^{3})}^{1/2}$ (see Equation (13)), leaving only the dimensionless constant C undetermined.

The comparison between the instability growth rates obtained from our simulations with the ones arising from Equation (8) esentially confirms this theoretical scenario, while providing an empirical determination for the dimensionless constant C, which amounts to $C\approx 18.8$. Perhaps more importantly, since ${\nu }_{\mathrm{turb}}\propto {u}_{\mathrm{turb}}\ {l}_{\mathrm{turb}}$ and given the fact that the instability has not been completely quenched by the turbulence (otherwise it would not have been observed), observational determinations of u_turb from nonthermal broadenings pose an upper limit to the correlation length of the turbulence l_turb. For observational values of ${u}_{\mathrm{turb}}\quad \approx$ 20–60 km.s⁻¹, the correlation length of turbulence is expected to be smaller than about ${l}_{\mathrm{turb}}\quad \approx$ 170–510 km, which is consistent with not having been spatially resolved by current coronal imaging spectrometers such as EIS on board Hinode.

D.G. and E.E.D. acknowledge financial support from grant SP02H1701R from Lockheed-Martin to SAO. E.E.D. was also supported by NASA contract NNM07AB07C. D.G. also acknowledges support from PICT grant 0454/2011 from ANPCyT to IAFE and P.D.M. acknowledges support from PICTs 2011-1529 and 2011-1626 from ANPCyT to IFIBA (Argentina).