© ESO 2022

1. Introduction

Recent studies have revealed that the turbulence in magnetized plasmas is greatly affected by the Alfvén wave effects. The well-documented example is the solar-wind turbulence, whose nature is essentially Alfvénic, and turbulent fluctuations can be approximately described as Alfvén waves (AWs) (Belcher & Davis 1971; Bruno & Carbone 2013). The standard magnetohydrodynamic (MHD) description of Alfvénic turbulence in astrophysical and laboratory plasmas is based on the interaction of oppositely propagating incompressible wave packets (Iroshnikov 1963; Kraichnan 1965).

Following significant previous work on the weak turbulence in incompressible MHD (Sridhar & Goldreich 1994; Montgomery & Matthaeus 1995; Ng & Bhattacharjee 1996; Galtier et al. 2000), more recent work (Howes & Nielson 2013) has described the mechanism of turbulent energy transfer via AW collisions in more detail. The authors showed analytically that two colliding counter-propagating AWs with wave vectors and first produce a specific intermediate wave with , and then its interaction with the initial waves produces the tertiary waves with the following wave vectors: and . Here, , , and are the unit Cartesian vectors, such that is parallel to the background magnetic field B₀. These analytical results have been confirmed by both gyrokinetic simulations in the MHD limit (Nielson et al. 2013) and experimentally in the laboratory (Drake et al. 2013, 2014, 2016). Since the energy is transferred to AWs with higher perpendicular wave numbers, this process represents an elementary step of the direct turbulent cascade in which energy is transferred from larger to smaller scales.

Goldreich & Sridhar (1995) introduced the critical balance conjecture and developed their famous model of strong anisotropic MHD turbulence. The critical balance assumes that the linear (wave-crossing) and nonlinear (eddy turnover) times are equal at each scale. Whereas the critical balance remains a physically reliable hypothesis not strictly derived from basic principles, it allows for a phenomenological prediction of turbulence properties, in particular the energy spectrum and anisotropy of turbulent fluctuations. The Goldreich & Sridhar model gave rise to many important insights in the turbulence nature and resulted in many theoretical, numerical, and experimental studies (see, e.g., Verniero & Howes 2018; Verniero et al. 2018; Mallet et al. 2015, and references therein). It is worth noting that the critical balance conjecture is essentially a statement implying persistence of linear wave physics in the strongly turbulent plasma.

Despite extended investigations of the critically balanced turbulence, many current problems remain, such as the nonzero cross-helicity effects in the presence of shear plasma flows (Gogoberidze & Voitenko 2016), or nonlocal effects in AW collisions (Beresnyak & Lazarian 2008). Also, the plasma compressibility can introduce surprising effects to the behavior of MHD waves (Magyar et al. 2019).

Numerical simulations of turbulence are usually done either via numerical codes for reduced MHD or using analytical (Beresnyak 2014, 2015; Mallet et al. 2015; Perez et al. 2020), pseudo-spectral (Chandran & Perez 2019) and gyrokinetic (Verniero et al. 2018) frameworks. Pezzi et al. (2017a,b) performed simulations using compressible MHD, Hall MHD, and hybrid Vlasov-Maxwell codes; the 2.5D geometry used in these works does not allow to take into account nonlinear terms ∼​(v^±⋅∇)v^∓ and ∼​(b^±⋅∇)b^∓ (v^± and b^± are velocity and magnetic fluctuations in ± waves) for AWs with .

Using a compressible MHD model in 3D, we numerically studied the commonly accepted presumption that the AW turbulence is generated by the collisions between counter-propagating AWs, particularly the wave-number dependence of the amplitudes of induced waves. Our simulations reveal that the AW collisions can occur in two regimes, the first one corresponding to the case of strong turbulence, which follows theoretical explanation, and the second one corresponding to larger scales, which obviously is governed by a different mechanism.

2. Physical and numerical setup

The simulations were performed in 3D using the numerical code MPI-AMRVAC (Porth et al. 2014). The code applies the Eulerian approach for solving the compressible resistive MHD equations:

(1)

(2)

(3)

(4)

where e, ρ, v, and B are the total energy density, mass density, velocity, and magnetic field, p = (γ − 1)(e − ρv²/2 − B²/2) is the thermal pressure, p_tot = p + B²/2 is the total pressure, J = ∇ × B is the electric current density, η is the electrical resistivity, and γ is the ratio of specific heats. The magnetic field is measured in units for which the magnetic permeability is 1. Since in this study we are not interested in dissipative processes, we take γ = 5/3, and η = 0. We used the following three normalization constants: the length L_N = 1 Mm, the magnetic field B_N = 20 G, and the density ρ_N = 1.67 × 10⁻¹⁵ g cm⁻³. This determined normalization for other physical quantities: electron concentration n_N = 10⁹ cm⁻³, speed km s⁻¹, and time t_N = L_N/v_N = 0.7246 s.

The simulations were performed in 3D in Cartesian geometry with a rectangular numerical box. The background magnetic field B₀ = 20 G is directed along the z-axis. Equilibrium plasma parameters were set typical for the solar coronal base: n_e = 10⁹ cm⁻³ (ρ₀ = 1.67 × 10⁻¹⁵ g cm⁻³) and temperature T = 1 MK, which determines the plasma beta parameter β = 0.017. The Alfvén speed in equilibrium plasma is km s⁻¹ or v_A = 1 in normalized units, and the sound speed is .

In order to induce counter-propagating Alfvén waves, we set the components of magnetic field and velocity at the z boundaries of the simulation volume. The forward wave propagating in the +z direction along B₀ is initiated at z = 0 by the following forcing:

(5)

(6)

(7)

and the backward wave propagating in −z direction is initiated at z = z_max:

(8)

(9)

(10)

where ω = k_zv_A = 2π(λ_z)⁻¹v_A is the angular frequency, the parallel wavelength λ_z = 10 Mm (always the same constant in all setups), the initial amplitudes of magnetic field b, and velocity were set either 3.33% or 10% of B₀ and v_A, respectively, and represents the ponderomotive component of the speed (its order is 10⁻³). Boundary conditions at other boundaries were periodic. Introduction of A_p≠ 0 insures a smooth solution of the MHD equations at the boundaries; its influence is studied in Sect. 3.1. The physical configuration is shown in Fig. 1.

Fig. 1. General picture of the numerical setup for AW collision. Top: physical setup; early phase. The green rectangle denotes the numerical box with equilibrium plasma, the red and blue areas represent velocity perturbations (positive and negative) of the z− wave (far boundary) and z+ wave (near boundary). Polarization planes and nonuniformity directions are annotated. Bottom: longitudinal sketches; late phase after AW collision. Main setups with one period in z− and three periods in the z+ wave (left), and complimentary setups with one period in each wave (right). The gray areas denote cross-sections taken for further analysis.

The forcing described above was applied during one period for the forward wave and during three periods for the backward wave, which we call the main setups hereafter (see Table 1 for setup parameters). Beside the main setups, we ran several complementary simulations without backward wave, or with different amplitudes of counter-propagating waves, or setups with a single period in both waves.

As suggested by the nonlinear term (z^±⋅∇)z^∓ in the Elsässer form of MHD equations, in order to allow for effective interactions, the counter-propagating AWs should have different polarizations. In our setups, the forward wave was polarized along x and its wave vector ; the backward wave was polarized along the y-axis and (see Fig. 1).

The numerical box has physical z length L_z of either 50 Mm (main setups) or 30 Mm (complementary setups). The sizes along x and y are set as equal to the perpendicular wavelength λ_⊥ (hence, they change from setup to setup). The numerical box for the main setups has either 256  ×  256  ×  512 pixels (high-resolution) or 128  ×  128  ×  256 pixels (low-resolution). We verified that the decrease of numerical resolution does affect the results: the waves start to decay during their propagation and the wave profiles become distorted. However, this effect is small even for the case of low-resolution setups. In complimentary setups, the numerical box always has 256  ×  256  ×  384 pixels; thus, its spatial resolution coincides with that of the high-resolution setups. We compared various numerical schemes and parameters of MPI-AMRVAC and chose the best settings (powel scheme for the ∇ ⋅ B corrector, high-resolution numerical box, etc.). We also paid special attention to distinguish the physical phenomena from numerical artifacts.

3. Results

3.1. Nonlinear effects in a single AW

First, we verified the effect of nonlinear self-interaction within a single Alfvén wave. In Fig. 2, we show the longitudinal (along z) profiles of v_x, v_y, and v_z of the forward Alfvén wave, initiated via boundary conditions described in Sect. 2. The instant shown corresponds to the developed phase of the wave but before the collision with the backward wave. For visualization, the quantities are normalized by the following constants: the mother wave v_x by the initial amplitude u = 0.10, the horizontal component v_y and the ponderomotive component v_z by A_p = 2.53 ⋅ 10⁻³.

Fig. 2. Velocities vx (left), vy (middle), and vz (right) in a single Alfvén wave. The color table for each quantity have different amplitudes to reflect the range of velocities. Initially, only vx and Bx (not shown) and vz are driven in the Alfvén wave; the vy component of the velocity is generated self-consistently due to nonlinear self-interaction within the Alfvén wave.

The amplitude and spatial structure of the ponderomotive component v_z perfectly reproduce theoretical predictions: its wave numbers are two times larger than in the mother wave, and its amplitude varies from 0 to 2 (McLaughlin et al. 2011; Zheng et al. 2016). We also observed a self-consistent generation of v_y that appears only in oblique waves with λ_⊥ ≠ 0 (v_y = 0 at λ_⊥ = 0). Our preliminary simulations (2D setups were sufficient there) show the following trend in the variation of v_y with varying cross-field wavelength: the amplitude of v_y grows proportionally to 1/λ_⊥ at the larger scales λ_⊥ > λ_z, this growth slows down at λ_⊥ ∼ λ_z, and eventually v_y becomes a constant independent on λ_⊥ at smaller scales λ_⊥ ≪ λ_z. The spatial extension of v_y in both parallel and perpendicular directions is two times shorter than that of the mother wave. The amplitude of v_y is always smaller than that of v_z. Similar perturbations of perpendicular velocity were also observed in torsional waves (Shestov et al. 2017).

The observed perturbations of v_y and v_z propagate along the magnetic field with the Alfvén speed v_A and are natural companions of AWs not caused by the numerical effects or boundary conditions for v_z. The perturbations always develop in AWs regardless of the how the waves are initiated – that is, by boundary or initial conditions, with or without boundary perturbations given by Eqs. (7) and (10). In other words, the observed propagating wave is the eigenmode of the compressible nonlinear MHD. We thus observe typical characteristics of AWs before they collide.

3.2. AW collision

To study effects of the AW collisions, we allowed the two counter-propagating waves to fully propagate through each other, and we analyzed perpendicular profiles of v_x of the forward-propagating z⁻ wave in its leading maximum – x − y plane with z = 42.5 at instant t = 45 (see Fig. 1, bottom left panel; main setups with u = 0.1 are used). In Fig. 3, the panels show v_x of three different setups with λ_⊥ = 0.5, 0.8, and 3.0 Mm. The perturbations of the wave profiles depend on the perpendicular scale: they are significant for the smallest λ_⊥ = 0.5, moderate for λ_⊥ = 0.8, and weak for the largest λ_⊥ = 3.0. Similar perturbations are also observed in the z⁺ wave. However, in setups with only one wave present, such perturbations do not appear, and hence their development can be attributed to AW collision. The appearance of such small-scale perturbations propagating with Alfvén velocity can be treated as a generation of new AWs at smaller perpendicular scales .

Fig. 3. Comparison of perpendicular profiles of vx in the leading maximum of the z− wave in setups with λ⊥ = 0.5 (left), 0.8 (middle), and 3.0 Mm (right). The axes are in Mm and are different in every panel.

3.3. Dependence on perpendicular scales

In order to distinguish the nonlinearly generated waves from the mother wave, we further analyzed the wave profiles in the perpendicular cross-section of z⁻ wave. We extracted the perturbed velocity Δv_x by subtracting the initial harmonic profile of v_x, , where the amplitude A is adjusted to cancel the perturbation in the wave maximum. The results are shown in Fig. 4 for the setups with λ_⊥ = 0.8 (top panels) and λ_⊥ = 3.0 (bottom panels). The left panels show v_x, middle Δv_x, and right v_y. The induced velocities Δv_x and v_y have amplitudes ∼(0.02 ÷ 0.03)u and are nonuniform in both x and y directions.

Fig. 4. Perpendicular profiles of velocities in setups with λ⊥ = 0.8 (top) and λ⊥ = 3 (bottom). Left: measured vx; middle: difference Δvx between the measured vx and a harmonic function; right: measured vy. Perturbations of Δvx and vy are produced as a result of AW collisions. In each panel the color table matches the maximum amplitude of the measured quantity.

To evaluate numerical effects, we made the similar analysis for z⁻ wave in the absence of z⁺ waves. Here, the perturbations Δv_x are observed as well; but they have an amplitude at least a factor of 10 smaller and are uniform along x. It means that numerical effects produce significantly weaker perturbations with different spatial profiles. On the contrary, after collisions with counter-propagating z⁺ waves, the perturbations co-propagating with z⁻ waves have both v_x and v_y components, larger amplitudes, and profiles that are nonuniform along both y and x, which cannot be ascribed to numerical effects. Furthermore, the perturbations of v_y generated by the AW collisions can not be attributed solely to the single AW self-interaction where perturbations of the v_y are zero at the original wave maximum.

The spatial patterns of the induced velocities fall in two distinct groups: all spatial patterns at are similar to that shown in the top panels of Fig. 4, and all patterns at are similar to that shown in the bottom panels (the transition scale for u = 0.1 used in this figure). The perturbations in the former group have a current-sheet structuring, similar to that reported by Verniero et al. (2018) for the strong turbulence regime. The perturbations in the second group have a symmetric structure. The same two groups of spatial structures are also observed in the setups with different amplitudes u, but with different transition scales, such that is larger for smaller u (for example, for u = 0.033).

The dependence of the amplitudes of induced waves on the perpendicular scales is shown in Fig. 5. The diamonds correspond to Δv_x and asterisks correspond to v_y. For u = 0.1 the symbols are blue and green, for u = 0.033 they are orange and red. Gray and pink regions indicate the wave number ranges where the wave collisions should generate the strong (critically balanced) turbulence with k_⊥/k_∥ ∼ v_A/u for u = 0.1 and u = 0.033, respectively. In both these cases, the amplitude behavior is similar. At the largest λ_⊥, the amplitudes of the induced waves are much smaller than the amplitudes of the original waves, and the resulting AW turbulence should be weak. As λ_⊥ decreases, the induced amplitudes first increase slowly and reach a maximum. This maximum is still much smaller than the initial AW amplitude and is reached at λ_⊥ = λ_⊥max that is still much larger than the perpendicular scale given by the critically balance condition, λ_⊥max≫ λ_⊥* = λ_∥u/v_A (λ_⊥max, λ_⊥* and other characteristic perpendicular scales are shown in Fig. 5). When λ_⊥ decreases further beyond λ_⊥max, the induced amplitudes decrease and reach a minimum at λ_⊥= λ_⊥min, which is still larger than λ_⊥*. After this minimum, a strong increase of induced perturbations occurs in the region where λ_⊥ becomes several times shorter than λ_⊥*. Amplitudes of generated perturbations become there comparable to the amplitudes of initial waves and such collisions can generate strong turbulence.

Fig. 5. Dependence of amplitudes of induced waves on k⊥/k∥ of the original wave for the u = 0.1 (top) and u = 0.033 (bottom) setups. Diamonds and asterisks correspond to the Δvx and vy, respectively. Filled areas correspond to the regions of classic critical balance in incompressible MHD calculated for a particular u. In the top panel, the values of λ⊥max, λ⊥min, , and are shown.

While the observed strengthening of the nonlinear interaction with decreasing λ_⊥ is expected when taking into account the fact that the responsible nonlinear term is ∼(z⋅∇)z∼ , the depression observed at λ_⊥≳ λ_⊥* and its influence on the transition from weak to strong turbulence need further investigation. At present, we can only state that this depression should result in a shift of the weak-strong turbulence transition to the perpendicular scales significantly shorter than that prescribed by the standard critical balance condition.

3.4. Influence of several collisions

Since the initiated z⁻ and z⁺ waves contain one and three periods, respectively, the z⁻ wave can interact with three periods of the counter-propagating wave, whereas each period of the z⁺ wave can only interact with one period of z⁻. We thus expect different amplitudes of the induced perturbation propagating in z⁻ and z⁺ directions. To verify this, we measured the perturbations accompanying the z⁺ wave using the same technique as for z⁻ wave (remember that in z⁺ wave the roles of v_x and v_y are exchanged). A comparison of the corresponding perturbations in the z⁻ and z⁺ waves is given in Fig. 6. The top panel shows the amplitudes of perturbations accompanying z⁻ (blue and green symbols) and z⁺ (orange and red symbols) and the bottom panel shows the ratio z⁻/z⁺ of the perturbations with the corresponding (orthogonal) polarizations. In both panels, the diamonds denote perturbations with the same polarization as in the original waves (Δv_x in z⁻, Δv_y in z⁺), and the asterisks denote the complimentary polarization.

Fig. 6. Dependence of amplitudes of induced waves on the wave-number ratio for different numbers of wave collisions (z− wave collides 3 times and z+ wave collides once). Top: amplitudes of the perturbations. Bottom: ratio of the amplitudes z−/z+, diamonds correspond to original polarization (vx in z−, and vy in z+), and asterisks correspond to perpendicular polarization.

The behavior of z⁺ perturbations as a function of λ_⊥ is qualitatively similar to that of z⁻ perturbations. At the smallest λ_⊥, the ratio of the −/+ perturbations is about 2, and it then approaches 3 with the scale increase, before increasing significantly at λ_⊥/λ_∥ ∼ 1, and it finally drops again to 2 at large perpendicular scales λ_⊥/λ_∥ > 1. In the region of (super-) strong turbulence, the observed ratio of z⁻/z⁺ < 3 denotes the inapplicability of the perturbation theory: after the first interaction, the wave profiles are already significantly distorted, and the following collisions do not add much.

3.5. Nonzero cross-helicity case

In this section, we analyze the effects of nonzero cross-helicity (imbalance) when the counter-propagating initial waves have different amplitudes. This situation is common in the fast solar wind (Tu et al. 1990; Lucek & Balogh 1998) and also occurs in numerical simulations in local subdomains of the simulation box (Perez & Boldyrev 2009).

We ran dedicated setups with initial amplitudes of u⁻ = 0.1 in the z⁻ wave and u⁺ = 0.033 in the z⁺ wave; both waves have one period. We compare measured perturbations with our main setups in Fig. 7. The black symbols denote z⁻ perturbation, and orange and red symbols denote z⁺ perturbations in imbalanced setups, and blue and green symbols denote main setups (u = 0.1, 1 period in z⁻ wave and 3 periods in z⁺ wave).

Fig. 7. Amplitude of induced waves for the case of nonzero cross-helicity: u− = 0.1, u+ = 0.033, each wave has a single period. Black symbols denote perturbations in the z− wave, and orange and red symbols denote perturbations in the z+ wave. Blue and green symbols represent the main setups (u = 0.1, 1 and 3 periods, respectively).

The perturbations observed in imbalanced cases are smaller than in the main setups. At the same time, the perturbations (expressed in initial amplitudes u) in the z⁻ wave are ∼3 times smaller than in the z⁺ wave.

3.6. Perpendicular Fourier spectra

In order to understand the spectral transport generated by the AW collisions, we analyzed the spatial Fourier spectra of the induced waves. The spectra of the v_x and v_y velocities at the leading maximum of z⁻ are given in Fig. 8 for λ_⊥ = 0.8 in the top row, λ_⊥ = 2.0 in the middle row, and λ_⊥ = 3.0 in the bottom row. In the left panels, the spectra of v_x in a single z⁻ wave are shown, in the middle and right panels the spectra of v_x and v_y after the AW collision are shown. In each panel, the (x, y)-coordinates represent corresponding Fourier wave numbers, and the color shows the intensity of a given spectral component. The quasi-logarithmic color scale is normalized to the intensity of an ideal harmonic function . This function would only have two peaks with spectral coordinates (0, ±1) that correspond to the brightest components in the left and middle panels. In what follows, we drop the ± sign, since the spectra are symmetrical in y.

Fig. 8. Perpendicular Fourier spectra of vx and vy velocities measured in the leading maximum of the z− wave. Top row: λ⊥ = 0.8; middle row: λ⊥ = 2; bottom row: λ⊥ = 3. Left: vx of the setup with the z− wave only; middle and right denote vx and vy of the setups after AWs collision. In each panel, the (x, y)-coordinates represent corresponding Fourier wave numbers, and the color shows intensity of a given spectral component. The quasi-logarithmic color table is normalized to the amplitude of a harmonic wave.

The higher wave number spectral components (0, |y| > 1) accompanying the single z⁻ wave without collisions (left panels) are due to numerical effects. We note the low level of these components and their uniform distribution. On the contrary, the real spectral components with higher wave numbers are generated by the AW collisions (middle and right panels).

The strongest induced components at λ_⊥ = 0.8 have spectral coordinates (1, 2) corresponding to the perpendicular wave vector . The generation of waves with such wave vectors supports the mechanism proposed by Howes & Nielson (2013) (see their Fig. 2 explaining the appearance of such “tertiary” waves). This mechanism is summarized in the Introduction.

The same Fourier components (1, 2) of v_x are also seen in the middle row of Fig. 8 in the case of intermediate scale λ_⊥ = 2; in addition, the spectral components of v_x with coordinates (1, 1) corresponding to are also significant.

The spatial spectra of v_x and v_y at the largest scale λ_⊥ = 3 are qualitatively different: the strongest induced components have coordinates (1, 1), while the others are negligible. These spectral components might be formed by a different mechanism than in the λ_⊥ = 0.8 case. In general, the spectral dynamics observed in our simulations, that is generation of higher wave number spectral components, supports scenarios with direct turbulent cascades generated by AW collisions.

3.7. Field-aligned structure of the induced Alfvén waves

Longitudinal behavior of the Fourier components of the induced Alfvén waves is studied using the following approach: we Fourier-analyzed perpendicular cross-sections at multiple z coordinates, covering the distance of slightly more then one full wavelength λ_∥ along z (see bottom right sketch in Fig. 1). In Fig. 9, we show longitudinal behaviour of the spectral components of v_x with coordinates (0, 1), (1, 1), (1, 2), (2, 1), and (2, 2) with different colors. The mother wave with spectral coordinates (0, 1) is shown in black. The top panel shows the setups with λ_⊥ = 0.5, and the bottom panel shows the setup with λ_⊥ = 2.0. The intensity of the spectral components is multiplied by a factor of 10 in the top panel, and by a factor of 100 in the bottom panel.

Fig. 9. Longitudinal dependence of amplitudes of spatial Fourier components of vx. Top: setup with λ⊥ = 0.5. Bottom: setup with λ⊥ = 2 Mm. Different colors correspond to particular spectral coordinates.

We observe a drastically different behavior of the spectral components in different setups. While we do not see any regularity in the larger scale setup, in the setup with λ_⊥ = 0.5 the growth of (1, 1) and (1, 2) components is highly correlated and their parallel scales are somehow shorter than in initial AWs. In addition, the energy of the induced waves tend to concentrate near the center of the mother wave.

4. Discussion and application

Results of our simulations revealed several new properties of AW collisions in compressional plasmas, which can affect Alfvénic turbulence and anisotropic energy deposition in plasma species. The most striking new property is the modified relation between the parallel and perpendicular scales in the strong turbulence regime where the energy is efficiently transferred to the smaller scale during one collision.

The turbulence strength is usually characterized by the nonlinearity parameter , where is the nonlinear mixing time, is the linear (correlation) crossing time of colliding AWs, and v_k is the velocity amplitude of the colliding AWs. δv_k denotes the velocity amplitude of generated waves. When the classic critical balance condition of incompressible MHD is satisfied,

(11)

the nonlinear mixing becomes as fast as the linear crossing, and the turbulence is believed to be strong, δv_k/v_k ∼ 1 (Goldreich & Sridhar 1995).

However, as follows from our simulations (see, e.g., Fig. 5 showing δv_k/v_k as a function of k_⊥/k_z for two fixed amplitudes, u ≡ v_k/v_A = 0.1 and 0.033, and ), the spectral transport in compressible MHD is strongly (about one order of magnitude) suppressed at k_⊥/k_z, satisfying Eq. (11). Namely, δv_k/v_k ≪ 1 at k_⊥/k_z = 10 for u = 0.1 and at k_⊥/k_z = 30 for u = 0.033. At k_⊥/k_z increasing further, the spectral transport eventually becomes fast and the turbulence strong, δv_k/v_k ∼ 1, which happens at k_⊥/k_z obeying the modified critical balance condition,

(12)

where α < 1 is the factor reducing efficiency of the nonlinear mixing (in other words, the effective nonlinear time increases by the factor 1/α). Consequently, the turbulence becomes strong at perpendicular wave numbers that are larger than in the classic, critically balanced case.

The origin and nature of α need further clarification. Since α ≠ 1 arises when the plasma compressibility is taken into account, it should depend on the relative content of thermal energy, that is on the plasma β. For parameters adopted in our simulations, α ≈ 0.3 in the critically-balanced state where the scale ratio k_⊥/k_z obeys . Such departure from the classic critical balance affects dynamics of the strong AW turbulence (see below). In the general case of arbitrary scales, the functional dependence α = α(β,v_k,k_⊥/k_z) is complex; in particular, δv_k/v_k becomes a decreasing function of k_⊥/k_z in some intervals (in Fig. 5, it happens at k_⊥/k_z ≲ 10), which should greatly affect the weak AW turbulence. We do not exclude that α may also depend on other parameters of plasma and waves.

Let us consider the Alfvénic turbulence driven by the fluctuating velocity v_k0 at the wave number ratio k_⊥0/k_z0, obeying the critical balance condition , in which case the turbulence is already strong at the driving scales. The spectral energy flux in the inertial range is

(13)

where is the AW collision time and ρ ≈ n₀m_i is the mass density. We note that the spectral flux ϵ_s from Eq. (13) is reduced as compared to the incompressible strong turbulence driven at the same perpendicular scale, but it remains the same for the turbulence driven at the same parallel scale.

We assume that there is a weak dependence α = α₀(v_k/v_k0)^δ, where 0 < δ < 3/4. Such dependence is suggested by the following semi-empirical considerations. As the observed spectra are power laws, the scaling of α with v_k should be a power law as well. Furthermore, the index δ of the power-law dependence should be small and positive to reproduce the observed in simulation mismatches between the classic and real critical balances (which is larger for larger wave amplitude). Moreover, such positive values of δ appear to be compatible with the observed spectral indexes of turbulence in the quasi-stationary solar wind, which are slightly larger than −5/3 (up to −3/2).

The kinetic energy spectrum is then flatter than the Kolmogorov one:

(14)

and its spectral index varies between −5/3 and −3/2, as is typically observed in the solar wind turbulence. In the case of α constant along the critical balance path, δ = 0, the spectrum reduces to the classic Kolmogorov The parallel wave number spectrum is, as usual, .

If the turbulence is weak at injection, , the cascade time increases from the strong turbulence value to the weak turbulence value . The resulting weakly turbulent energy flux ϵ_w decreases as compared to the strongly turbulent energy flux (13), :

(15)

The weakly turbulent spectrum is problematic to calculate because of a complex dependence of δv_k/v_k upon k_⊥ and v_k (see Fig. 5), which is unknown and difficult to estimate. At present, we can only note that the strength of the compressional weak turbulence is much smaller (about one order, as is demonstrated by Fig. 5) than the incompressional one, , which drastically decreases the not very turbulent energy flux.

Although the large-scale MHD AWs do not dissipate directly, the turbulent cascade transfers their energy to small scales where dissipative effects come into play in heating plasma. MHD Alfvénic turbulence has been employed as the mechanism for plasma heating in the solar corona and solar wind, both from the theoretical/modeling perspective (Van Ballegooijen et al. 2011; Verdini et al. 2012) and based on experimental observations of quiescent (Morton et al. 2016; De Moortel et al. 2014; Xie et al. 2017) and flaring loops (Doschek et al. 2014; Kontar et al. 2017). Here, we discuss how the new properties of AW collisions observed in our simulations can affect models of quasi-steady turbulent plasma heating in coronal loops.

Recently, Xie et al. (2017) analyzed as many as 50 loops in active regions using observations of the Extreme-ultraviolet Imaging Spectrometer (EIS) (Culhane et al. 2007) on board the Hinode satellite. They observed nonthermal widths of spectral lines and found corresponding nonthermal velocities in the v_nt = 30 ÷ 40 km s⁻¹ range, magnetic field in the loop apexes up to 30 G, loop widths of L_⊥ ∼ 2 ÷ 4 Mm, and loop lengths of L_z ∼ 100 Mm. Brooks & Warren (2016) also used spectroscopic data from EIS and evaluated nonthermal velocities in loops in 15 active regions. The typical values were somewhat smaller, with typical values of v_nt ∼ 20 km s⁻¹; the authors did not, however, provide any other parameters. Furthermore, Gupta et al. (2019) analyzed nonthermal widths of spectral lines in high coronal loops (with heights up to 1.4 R_⊙) measured by EIS and found the nonthermal velocities to be in the 20 ÷ 30 km s⁻¹ range. The above values can be used to evaluate the turbulent heating of coronal loops.

We assume that there are AW sources at the loop footpoints. These sources can be due to magnetic reconnection and/or photospheric motion (we do not specify their origin in more detail here). The perpendicular AW wavelengths λ_⊥0 are limited by the cross-B₀ scale l_⊥ of density filaments comprising the loops, λ_⊥0 ≲ l_⊥ (wave number k_⊥0 ≳ 2π/l_⊥). We note that l_⊥ can be significantly smaller than the visible loop width r. On the contrary, the coronal plasma is quite homogeneous along B₀ and the possible parallel wavelength λ_z0 are restricted by the loop length L_z, λ_z0 = 2π/k_z0 ≤ L_z (wave number k_z0 ≳ 2π/L).

For the wavelengths within the limits mentioned above, a large spectral flux, and hence a strong plasma heating, can be established by the strong turbulence driven at the critically balanced anisotropy αk_⊥0/k_z0= αλ_z0/λ_⊥0= v_A/v_k0. The corresponding energy flux injected in the unit volume is . Assuming that the turbulent velocity at injection v_k0 is observed as the nonthermal velocity, v_k0 ≈ v_nt, and taking v_nt = 30 km s⁻¹, magnetic field B = 30 G, and density n_e = 2 ×  10⁹ cm⁻³ from Xie et al. (2017), we obtain the energy flux 3 × 10⁻⁴(l_⊥/L_⊥)⁻¹ erg cm⁻³ s⁻¹. For sufficiently small l_⊥ ≲ 0.1L_⊥, the energy flux ϵ_cor≳ 3 × 10⁻³ erg cm⁻³ s⁻¹ is enough to heat typical coronal loops. The corresponding parallel wavelengths at injection are λ_z0 ∼ (αv_nt/v_A)⁻¹λ_⊥0 ≲ 0.5L_z. Therefore, the turbulent cascade and related plasma heating can be effective if the perpendicular length scales of the loop substructures are about ten times smaller than the loop width, which implies that the loops should be structured more than was required by previous turbulent heating models.

5. Conclusions

In the framework of compressional MHD, we studied the spectral transport produced by the collisions between counter-propagating Alfvén waves numerically. The initial two waves are linearly polarized in two orthogonal planes, and their cross-field profiles vary normally from their polarization planes. Polarization and spectral characteristics of the perturbations generated after single and multiple collisions between such AWs were analyzed in detail. The main properties of the resulting spectral transfer are as follows:

• The perturbations generated by AW collisions have smaller scales than the original waves, which supports turbulence scenarios based on the direct turbulent cascade generated by AW collisions;

• We observed two regimes of the AW interaction: the first one is typical for the case of strong turbulence, and the second one is governed by a different mechanism;

• The spectral transfer generated by the AW collisions is strongly suppressed at the scales satisfying the classic critical balance condition (11) of incompressional MHD, which makes the turbulence weak at these scales;

• The strong turbulence is reestablished at significantly smaller perpendicular scales satisfying the modified critical balance condition (12).

We used these properties to reevaluate the turbulent heating of the solar coronal loops. The main conclusion is that the turbulent cascade can heat the loop plasma, provided the loop is structured and the characteristic widths of the loop substructures are more than ten times smaller than the loop width.

References

All Tables

All Figures

Fig. 1. General picture of the numerical setup for AW collision. Top: physical setup; early phase. The green rectangle denotes the numerical box with equilibrium plasma, the red and blue areas represent velocity perturbations (positive and negative) of the z− wave (far boundary) and z+ wave (near boundary). Polarization planes and nonuniformity directions are annotated. Bottom: longitudinal sketches; late phase after AW collision. Main setups with one period in z− and three periods in the z+ wave (left), and complimentary setups with one period in each wave (right). The gray areas denote cross-sections taken for further analysis.

In the text

	Fig. 2. Velocities vx (left), vy (middle), and vz (right) in a single Alfvén wave. The color table for each quantity have different amplitudes to reflect the range of velocities. Initially, only vx and Bx (not shown) and vz are driven in the Alfvén wave; the vy component of the velocity is generated self-consistently due to nonlinear self-interaction within the Alfvén wave.
In the text

	Fig. 3. Comparison of perpendicular profiles of vx in the leading maximum of the z− wave in setups with λ⊥ = 0.5 (left), 0.8 (middle), and 3.0 Mm (right). The axes are in Mm and are different in every panel.
In the text

	Fig. 4. Perpendicular profiles of velocities in setups with λ⊥ = 0.8 (top) and λ⊥ = 3 (bottom). Left: measured vx; middle: difference Δvx between the measured vx and a harmonic function; right: measured vy. Perturbations of Δvx and vy are produced as a result of AW collisions. In each panel the color table matches the maximum amplitude of the measured quantity.
In the text

	Fig. 5. Dependence of amplitudes of induced waves on k⊥/k∥ of the original wave for the u = 0.1 (top) and u = 0.033 (bottom) setups. Diamonds and asterisks correspond to the Δvx and vy, respectively. Filled areas correspond to the regions of classic critical balance in incompressible MHD calculated for a particular u. In the top panel, the values of λ⊥max, λ⊥min, , and are shown.
In the text

	Fig. 6. Dependence of amplitudes of induced waves on the wave-number ratio for different numbers of wave collisions (z− wave collides 3 times and z+ wave collides once). Top: amplitudes of the perturbations. Bottom: ratio of the amplitudes z−/z+, diamonds correspond to original polarization (vx in z−, and vy in z+), and asterisks correspond to perpendicular polarization.
In the text

	Fig. 7. Amplitude of induced waves for the case of nonzero cross-helicity: u− = 0.1, u+ = 0.033, each wave has a single period. Black symbols denote perturbations in the z− wave, and orange and red symbols denote perturbations in the z+ wave. Blue and green symbols represent the main setups (u = 0.1, 1 and 3 periods, respectively).
In the text

Fig. 8. Perpendicular Fourier spectra of vx and vy velocities measured in the leading maximum of the z− wave. Top row: λ⊥ = 0.8; middle row: λ⊥ = 2; bottom row: λ⊥ = 3. Left: vx of the setup with the z− wave only; middle and right denote vx and vy of the setups after AWs collision. In each panel, the (x, y)-coordinates represent corresponding Fourier wave numbers, and the color shows intensity of a given spectral component. The quasi-logarithmic color table is normalized to the amplitude of a harmonic wave.

In the text

	Fig. 9. Longitudinal dependence of amplitudes of spatial Fourier components of vx. Top: setup with λ⊥ = 0.5. Bottom: setup with λ⊥ = 2 Mm. Different colors correspond to particular spectral coordinates.
In the text