4.1. Global

In § 2.2 we saw that the collisional strength $\tilde {\alpha }$ sets the degree of coupling between the two fluids: in the small-$\tilde {\alpha }$ limit we expect the two fluids to move independently, whereas, in the large-$\tilde {\alpha }$ limit, they should be almost identical, following the MHD equations (2.8a) and (2.8b). In figure 1 we see three snapshots of the vorticity for each species, $\omega _s \equiv \boldsymbol {\nabla } \times \boldsymbol {v}_s$, representing, from left to right, $\tilde {\alpha } \ll 1$, $\tilde {\alpha } \sim 1$ and $\tilde {\alpha } \gg 1$, all from the K4 runs. The top row shows the neutral vorticity for each run, whereas the bottom row shows the ion vorticity for those runs. By visual inspection of figure 1(a) and (d), with $\tilde {\alpha } = 0.003$, one can indeed see that for small $\tilde {\alpha }$ the two species behave as they would in a completely uncoupled regime. The neutral vorticity shows clear signs of the 2-D HD inverse cascade of kinetic energy given by the two large-scale vortices, whereas the ion vorticity shows many filamented structures, typical of 2-D MHD turbulence. In this regime, we can approach the question of whether or not a one-fluid description is adequate for the proper determination of the dynamics by going into the centre-of-mass frame, typically done in other plasma settings (e.g. when deriving MHD itself). We define $\boldsymbol {V} \equiv \chi \boldsymbol {v}_i + (1-\chi ) \boldsymbol {v}_n$ and $\boldsymbol {D} = \boldsymbol {v}_i - \boldsymbol {v}_n$. The question now becomes: Is it possible only to account for the dynamics of $\boldsymbol {V}$ without knowing or integrating the dynamics of $\boldsymbol {D}$? Although it is not shown here, the simulations reveal that in the $\tilde {\alpha } \ll 1$ regime this is generally not possible – that is, the dynamics of $\boldsymbol {V}$ is partially determined by $\boldsymbol {D}$ and vice versa. However, as $\chi \rightarrow 0$ or $\chi \rightarrow 1$, we find that a one-fluid description ($\boldsymbol {V}$ only) is sufficient, as one might expect. This was shown by comparing the relative magnitudes of $\boldsymbol {V}$ and $\boldsymbol {D}$ and noting when $|\boldsymbol {D}|\ll |\boldsymbol {V}|$.

Figure 1.

Snapshots at steady state from K4 runs of the vorticity for each species, $\omega _s \equiv \boldsymbol {\nabla } \times \boldsymbol {v}_s$, representing, from left to right, $\tilde {\alpha } \ll 1$, $\tilde {\alpha } \sim 1$ and $\tilde {\alpha } \gg 1$. Panels (a)–(c) show the neutral vorticity for the three values of $\tilde {\alpha }$ whereas panels (d)–(f) show the ion vorticity for those runs.

As we increase $\tilde {\alpha }$ so that it is of order one, we see the neutral species begin to lose the large-scale vortices, which are dissipated away by collisional heating. In figure 1(b) and (e), at $\tilde {\alpha } = 3.15$, we no longer see obvious 2-D HD behaviour from the neutral species, but it also does not appear to be of a similar nature to the ion vorticity. Later scale-by-scale analysis will reveal that this regime is where the highest collisional heating is found and that most of the energy injected into the neutrals will be transferred to the ions or simply dissipated away. In figure 1(c) and (f), we see the case of $\tilde {\alpha } = 314.98$, and note immediately that the two fluids look identical from visual inspection. The two fluids are coupled and so the neutral species is behaving like an MHD fluid, as was predicted.

Although our predictions seem to agree qualitatively given the snapshots in figure 1, we now aim to confirm our results from § 2.2 in a quantitative way. In figure 2 we see two panels comparing global (time- and volume-averaged) quantities, where each data point is a single simulation. All runs performed in this study are included. Each panel aims to test the predictions made for each extreme of $\tilde {\alpha }$.

Figure 2.

(a) Non-dimensional total energy versus the collision strength $\tilde {\alpha }$, rescaled by the energy injection rate of the neutrals over the ionization fraction. The red dashed line shows the prediction from (2.6), valid in the $\tilde {\alpha } \ll 1$ limit. (b) Non-dimensional collisional heating, rescaled by the square of the averaged Lorentz force, versus collision strength $\tilde {\alpha }$. The red dashed line shows the prediction from (2.8c), valid in the $\tilde {\alpha } \gg 1$ limit.

Figure 2(a) compares total energy with the collision strength rescaled by the energy injection rate of the neutral species over the ionization fraction, based on the right-hand side of (2.6), and with tildes representing suitable non-dimensionalization by a combination of $u_f$ and $k_f$. The red dashed line shows the prediction from (2.6), which is valid for small values of $\tilde {\alpha }$, and represents the balance between collisional heating and energy injection rate into the neutral species. The collapse of the data – for various values of ionization fraction, $Re$ and forcing wavenumber – on a single line, whose slope agrees with the predicted relationship over various orders of $\tilde {\alpha }$, confirms our claims.

Figure 2(b) focuses on the highly coupled regime, where the prediction for collisional heating was based on an expansion over $\tilde {\alpha }^{-1}$, revealing that $CH$ was given by the square of the difference in the forces acting on each species – given by (2.8c). This prediction is quite general but simplifies significantly for our simulations due to the fact that densities are uniform, $\nu _i = \nu _n$, and $\boldsymbol {F}_i = \boldsymbol {F}_n$. After non-dimensionalizing using $u_f$ and $k_f$, we are left with $\widetilde {CH} = (1-\chi )\chi ^{-1} \tilde {\alpha }^{-1} \widetilde {LF}$, where $\widetilde {LF} = \langle | \tilde {\boldsymbol {J}} \times \tilde {\boldsymbol {B}} |^2 \rangle$ and $\widetilde {(\cdot)}$ denotes the dimensionless version. Moving everything except $\tilde {\alpha }$ to the left-hand side, we get that, in the high collisional limit, the rescaled collisional heating ($y$-axis of figure 2b) should be proportional to $\tilde {\alpha }^{-1}$, denoted by a red dashed line. We see that indeed the data collapse onto a single line, once again for various ionization fractions, $Re$ and forcing wavenumbers. The slope of the rescaled collisional heating seems to approach the theoretical one for large values of $\tilde {\alpha }$, confirming our predictions for the highly collisional limit.

4.2. Spatial

Our spatially averaged (‘global’) predictions for the extreme limits of $\tilde {\alpha }$, taken from § 2.2, have been confirmed. However, turbulence is an out-of-equilibrium, multi-scale process whose scale-by-scale analysis can reveal further interesting phenomena that are difficult to identify or predict otherwise. This is the purpose of the current subsection.

In our scale-by-scale analysis we will look at two general quantities – spectra, which tell how a quadratic quantity is distributed over scales, and fluxes, which tell us how that quadratic quantity is flowing through scales (Alexakis & Biferale Reference Alexakis and Biferale2018). The one-dimensional spectrum $KE_s(k)$ of the kinetic energy of a species $s$ is

(4.1)\begin{equation} KE_s(k) = \frac{1}{2} \sum_{| \boldsymbol{k}|=k} | \widehat{\boldsymbol{v}_s}|^2( \boldsymbol{k}), \end{equation}

where $\widehat {(\,\cdot \,)}$ denotes the Fourier transform of $\boldsymbol {v}_s$. We will also be looking at the dimensionless collisional heating spectra,

(4.2)\begin{equation} \widetilde{CH}(k) = \frac{(1-\chi) \chi \alpha}{u^3_f k_f} \sum_{| \boldsymbol{k}|=k} | \widehat{\boldsymbol{v}_i - \boldsymbol{v}_n}|^2( \boldsymbol{k}), \end{equation}

and the dimensionless spectra of the square of the Lorentz force

(4.3)\begin{equation} \widetilde{LF}(k) = \frac{1}{k_f^2 u_f^4} \sum_{| \boldsymbol{k}|=k} | \widehat{\boldsymbol{J} \times \boldsymbol{B}}|^2( \boldsymbol{k}). \end{equation}

We have seen that $\tilde {\alpha }$ measures the relative strength between collisions and the nonlinear term at a typical scale $L$ with velocity $U$ (or, in the case of our simulations, $k^{-1}_f$ and $u_f$). Therefore, it measures the degree of coupling between the two fluids at $L$. However, due to the multi-time-scale nature of turbulence, different length scales couple at different values of $\tilde {\alpha }$. The global results from the last subsection showed us how to predict the average coupling of the fluid, given some information about the collision strength at the forcing scale. However, a scale-by-scale analysis allows us to see exactly how the two fluids gradually couple together among scales and to better understand non-extreme cases for $\tilde {\alpha }$. In figure 3 we see the steady-state average kinetic energy spectra of both species, for three different values of $\tilde {\alpha }$, as in figure 1. These spectra have ionization fraction $\chi = 0.5$ and were taken from the K8 set of runs, forced at intermediate wavenumbers so as to be able to resolve approximate inertial ranges in both large and smaller scales.

Figure 3.

Kinetic energy spectra for the neutral species (orange, solid line) and ion species (green, dashed line) at (a) $\tilde {\alpha } = 0.002$, (b) $\tilde {\alpha } = 2.520$ and (c) $\tilde {\alpha } = 29.923$, all with ionization fraction $\chi =0.5$. As $\tilde {\alpha }$ increases, the two species become more and more coupled, starting from the largest scales which couple first.

In figure 3(a) we have the uncoupled limit, evident by the two spectra and their distinct shapes. The black dot-dashed lines show a $-5/3$ slope, the prediction for the large scales of $KE_n$ and the small scales of $KE_i$. The slope of the neutral kinetic energy is a bit steeper at large scales due to the presence of a condensate. Despite this, as well as a poorly resolved inertial range for the ion kinetic energy, we see reasonable agreement between expected and observed behaviour for each individual species. Figure 3(b) shows $\tilde {\alpha } = 2.52$, a moderately coupled case. In this figure we see that, indeed, the energies at scales down to the forcing scale seem to be lying on top of each other, implying a coupling at those scales. The collisions are strong enough at these scales to completely remove the inverse cascade of kinetic energy in the neutral species, either by dissipation or by transfer of energy to the ions (to be discussed further when we look at the fluxes). However, at scales smaller than the forcing there is a clear distinction between the two energies, implying a lack of coupling. It is this remaining ‘slippage’ between the ions and the neutrals in the small scales, along with the fact that $\tilde {\alpha } \sim \textit{O}(1)$, that maintains a large collisional heating, despite the small value of $|\boldsymbol {v}_i-\boldsymbol {v}_n|$ at large scales. In fact, collisional heating is largest at these values of $\tilde {\alpha }$, as seen in figure 2(b).

Finally, figure 3(c) shows a higher collisional case with $\tilde {\alpha } = 29.923$. In this figure we observe that scales smaller than the forcing scale are now beginning to couple; however, $\tilde {\alpha }$ is small enough so that the smallest scales are still not coupled. The two fluids look almost identical, like the same MHD fluid. At this point a natural question arises: At what scales do the fluids decouple and how does that depend on $\tilde {\alpha }$? We have attempted to numerically investigate the dependence of the decoupling wavenumber, which we are calling $k_{coll}$, on $\tilde {\alpha }$, but were not able to reach a definitive conclusion. Possible issues include: a small inertial range (limited by the numerical capabilities); a failure of time-scale assumptions common in homogeneous, isotropic turbulence which only are valid in a true inertial range (here, dissipation due to collisions might render that irrelevant); and more. We do, however, think that the results are worth showing, and so we have included this work in appendix B.

A quantity intimately related to the degree of coupling between the two species is the collisional heating, which is proportional to $|\boldsymbol {v}_i - \boldsymbol {v}_n|^2$. A scale-by-scale decomposition of $CH$ can tell us where kinetic energy dissipation (heating) due to collisions is most active, which might be of interest in astrophysical applications. In figure 4 we plot the non-dimensional spectra of $CH$, the black dashed curves, for the usual three cases of $\tilde {\alpha }$. These are taken from steady-state averages in the K32 set of runs, which were forced at small scales.

Figure 4.

The spectra of dimensionless collisional heating $\widetilde {CH}$ (bold, black, dashed lines) and of the dimensionless right-hand side of (2.8c), proportional to the spectrum of the Lorenz force squared, $\widetilde {LF}$ (bold, light blue, solid lines). These two curves are expected to be equal when the two species are highly coupled. Shown are three limits of $\tilde {\alpha }$, namely (a) weak coupling, (b) moderate coupling and (c) strong coupling, allowing one to see the gradual coupling of the species and the change in form of the $CH$ spectrum.

Focusing only on the black dashed curves for now, we look first at figure 4(a), the low coupling case. Since the two species are not coupled, we can approximately say that, at scales larger than the forcing, the kinetic energy of the neutrals dominates and hence $CH \propto |\boldsymbol {v}_i-\boldsymbol {v}_n|^2 \approx |\boldsymbol {v}_n|^2 \propto KE_n$. This is confirmed by the $-5/3$ slope shown by the thin dot-dashed black line. Hence, the collisional heating is acting like a friction term for the neutral species, causing significant dissipation at the largest scales where the condensate lies. Although not shown here, apart from directly dissipating energy to heat, collisions also act to transfer energy from neutrals to ions at large scales (particularly for low ionization fractions). This possibility is evident by looking at (2.2a) – when $|\boldsymbol {v}_n|\gg |\boldsymbol {v}_i|$, the sign-indefinite term will dominate $|\boldsymbol {v}_i|^2$ for the ion species kinetic energy equation, thus providing the possibility to gain kinetic energy via collisions with the highly energetic neutrals at those scales. The species first begin to couple at the largest scales, as we saw when looking at the kinetic energy spectra. This coupling reduces $CH$, causing the length scale of maximal collisional heating to become smaller and smaller as we increase $\tilde {\alpha }$.

This is exemplified in figure 4(b), where we already see a significant change in the shape of the collisional heating spectrum, its peak around $k = 15$. This happens until $\tilde {\alpha } > 1$, after which the maximal heating is at the forcing scale and the shape of the $CH$ spectrum remains practically unchanged, with only its magnitude being reduced (proportional to $\tilde {\alpha }^{-1}$). The shape of the collisional heating spectrum at scales that are already coupled is set by the spectrum of $\langle |\boldsymbol {J} \times \boldsymbol {B}|^2 \rangle$, as predicted by (2.8c). The light blue, solid lines in figure 4 represent the spectra of the right-hand side of this equation, and indeed the two curves lie on top of each other at the scales that are coupled. Thus, understanding the spectrum of the square of the Lorentz force is key to understanding the scale-by-scale collisional heating. The small dark blue dashed line depicting a slope of $2.4$ (found by fitting that region of the spectrum) can be seen in figure 4(a) and (c), showing that the spectrum of the Lorentz force has not changed and is not affected by the collisions. One would still like to know what sets the spectrum of the Lorentz force (and hence of $CH$). A very clear power law with a positive exponent is seen in the figure, so one is tempted to speculate about its origins, particularly using equilibrium spectra arguments, as is done in 3-D homogeneous and isotropic turbulence. However, such equilibrium spectra are not necessarily universal if the forcing in the spectral shell around $k_f$ is dense, as is ours (Alexakis & Brachet Reference Alexakis and Brachet2019). A prediction of the spectral slope of the Lorentz force, although interesting and relevant for our applications, is beyond the scope of this paper.

The spectra have revealed the length scales at which coupling between the species occurs, as well as allowed us to observe at what scales collisional heating dominates and confirm our predictions about the spectral shape of $CH$. While these spectra have informed us of the distribution of various quantities among length scales, also of interest in turbulence is how certain conserved quantities, such as energy, move across scales. This is measured by the (spectral) flux. In this study, we will only look at two components of energy flux, what we call the kinetic component of energy flux for a species $s$, $\Pi _{U_s}$, and the magnetic component of energy flux for the ions, denoted by $\Pi _{B}$. The former is defined to be

(4.4)\begin{equation} \Pi_{U_s}(k) = \langle \boldsymbol{v}^{{<}k}_s \boldsymbol{\cdot} (\boldsymbol{v}_s \boldsymbol{\cdot} \boldsymbol{\nabla} \boldsymbol{v}_s) \rangle, \end{equation}

where $\boldsymbol {v}^{<k}_s$ stands for a filtering of the velocity $\boldsymbol {v}_s$ in Fourier space so that only the wavenumbers with modulus smaller than $k$ are kept. The flux $\Pi _{U_s}$ expresses the rate at which kinetic energy of a species $s$ is flowing out of scales larger than $\ell = 2 {\rm \pi}/ k$ due to nonlinear interactions only in velocity. Therefore, if energy is going from large to small scales, the energy flux will be positive, and vice versa. The magnetic component of energy flux for the ion species is defined to be

(4.5)\begin{equation} \Pi_{B}(k) ={-} \langle \boldsymbol{v}_i^{{<}k} \boldsymbol{\cdot} (\boldsymbol{B} \boldsymbol{\cdot} \boldsymbol{\nabla} \boldsymbol{B} ) \rangle + \langle \boldsymbol{B}^{{<}k} \boldsymbol{\cdot} ( \boldsymbol{v}_i \boldsymbol{\cdot} \boldsymbol{\nabla} \boldsymbol{B} - \boldsymbol{B} \boldsymbol{\cdot} \boldsymbol{\nabla} \boldsymbol{v}_i ) \rangle. \end{equation}

As a reminder, in 2-D HD turbulence, energy flows to larger scales whereas 2-D MHD turbulence has been shown to cascade total energy forward to smaller scales. Thus, in an uncoupled system with $\tilde {\alpha } = 0$, we expect the neutral kinetic energy flux $\Pi _{U_n}$ to be negative at scales larger than the forcing and zero at scales smaller than the forcing. On the other hand, $\Pi _{U_i} + \Pi _{B}$ ought to be zero at large scales and positive at scales smaller than the forcing. Although $KE_i$ and $E_B$ are not individually conserved, the nonlinear interactions that cause the respective fluxes $\Pi _{U_i}$ and $\Pi _{B}$ each conserve energy and can thus also be analysed.

In figure 5 we look at all four fluxes, $\Pi _{U_n}$, $\Pi _{U_i}$, $\Pi _{B}$ and $\Pi _{U_i}+\Pi _{B}$, normalized by the total energy injection rate $I_{tot} = I_n + I_i + I_B$. We look at three cases of $\tilde {\alpha }$, the first two taken from the K8 set of runs, and the large $\tilde {\alpha }$ limit from the K8R run (all with $\chi = 0.5$). These fluxes are averaged over the steady state. The $\tilde {\alpha } \ll 1$ case, figure 5(a), confirms the uncoupled predictions. The neutral kinetic energy flux (orange solid line) is negative and roughly equal to the neutral injection rate: $I_n/ I_{tot} = 4/9$ (recall that $I_n = I_i = 4 I_B$). The total ion energy flux (black solid line) is positive and seen to be close to $(I_i+I_B)/I_{tot} = 5/9$. Since $\tilde {\alpha } \neq 0$, we do see a hint of some large-scale energy dissipation in the neutrals, given by the slight negative slope of the flux at large wavenumbers.

Figure 5.

Kinetic component of energy flux for the neutral species, $\Pi _{U_n}$ (orange, solid line), and ion species, $\Pi _{U_i}$ (green, dashed line), as well as the magnetic component of energy flux for the ions, $\Pi _B$ (blue, dot-dashed line). Total ion energy flux is denoted by the solid, black line, and grand total energy flux $(\Pi _{U_n}+\Pi _{U_i}+\Pi _{B})$ by the thin, solid, red line, seen only in panel (c). The three panels depict the three limiting cases for $\tilde {\alpha }$: (a) uncoupled, where we see the HD and MHD predictions for fluxes; (b) moderate coupling, which is highly dissipative; and (c) high coupling, in which the neutrals are fully coupled to the ions and we recover the prediction from § 2.2 that the system behaves as a single MHD fluid with $\boldsymbol {F} = \boldsymbol {F}_i + \boldsymbol {F}_n$.

Turning now to the case where $\tilde {\alpha } \approx 1$, figure 5(b) shows that the neutral kinetic energy flux is practically zero for all wavenumbers, meaning that collisions have almost completely stopped any nonlinear transfer of kinetic energy among the neutrals. At the forcing scale and larger, the two species are coupled and behave like MHD, as we have seen. Since MHD has no inverse cascade of kinetic energy, the neutrals have no inverse cascade at these scales. At scales smaller than the forcing, the species are not yet coupled, and since the neutral species does not have a forward cascade of energy, the flux remains zero. However, there is still constant injection rate of kinetic energy into the neutrals of magnitude $I_n$, and this energy must go somewhere. That energy is either dissipated away by collisions, or transferred to the ion species. Seeing as, at best, approximately half of the total energy injection rate is being transferred spectrally by the ions, it seems that in this case approximately half of the energy injected at the forcing scale is immediately dissipated by collisions, the rest being transferred away by the ion species. For the case of uncoupled 2-D MHD turbulence, one would expect this flux of energy by the ion species and magnetic field to be constant in the inertial range down to the dissipation scales. However, in figure 5(b) we see a gradual drop-off of energy flux. This non-constant flux is a sign of strong dissipation of energy in the would-be inertial range due to collisional heating, which is largest in the runs where $\tilde {\alpha } \sim \textit{O}(1)$, as we have seen.

Moving on to figure 5(c), we notice that the ion and neutral species are fully coupled, even at the small scales, as seen by the fact that $\Pi _{U_n} = \Pi _{U_i}$. Now that the two species are coupled throughout most of the scales, the collisional heating is very weak and does not dissipate practically any energy that is being injected at the forcing scale. Figure 5(c) indicate that all of the energy being injected by the forcing – including the neutral energy injection – is going to small scales. This is shown by the total energy flux from all fields, denoted by the thin, red, solid line, which we see is positive and very close to 1. This is consistent with our $\tilde {\alpha } \rightarrow \infty$ limit studied in § 2.2, where we showed that the system behaves as a single MHD fluid with $\boldsymbol {F} = \boldsymbol {F}_i + \boldsymbol {F}_n$, as seen in (2.8a). Sticking to $\boldsymbol {v}_i$ and $\boldsymbol {v}_n$, what seems to be happening is that the two species exchange energy between each other very efficiently due to collisions, but, since ions also exchange energy with the magnetic field, part of the energy being injected into the neutral species ends up as magnetic energy, which is transferred to smaller scales. Although it is not shown here, we ran a few decaying turbulence experiments with no forcing to see if this description holds true. Indeed, after initializing the two species with identical random initial conditions (and $KE_i = KE_n = E_B$ at $t = 0$), the runs with very large $\tilde {\alpha }$ showed larger initial magnetic field growth, implying an indirect transfer of energy from the neutrals.