# Impact of the ambipolar diffusion in the structuration of the magnetic Rayleigh Taylor instability with oblique magnetic field

Edouard Callies, V. Guillet, A. Marcowith, Z. Meliani, P. Lesaffre

► To cite this version:

Edouard Callies, V. Guillet, A. Marcowith, Z. Meliani, P. Lesaffre. Impact of the ambipolar diffusion in the structuration of the magnetic Rayleigh Taylor instability with oblique magnetic field. Astronomy & Astrophysics - A&A, 2025, 698, pp.A147. 10.1051/0004-6361/202554074. hal-05128070

HAL Id: hal-05128070

https://hal.science/hal-05128070v1

Submitted on 24 Jun 2025

HAL is a multi-disciplinary open access
archive for the deposit and dissemination of sci-
entific research documents, whether they are pub-
lished or not. The documents may come from
teaching and research institutions in France or
abroad, or from public or private research centers.

L'archive ouverte pluridisciplinaire HAL, est
destinée au dépôt et à la diffusion de documents
scientifiques de niveau recherche, publiés ou non,
émanant des établissements d'enseignement et de
recherche français ou étrangers, des laboratoires
publics ou privés.

<figure><img src="image_1.png" /></figure>

Distributed under a Creative Commons Attribution 4.0 International License
# Impact of the ambipolar diffusion in the structuration of the magnetic Rayleigh Taylor instability with oblique magnetic field

E. Callies<fnref n="1" />, V. Guillet<fnref n="1" />,<fnref n="2" />, A. Marcowith<fnref n="1" />, Z. Meliani<fnref n="3" /> and P. Lesaffre<fnref n="4" />

June 24, 2025

## ABSTRACT

**Aims.** We investigate the impact of ambipolar diffusion on the development of the Rayleigh-Taylor instability (RTI) with an oblique magnetic field in the incompressible limit.

**Methods.** We developed a general bi-fluid framework comprising charges and neutrals with the specific feature of differing gravity between charges and neutrals. We derived the perturbed magnetohydrodynamic (MHD) equations and obtained an analytic dispersion relation for an oblique magnetic field. The growth rate was then evaluated using a numerical integration of the dispersion relation. In particular, we focussed on the anisotropy in the mode growth induced ambipolar diffusion.

**Results.** In contrast to the case of a magnetic field within the interface, an oblique magnetic field is much less restrictive with respect to the wavenumbers that can develop; rather than a sharp cut-off, it results in a selection of preferred scales for mode growth. The presence of a second neutral fluid that is insensitive to gravity tends to amplify the anisotropy in the possible direction of the instability development. In particular, we show that this effect is maximal when the coupling is in an intermediate range between full and no charge-neutral coupling, specifically in the region where the ambipolar diffusion is highest.

## Conclusions.

**Key words.** Magnetohydrodynamics (MHD) – instabilities – Ambipolar Diffusion–Interstellar Medium

## 1. Introduction

1 The Rayleigh-Taylor instability (RTI) is one of the main hy-

2 drodynamic instability that develops in the interstellar medium

3 (ISM). This instability has, for instance, been thoroughly investi-

4 gated in the context of pulsar wind nebula; for instance, by Hes-

5 ter et al. (1996); Bucciantini et al. (2004), Porth et al. (2014),

6 as well as in the context of supernova remnant expansion (e.g.

7 Gull (1975) ) and their interaction with molecular clouds (Jun

8 et al. 1996). Hillier (2016) presented an overview of analytical

9 and numerical applications of RT in astrophysical systems. The

10 onset of the RT instability has also been invoked in the context

11 of molecular cloud formation (Mouschovias et al. 1974) and in

12 gas dynamics structure developments in galactic discs (Rosen &

13 Bregman 1995).

14 Generally, the RT instability is triggered when one fluid is

15 resting on top of a less dense fluid and an unstable equilibrium

16 is established. The resulting dispersion relation is set as follows

17 (Chandrasekhar 1961):

$$ \omega^2 = -\mathcal{A}gk, \quad \mathcal{A} = \frac{\rho_2 - \rho_1}{\rho_1 + \rho_2}, \qquad (1) $$

19 where g is the gravity, k is the wavenumber, and A is the Atwood number, with $\rho_1$ (resp. $\rho_2$) being the density of the light (resp. heavy) fluid. This initial dispersion relation is simple to interpret:

20 whatever the scale, the instability develops and the smaller the

21 scale, the faster the instability develops.

Chandrasekhar (1961) also investigated the impact of the presence of a magnetic field on the development of this instability. In particular, he explored the case of a magnetic field tangent to the interface (i.e. with vanishing obliquity). By deriving the equations within a magnetohydrodynamic (MHD) framework, we can obtain the following dispersion relation (see also the analysis proposed in Priest (1987), Jun et al. (1995)):

$$ \omega^2 = -\mathcal{A}gk + \frac{2(\mathbf{B} \cdot \mathbf{k})^2}{\mu(\rho_1 + \rho_2)}, \qquad (2) $$

with **B** the magnetic field, **k** the wave vector, and $\mu$ the vacuum permeability. The effect of magnetic fields increases at small scales up to the magnetic tension to quench the instability.

Many theoretical studies have been conducted on the magnetic RTI (MRTI). The vast majority assume a magnetic field tangent to the discontinuity. This restriction does not allow for a more general behaviour of the magnetic field to be considered. Vickers et al. (2020) studied the impact of the direction of the magnetic field on the growth rate of the MRTI. The authors proposed a model for an incompressible fluid where the magnetic field has an oblique direction relatively to the interface where the instability develops. Vickers et al. (2020) concluded that the magnetic field's cutoff is not as abrupt as the one given by Eq. 2 when the magnetic field leaves the interface.

In the context of solar physics, the plasma is not totally ionised. Díaz et al. (2012) developed a model with two fluids

<footnote n="1">¹ Laboratoire Univers et Particules de Montpellier (LUPM) Université Montpellier, CNRS/IN2P3, CC72, Place Eugène Bataillon, F-34095 Montpellier Cedex5, France</footnote>

<footnote n="2">² Université Paris-Saclay, CNRS, Institut d'Astrophysique Spatiale, 91405 Orsay, France</footnote>

<footnote n="3">³ Laboratoire Univers et Théories, Observatoire de Paris, Universite PSL, Université de Paris, CNRS, Meudon, France</footnote>

<footnote n="4">⁴ Laboratoire de Physique de l'École Normale Supérieure, ENS, Université PSL, CNRS, Sorbonne Université, Université Paris Cité, 75005, Paris, France</footnote>

<footnote>e-mail: edouard.callies@umontpellier.fr</footnote>
(one charged, the other neutral) both compressible, with a magnetic field lying in the interface. The authors investigated the impact of the magnetically constrained charge motion over the development of instability in the neutral fluid. They found that the impact of the charges and their coupling with the neutrals mainly results in an overall slowing down in the development of the instability such that the magnetic field cutoff was almost no longer present (see also Díaz et al. (2014), Soler & Ballester (2022)).

In the context of ISM physics, the plasma is known to be magnetised and in the galactic disc to be primarily composed of neutrals (H, H₂, or CO depending on the its phase), with a small fraction of charged particles, including ions and grains. ISM is known to often be structured into filaments and is subject to both HD and MHD instabilities and, in particular, the RTI (Jacquet & Krumholz 2011). In the diffuse ISM, gravity from surrounding stars and self-gravity can be safely ignored. However, the medium can be subjected to radiation pressure exerted by a star on the background medium. The impact of this pressure on the medium will be two-fold. On the one hand, the pressure that could push the neutrals is absorbed through the ionisation of the neutrals. Beyond this distance, where the neutrals are ionised, they will no longer be subject to radiation pressure (see chapters 6-7 Osterbrock & Ferland (2006)). On the other hand Krolik (1977) and Jacquet & Krumholz (2011) found that a radiative pressure can act similarly to a steady acceleration on the grains, which (through the equivalent principle) will create an gravity field and trigger the RTI.

The difference in behaviour between the grains and the neutrals is therefore expected to be significant, as only the grains are subject to both the effective gravity and the magnetic field. With this perspective in mind, we aim to study the MRTI by focusing on three key aspects: the oblique nature of the magnetic field (since an interstellar magnetic field is unlikely to be purely contained within the interface), the fact that gravity affects only the grains, and the inclusion of the neutral fluid, which is dominant in this environment but is subject to neither gravity nor the magnetic field. Our final aim in this paper is to examine to what extent ambipolar diffusion induced by the presence of neutrals can impact the structuring of dust filaments created by the RTI.

In Section 2, we develop a theoretical model of an incompressible two-fluid system subject to RTI in the presence of an oblique magnetic field, and we derive the dispersion relation that gives the growth rate of the charges. We then highlight in Section 3, the main differences this relation introduces compared to classical models, namely those of a single-fluid and a two-fluid system, where both fluids experience gravity as well as the impact of the ambipolar diffusion. In Section 4, we focus on the application to the interstellar medium, where in our context, only the charges are subject to gravity. We thus study the impact of the magnetic field on the structuring of anisotropy in the development of the RTI to examine the effect of ambipolar diffusion on this anisotropy. We summarise our results and present our conclusions in Section 5.

## 2. Problem formulation and dispersion relation

### 2.1. Plasma configuration

To keep our approach as general as possible, we chose to consider a bi-fluid system composed of a neutral fluid and a charged fluid. We consider the situation schemed in Fig. 1. The charged (resp. neutral) species have a volumic mass, $\rho_c$ (resp. $\rho_n$). The plane $z = 0$ defines the interface between two regions with dif-

<figure><img src="image_2.png" /><figcaption>Fig. 1: Geometry of the problem. We assume that the x, y plane separates two regions, this is the interface. Each region contains a density $\rho_1$ (for the lower region) and $\rho_2$ (for the upper region) with $\rho_1 < \rho_2$. The entire system is immersed into a magnetic field $\mathbf{B}$ lying in the x, z plane, forming an angle $\theta$ with the x-axis. We assume that the perturbation develops along the wave vector $\mathbf{k}$, which lies in the x, y plane and forms an angle $\alpha$ with the x-axis. The charged (resp. neutral) fluid is subject to the gravity field $\mathbf{g}_c$ (resp. $\mathbf{g}_n$), directed along $-\mathbf{e}_z$.</figcaption></figure>

ferent densities: $\rho_{c2}$ (resp. $\rho_{n2}$) for the upper region, and $\rho_{c1}$ (resp. $\rho_{n1}$) for the lower region. We assume that $\rho_{0c2}/\rho_{0c1} = \rho_{0n2}/\rho_{0n1} = d \ge 1$ is constant and uniform. To better encompass all possible physical situations, including real gravity or radiation pressure where the acceleration generated on each species may differ, we assume that each species is subject to a gravitational force directed along $-\mathbf{e}_z$, denoted as $g_n$ (for the neutrals) and $g_c$ (for the charges). We also consider that this medium is immersed into a magnetic field $\mathbf{B}$, which forms an angle $\theta$ with $\mathbf{e}_x$. The instability develops in the interface ported by the axis $\mathbf{e}_x$ and $\mathbf{e}_y$.

### 2.2. Linear analysis

Our first goal is to port the relations presented in Vickers et al. (2020) (hereafter V20) to a bi-fluid MHD model. Díaz et al. (2012) proposed such kind of analysis but with a magnetic field tangential to the interface. Since we are studying a bi-fluid, we will start with the following equation, using the same framework as Díaz et al. (2012) :

$$ \rho_c \left( \frac{\partial \mathbf{u}_c}{\partial t} + \mathbf{u}_c \cdot \nabla \mathbf{u}_c \right) = -\nabla p_{ce} + \frac{1}{\mu} (\nabla \times \mathbf{B}) \times \mathbf{B} + \rho_c \mathbf{g}_c - \rho_c \nu_{nc}(\mathbf{u}_c - \mathbf{u}_n), \quad (3) $$

$$ \rho_n \left( \frac{\partial \mathbf{u}_n}{\partial t} + \mathbf{u}_n \cdot \nabla \mathbf{u}_n \right) = -\nabla p_n + \rho_n \mathbf{g}_n - \rho_n \nu_{cn}(\mathbf{u}_n - \mathbf{u}_c), \quad (4) $$

where $\mathbf{u}_c$ (resp. $\mathbf{u}_n$) is the velocity of the charged (resp. neutral) fluid, $p_{ce}$ (resp. $p_n$), the pressure of the ion-electron (resp. neutral) fluid, $\nu_{nc}$ (resp. $\nu_{cn}$), the frequency of collision for a given charged (resp. neutral) particle with any neutral (resp. charged) particle. Conservation of momentum yields $\rho_c \nu_{nc} = \rho_n \nu_{cn}$.
132

The induction equation is expressed as:

$$ \frac{\partial \mathbf{B}}{\partial t} = \nabla \times (\mathbf{u}_c \times \mathbf{B}). \qquad (5) $$

133

134

135

136

137

138

To simplify the study, we choose, in a very general case, to assume incompressibility. The justification for the particular case we study in this article is given in Sect. 4.1. The compressible case is much more complex to investigate. It will deserve a specific future work. The incompressible hypothesis is written in general as :

$$ \nabla \cdot \mathbf{u}_c = 0, \qquad (6) $$

$$ \nabla \cdot \mathbf{u}_n = 0. \qquad (7) $$

139

140

141

142

We assume mass conservation and adiabatic energy equation for all the species and the hydrostatic equilibrium verified for the initial pressures (hereafter, the index 0 refers to equilibrium quantities):

$$ \frac{\partial p_{0\mathrm{ce}}}{\partial z} = -g_c \rho_{0c}, \qquad (8) $$

$$ \frac{\partial p_{0n}}{\partial z} = -g_n \rho_{0n}. \qquad (9) $$

143

144

With these considerations we obtain our system of basic equations:

$$ \rho_c \left( \frac{\partial \mathbf{u}_c}{\partial t} + \mathbf{u}_c \cdot \nabla \mathbf{u}_c \right) = -\nabla p_{ce} + \frac{1}{\mu} (\nabla \times \mathbf{B}) \times \mathbf{B} + \rho_c \mathbf{g}_c \\ - \rho_c \nu_{nc} (\mathbf{u}_c - \mathbf{u}_n), \qquad (10) $$

$$ \rho_n \left( \frac{\partial \mathbf{u}_n}{\partial t} + \mathbf{u}_n \cdot \nabla \mathbf{u}_n \right) = -\nabla p_n + \rho_n \mathbf{g}_n - \rho_n \nu_{cn} (\mathbf{u}_n - \mathbf{u}_c), \qquad (11) $$

$$ \frac{\partial p_{ce}}{\partial t} + \mathbf{u}_c \cdot \nabla p_{ce} = 0, \qquad (12) $$

$$ \frac{\partial p_n}{\partial t} + \mathbf{u}_n \cdot \nabla p_n = 0, \qquad (13) $$

$$ \frac{\partial \rho_c}{\partial t} + \mathbf{u}_c \cdot \nabla \rho_c = 0, \qquad (14) $$

$$ \frac{\partial \rho_n}{\partial t} + \mathbf{u}_n \cdot \nabla \rho_n = 0, \qquad (15) $$

$$ \frac{\partial \mathbf{B}}{\partial t} = \nabla \times (\mathbf{u}_c \times \mathbf{B}). \qquad (16) $$

145

### 2.2.1. System linearisation

146

147

148

149

150

151

152

153

154

155

156

To find the dispersion relation, we linearised the previous equations with respect to perturbations written in their Fourier form as $e^{i(k_x x+k_y y-\omega t)}$ in each region separately. Since our study is aimed at studying the growth rate of the instability, we introduced $\Omega = -i\omega$, a positive real that defines our perturbation as $e^{i(k_x x+k_y y)} e^{\Omega t}$. We connected the solution of each region with our boundary conditions. We called **b** the perturbation of the magnetic field and we normalised the magnetic field to $B_0$. Hence, we have the notations $\hat{\mathbf{b}} = \mathbf{b}/B_0$, $\hat{\mathbf{B}} = \mathbf{B}/B_0$. We retained all terms to the first order. With these considerations, we can obtain the following system:

$$ \Omega \rho_{0c} u_{cx} = -ik_x p_{ce} - \rho_{0c} \nu_{nc} (u_{cx} - u_{nx}) \\ + \rho_{0c} c_{Ac}^2 \sin \theta (\partial_z \hat{\mathbf{b}}_x - ik_x \hat{\mathbf{b}}_z), \qquad (17) $$

$$ \Omega \rho_{0c} u_{cy} = -ik_y p_{ce} - \rho_{0c} \nu_{nc} (u_{cy} - u_{ny}) \\ + \rho_{0c} c_{Ac}^2 [\sin \theta (\partial_z \hat{\mathbf{b}}_y - ik_y \hat{\mathbf{b}}_z) \\ + \cos \theta (ik_x \hat{\mathbf{b}}_y - ik_y \hat{\mathbf{b}}_x)], \qquad (18) $$

$$ \Omega \rho_{0c} u_{cz} = -\frac{\partial p_{ce}}{\partial z} - \rho_{0c} \nu_{nc} (u_{cz} - u_{nz}) \\ + \rho_{0c} c_{Ac}^2 \cos \theta (ik_x \hat{\mathbf{b}}_z - \partial_z \hat{\mathbf{b}}_x) - \rho_c g_c, \qquad (19) $$

$$ \Omega \rho_{0n} u_{nx} = -ik_x p_n - \rho_{0n} \nu_{cn} (u_{nx} - u_{cx}), \qquad (20) $$

$$ \Omega \rho_{0n} u_{ny} = -ik_y p_n - \rho_{0n} \nu_{cn} (u_{ny} - u_{cy}), \qquad (21) $$

$$ \Omega \rho_{0n} u_{nz} = -\frac{\partial p_n}{\partial z} - \rho_{0n} \nu_{cn} (u_{nz} - u_{cz}) - \rho_n g_n, \qquad (22) $$

$$ \Omega p_{ce} = -u_{cz} \partial_z p_{0c}, \qquad (23) $$

$$ \Omega p_n = -u_{nz} \partial_z p_{0n}, \qquad (24) $$

$$ \Omega \rho_c = -u_{cz} \partial_z \rho_{0c}, \qquad (25) $$

$$ \Omega \rho_n = -u_{nz} \partial_z \rho_{0n}, \qquad (26) $$

$$ \Omega \hat{\mathbf{b}} = \hat{\mathbf{B}} \cdot \nabla \mathbf{u}_c, \qquad (27) $$

with

157

$$ c_{Ac} \equiv \sqrt{\frac{B_0^2}{\mu \rho_{0c}}}, \qquad (28) $$

the Alfvén velocity in the charged fluid. At this point, we can express the components of the neutral velocity in Eqs. (20-22) as a function of the components of the charges velocity and substituted them into Eqs. (17-19). Then, deriving according to time, we can replace the expressions of the magnetic components and we obtain

158

159

160

161

162

163

$$ \frac{\Omega^2}{\phi_c^2} u_{cx} = -\frac{ik_x \Omega}{\rho_{0c}} (p_{ce} + p_n \psi_n) \\ + c_{Ac}^2 (\hat{\mathbf{B}} \cdot \nabla) [\sin \theta (\partial_z u_{cx} - ik_x u_{cz})], \qquad (29) $$

$$ \frac{\Omega^2}{\phi_c^2} u_{cy} = -\frac{ik_y \Omega}{\rho_{0c}} (p_{ce} + p_n \psi_n) \\ + c_{Ac}^2 (\hat{\mathbf{B}} \cdot \nabla) [\sin \theta (\partial_z u_{cy} - ik_y u_{cz}) \\ + \cos \theta (ik_x u_{cy} - ik_y u_{cz})], \qquad (30) $$

$$ \frac{\Omega^2}{\phi_c^2} u_{cz} = -\frac{\Omega}{\rho_{0c}} \frac{\partial}{\partial z} (p_{ce} + p_n \psi_n) \\ + c_{Ac}^2 (\hat{\mathbf{B}} \cdot \nabla) [\cos \theta (ik_x u_{cz} - \partial_z u_{cx})] \\ + \frac{\partial_z \rho_{0n}}{\rho_{0n}} \psi_n g_n u_{nz} + \frac{\partial_z \rho_{0c}}{\rho_{0c}} g_c u_{cz}, \qquad (31) $$

where we have introduced the adimensional quantities functions of the growth rate $\Omega$ and the collisional frequencies $\nu_{nc}$ and $\nu_{cn}$:

164

165

166

$$ \psi_n \equiv \frac{\nu_{cn}}{\Omega + \nu_{cn}}, \qquad (32) $$

$$ \phi_c \equiv \sqrt{\frac{\Omega + \nu_{cn}}{\Omega + \nu_{cn} + \nu_{nc}}}. \qquad (33) $$
167 To eliminate the pressure terms, we calculated $\partial_z(29) - ik_x \times (31)$ and $k_y(29) - k_x \times (30)$. We obtained the following system:

$$168 \qquad ik_x \phi_c^2 \left[ \frac{\partial_z \rho_{0n}}{\rho_{0n}} \psi_n g_n u_{nz} + \frac{\partial_z \rho_{0c}}{\rho_{0c}} g_c u_{cz} \right] + \Omega^2 \left[ \frac{\partial_z \rho_{0c}}{\rho_{0c}} u_{cx} \right], \quad (34)$$

$$169 \qquad \left[ c_{Ac}^2 \phi_c^2 (\hat{\mathbf{B}} \cdot \nabla)^2 - \Omega^2 \right] (\partial_z u_{cx} - ik_x u_{cz}) = 0. \quad (35)$$

comparing our equation and the one obtained in V20 (their Eq. 13), we can define the same function $f$ by replacing $v_A^2$ by $c_{Ac}^2 \phi_c^2$:

$$f = c_{Ac}^2 \phi_c^2 (\hat{\mathbf{B}} \cdot \nabla)^2 - \Omega^2. \quad (36)$$

At this stage, we can adopt the same assumption as V20, namely, $\partial_z \rho_{0c} = \rho'_{0c} = 0$ and $\partial_z \rho_{0n} = \rho'_{0n} = 0$. This assumption is motivated by two main aspects. First, in our model, we a priori assume that there is a uniform density on either side of the interface and we treat each side of the interface separately. Thus, even if we have a discontinuity of the density across the interface, the density of each side is constant. This simplifies the mathematical treatment and aligns with the physical intuition that the system is initially homogeneous at large scales. Second, as argued bin V20, the instability operates on a length scale much smaller than the characteristic scale over which the density varies. In other words, the density profile does not change significantly over the small spatial regions where the instability develops. This separation of scales justifies neglecting the spatial variation of the density when analysing the instability. With this assumption in place, we recovered the same differential equation for $f$ as V20. This equivalence further validates our approach and confirms that the underlying physics of the instability remains consistent under these assumptions. The charge velocity hence verifies the ordinary differential equation

$$\left(k^2 - \frac{\partial^2}{\partial z^2}\right) f(u_{cz}) = 0. \quad (37)$$

By developing this equation, we can obtain a linear differential equation in $u_{cz}$. The roots of its characteristic equation are then:

$$S = \left\{ \pm k, i \left( \frac{k_x}{\tan \theta} \right) \pm \frac{\Omega}{\sin \theta c_{Ac} \phi_c} \right\}. \quad (38)$$

From this point, we would generally require a complex growth rate, on the one hand, but also call for the phenomenon to be evanescent, as the growth starts at $z = 0$ and $t = 0$. With this in mind, we can explicitly write the expression for the $z$-component of the velocity as:

$$u_{cz} = \begin{cases} A_1 e^{kz} + B_1 e^{m_1 z}, & z < 0, \\ A_2 e^{-kz} + B_2 e^{m_2 z}, & z > 0. \end{cases} \quad (39)$$

We see that the $x$-component of the velocity is also a linear combination of the same exponential. We have:

$$u_{cx} = \begin{cases} i \frac{k_x}{k} A_1 e^{kz} + C_1 e^{m_1 z}, & z < 0, \\ -i \frac{k_x}{k} A_2 e^{-kz} + C_2 e^{m_2 z}, & z > 0. \end{cases} \quad (40)$$

To fully characterise our system, we also need to derive expressions for the $x$ and $z$ components of the velocity of the neutrals.

By using the equation of motion for the neutrals (Eqs. 20-22)

and applying the same transformation as for the charges equations we obtain, under the assumption that $\rho'_{0n} = 0$, the following relation:

$$\left( \partial_z u_{nx} - ik_x u_{nz} - ik^2 u_{nx} - k_x \partial_z u_{nz} \right) = \psi_n \left( -ik^2 u_{cx} - k_x \partial_z u_{cz} \right). \quad (41)$$

With this relation, it becomes clear that the $x$ and $z$ components of the velocity of the neutrals are also a combination of the two same exponential as the $x$ and $z$ components of the velocity of the charges. More precisely, we can write:

$$u_{nz} = \begin{cases} D_1 e^{kz} + \psi_n B_1 e^{m_1 z}, & z < 0, \\ D_2 e^{-kz} + \psi_n B_2 e^{m_2 z}, & z > 0, \end{cases} \quad (42)$$

$$u_{nx} = \begin{cases} i \frac{k_x}{k} D_1 e^{kz} + \psi_n C_1 e^{m_1 z}, & z < 0, \\ -i \frac{k_x}{k} D_2 e^{-kz} + \psi_n C_2 e^{m_2 z}, & z > 0. \end{cases} \quad (43)$$

### 2.2.2. Boundary conditions 210

Given the form of the $x$ and $z$ velocity components for the charges and neutrals, we required eight boundary conditions to determine our dispersion relation. The first condition is meant to consider the continuity of $u_c$. We can also add the incompressibility and the continuity of the magnetic field at the interface to obtain five boundary conditions:

$$[u_{cz}] = 0, \quad (44)$$

$$[u_{cx}] = 0, \quad (45)$$

$$\left[ \frac{\partial u_{cz}}{\partial z} \right] = 0, \quad (46)$$

$$\left[ \frac{\partial^2 u_{cz}}{\partial z^2} \right] = 0, \quad (47)$$

$$\left[ \frac{\partial u_{cx}}{\partial z} \right] = 0, \quad (48)$$

where $[X]$ stands for the jump of the quantity $X$ across the interface. We used the energy equations for charges assuming its continuity, which gives:

$$[g_c \rho_c u_{cz} - \frac{\partial p_{ce}}{\partial t}] = 0. \quad (49)$$

To obtain the eight mandatory boundary conditions, we added the continuity of the $z$ component of $u_n$ and the continuity of the energy equation,

$$[u_{nz}] = 0, \quad (50)$$

$$[g_n \rho_n u_{nz} - \frac{\partial p_n}{\partial t}] = 0. \quad (51)$$

To obtain a proper matrix for the boundary conditions, we replaced the expression of the pressure by using the equation of motion on the $x$ component for each species, which gives:

$$[\rho_{0c}(\phi_c^2 \Omega(\nu_{nc}) + \Omega) - c_{Ac}^2 \phi_c^2 \sin^2 \theta \partial_z^2] u_{cx} + \phi_c^2 \rho_{0c} \Omega \nu_{nc} u_{nx} + ik_x \phi_c^2 g_c \rho_{0c} u_{cz}] = 0, \quad (52)$$

$$[\Omega \rho_{0n}(\nu_{cn} + \Omega) u_{nx} + \rho_{0n} \Omega \nu_{cn} u_{cx} + ik_x g_n \rho_{0n} u_{nz}] = 0. \quad (53)$$
226

227

228

Using the expression of $u_{cx}$, $u_{cz}$, $u_{nx}$, $u_{nz}$, we can express the boundary conditions in a matrix form with each line corresponding respectively to Eqs. 44-51 as follows:

$$ \begin{pmatrix} 1 & 1 & 0 & -1 & -1 & 0 & 0 & 0 & 0 \\ ik_x/k & 0 & 1 & ik_x/k & 0 & -1 & 0 & 0 & 0 \\ k & m_{1-} & 0 & k & -m_{2+} & 0 & 0 & 0 & 0 \\ k^2 & m_{1-}^2 & 0 & -k^2 & -m_{2+}^2 & 0 & 0 & 0 & 0 \\ ik_x & 0 & m_{1-} & -ik_x & 0 & -m_{2+} & 0 & 0 & 0 \\ \alpha_{1c} & \beta_{1c} & \gamma_{1c} & -\alpha_{2c} & -\beta_{2c} & -\gamma_{2c} & \delta_{1c} & -\delta_{2c} & 0 \\ 0 & \psi_n & 0 & 0 & -\psi_n & 0 & 1 & -1 & 0 \\ \alpha_{1n} & \beta_{1n} & 0 & -\alpha_{2n} & -\beta_{2n} & 0 & \delta_{1n} & -\delta_{2n} & 0 \end{pmatrix} \begin{pmatrix} A_1 \\ B_1 \\ C_1 \\ A_2 \\ B_2 \\ C_2 \\ D_1 \\ D_2 \end{pmatrix} = \begin{pmatrix} 0 \\ 0 \\ 0 \\ 0 \\ 0 \\ 0 \\ 0 \end{pmatrix} \quad (54) $$

229

with

$$ \alpha_{1c} \equiv k_x \phi_c^2 \left[ g_c - \left( c_{Ac1}^2 k \sin^2 \theta - \frac{\Omega (v_{nc} + \Omega)}{k} \right) \right], \quad (55) $$

$$ \alpha_{2c} \equiv d^{-1} k_x \phi_c^2 \left[ g_c + \left( d c_{Ac1}^2 k \sin^2 \theta - \frac{\Omega (v_{nc} + \Omega)}{k} \right) \right], \quad (56) $$

$$ \beta_{1c} \equiv k_x \phi_c^2 g_c, \quad (57) $$

$$ \beta_{2c} \equiv d^{-1} k_x \phi_c^2 g_c, \quad (58) $$

$$ \gamma_{1c} \equiv i (\phi_c^2 c_{Ac1}^2 m_{1-}^2 \sin^2 \theta - \Omega^2), \quad (59) $$

$$ \gamma_{2c} \equiv d^{-1} i (d \phi_c^2 c_{Ac1}^2 m_{2+}^2 \sin^2 \theta - \Omega^2), \quad (60) $$

$$ \delta_{1c} \equiv -k_x \phi_c^2 \frac{\Omega v_{nc}}{k}, \quad (61) $$

$$ \delta_{2c} \equiv d^{-1} k_x \phi_c^2 \frac{\Omega v_{nc}}{k}, \quad (62) $$

$$ \alpha_{1n} \equiv -k_x \frac{\Omega v_{cn}}{k}, \quad (63) $$

$$ \alpha_{2n} \equiv d^{-1} k_x \frac{\Omega v_{cn}}{k}, \quad (64) $$

$$ \beta_{1n} \equiv k_x \psi_n g_n, \quad (65) $$

$$ \beta_{2n} \equiv d^{-1} k_x \psi_n g_n, \quad (66) $$

$$ \delta_{1n} \equiv k_x \left( g_n + \frac{\Omega(v_{cn} + \Omega)}{k} \right), \quad (67) $$

$$ \delta_{2n} \equiv d^{-1} k_x \left( g_n - \frac{\Omega(v_{cn} + \Omega)}{k} \right). \quad (68) $$

230

231

232

233

234

235

236

237

We defined $c_{Ac1} \equiv \sqrt{B_0^2/(\mu\rho_{0c1})}$ as the Alfvén velocity in the medium one (the light medium). We replace the analogous $c_{Ac2}$ by $d^{1/2}c_{Ac1}$. We can see that if we suppose that the two fluids are fully decoupled (i.e. $v_{cn} = v_{nc} = 0$) we will obtain a block matrix with the first 6 rows and column giving the matrix obtained by V20 and the last two giving the classic HD dispersion relation: $\Omega_n^2 = \mathcal{A}g_n k$. The next step is to further develop the determinant and after some algebraic manipulations, we obtain:

$$ \Delta_{bi}(k - m_{1-})(k + m_{2+})(m_{1-} + m_{2+}) = 0, \quad (69) $$

238

where

$$ \begin{aligned} \Delta_{bi} \equiv \Delta_{\text{MHD}} \left( 1 - \frac{\Omega(v_{cn} + \Omega)}{\mathcal{A}g_n k} \right) &+ 2k^2 \phi_c^3 g_n c_{Ac1} \frac{v_{nc} v_{cn}}{(v_{cn} + \Omega)^2} \sin\theta \frac{1-d}{1+d^{-1/2}} \\ &+ k \phi_c^2 g_n \frac{v_{nc} v_{cn}}{(v_{cn} + \Omega)^2} (1-d)\Omega, \end{aligned} \quad (70) $$

and

239

$$ \begin{aligned} \Delta_{\text{MHD}} \equiv (1+d)\Omega^3 &+ 2kd^{1/2} \phi_c c_{Ac1} \sin\theta (1+d^{1/2})\Omega^2 \quad (71) \\ &+ \phi_c^2 [-k\tilde{g}(1-d) + 2k^2 dc_{Ac1}^2 (1 - \sin^2 \alpha \cos^2 \theta)]\Omega \\ &- 2k^2 \phi_c^3 c_{Ac1} \tilde{g} \sin\theta \frac{d^{1/2}(1-d)}{1+d^{1/2}}, \end{aligned} $$

with

240

$$ \tilde{g} \equiv g_c + g_n \frac{v_{nc} v_{cn}}{(v_{cn} + \Omega)^2}. \quad (72) $$

The discriminant $\Delta_{\text{MHD}}$ corresponds to the expression of $S$ given by V20 in their Eq. 26, with a correction<fnref n="1" />. The effect of the second fluid on this expression could be summed up by the substitution of $g \to \phi_c^2 \tilde{g}$ and $v_A \to \phi_c c_{Ac1}$.

241

242

243

244

## 3. Results

245

### 3.1. Context and normalisation

246

We start by defining $\Omega_{\text{MHD}}$ as the solution to the dispersion relation obtained by setting $v_{cn} = 0$, namely, the solution obtained in V20. We use as a normalisation for wavenumbers:

247

248

249

$$ k_{\text{cut}} \equiv \frac{(1-d)g_c}{2dc_{Ac1}^2}. \quad (73) $$

This represents the maximum wavenumber value for which the instability can develop with $k$ along $B$ when the magnetic field is tangent ($\theta = 0^\circ$), namely, solving $\omega = 0$ for $\alpha = 0$ in Eq. 2. We define the associated pulsation:

250

251

252

253

$$ \Omega_{\text{cut}} \equiv \sqrt{\mathcal{A} g_c k_{\text{cut}}} = \frac{(1-d)g_c}{\sqrt{2d(1+d)}} c_{Ac1}. \quad (74) $$

All frequencies, if not specified otherwise, have been normalised by this reference frequency. Specifically, we used normalised coupling parameter $v_{nc}/\Omega_{\text{cut}}$. Physically, this normalised quantity is a way to compare the effectiveness of the coupling to the characteristic timescale for the growth of the magnetic RT instability. In order to solve the dispersion relation, we used Python's `fsolve` routine to extract the localisation of the roots. We then selected the one(s) with a positive real part.

254

255

256

257

258

259

260

261

### 3.2. Gravity on the two species

262

The dispersion relation $\Delta_{bi}(\Omega) = 0$ has two solutions for $g_n \neq 0$ and $g_c \neq 0$. In this section, to facilitate identifying in figures the different modes destabilised by the action of the gravity fields associated to $g_n$ and $g_c$ we set $g_n = 10g_c$. We also set $\rho_n = \rho_c$ so that $v_{cn} = v_{nc}$. We show the solution to the dispersion relation in Fig. 2. We normalise the growth rate to $\Omega_{\text{cut}}$. In the case of a vanishing coupling ($v_{nc} = 0$), there are two solutions to the dispersion relation. The first one corresponds to the growth of the instability in the neutral fluid ("HD" mode); it scales as $\sqrt{\mathcal{A}g_n k}$. The second is the one that is obtained by solving the V20 dispersion relation ("MHD" mode). As soon as ion-neutral coupling is activated, the two modes are modified ("HD-modified" and "MHD-modified" modes, respectively) and cannot directly be identified to development of the instability in the neutral or charged fluid, respectively.

263

264

265

266

267

268

269

270

271

272

273

274

275

276

277

278

<footnote n="1">¹ Vickers writes the $\omega^3$ term as $-(d^{1/2} + 1)(1-d)\omega^3$. The correct expression is $-(d^{1/2} + 1)(1+d)\omega^3$. This typo does not impact their results.</footnote>
<figure><img src="image_3.png" /><figcaption>Fig. 2: Solutions to the dispersion relation, normalised by $\Omega_{\text{cut}}$, for various values of the coupling parameter $v_{\text{nc}}/\Omega_{\text{cut}}$: 0 (black lines), $10^{-2}$ (red lines), 1 (blue lines), $10^3$ (green lines), and $10^7$ (yellow lines). We plot the two modes obtained: HD modified (plain lines) MHD modified (dashed lines). We set $\theta = 1^\circ$, $\alpha = 0^\circ$, $d = 0.25$, $g_n = 10g_c$, $\rho_c = \rho_n$.</figcaption></figure>

### 3.2.1. Behaviour of the MHD-modified mode

Without coupling ($v_{\text{nc}} = 0$), the behaviour of the solution for the MHD-modified mode corresponds to the behaviour described by V20: since the magnetic field is not purely tangent to the interface ($\theta = 1^\circ$), the cut-off effect at $k_{\text{cut}}$ disappears. The growth rate reaches a maximum before tending towards an asymptote given by equation 29 in V20. With coupling, the behaviour is different. The asymptote for large $k$ is still present, but its value varies. By expanding our dispersion relation and retaining only the leading orders in $k$, we find that:

$$ \frac{\Omega}{\phi_c}(k \to \infty) \to \Omega_{\text{MHD},k_{\text{co}}}, \qquad (75) $$

where

$$ \Omega_{\text{MHD},k_{\text{co}}} \equiv \frac{g_c}{c_{\text{Ac}1}} \times \frac{(1-d) \sin \theta}{d^{1/2} (1+d^{1/2}) (1 - \sin^2 \alpha \cos^2 \theta)}, \quad (76) $$

is the asymptotic value given by Eq. 29 of V20. Therefore in a general case, we can find the asymptotic value for large $k$ by solving this third-order equation:

$$ \Omega^3 + (v_{\text{nc}} + v_{\text{cn}})\Omega^2 - \Omega_{\text{MHD},k_{\text{co}}}^2 \Omega - v_{\text{cn}}\Omega_{\text{MHD},k_{\text{co}}}^2 = 0. \quad (77) $$

Overall, this equation has no simple analytical solutions. However, we can calculate the solution for extreme value of $v_{\text{nc}}$. Indeed, taking the limit $v_{\text{nc}}/\Omega \ll 1$, $\phi_c$ tends towards 1 and $\Omega$ tends towards $\Omega_{\text{MHD},k_{\text{co}}}$. At the opposite end, for $v_{\text{nc}}/\Omega \gg 1$, $\phi_c$ tends towards $\sqrt{\rho_c/(\rho_c + \rho_n)}$, which explains why in our case of $\rho_c = \rho_n$, the asymptote is reduced by a factor $\sqrt{2}$ with respect to

V20.

In the general case, the behaviour depends on whether $v_{\text{nc}} \ll \Omega$ or $v_{\text{nc}} \gg \Omega$. Indeed, when $v_{\text{nc}} \ll \Omega$, we observe that the growth rate corresponds almost exactly to the growth rate in the uncoupled case. For instance, we can plot a horizontal line for $v_{\text{nc}}/\Omega_{\text{cut}} = 10^{-2}$ (red dashed curves) at $\Omega/\Omega_{\text{cut}} = 10^{-2}$. Above this horizontal line, the growth rate corresponds almost to that of the uncoupled case. At small $k$, the growth rate transitions to a linear growth in $k$. This characteristic effect can be understood from a physical point of view. For this, we compare, on the one hand, the characteristic growth time of the instability and, on the other, the characteristic coupling time. When the former is much shorter than the latter, the instability develops too fast for the coupling to have any significant impact. When the latter is much longer than the former, the instability can only occur by considering a strong coupling. As $v_{\text{nc}}$ increases, the range of $k$ such that $\Omega > v_{\text{nc}}$ shrinks and the entire growth rate eventually transitions to a linear growth in $k$ before reaching the asymptote described above. For an infinite value of the coupling parameter, the MHD-modified mode completely disappears.

### 3.2.2. Behaviour of the HD-modified mode

Without coupling, the solution for the HD-modified mode behaves as predicted by Eq. 1, namely, a scaling as $\sqrt{A}g_n k$. We first considered small values for the coupling parameter ($v_{\text{cn}} \ll \Omega_{\text{cut}}$). At the largest scales ($k \ll k_{\text{cut}}$), the charged and neutral fluids are perfectly coupled and are subject to the effective gravity field $g_{\text{bf}}$ defined as:

$$ g_{\text{bi}} \equiv \frac{\rho_c g_c + \rho_n g_n}{\rho_c + \rho_n}, \qquad (78) $$

which is the weighted average of the gravity experienced by the two species and ultimately assumes the fluids are perfectly coupled. At smaller scales ($k \gg k_{\text{cut}}$), as soon as $\Omega > v_{\text{cn}}$, the neutral fluid decouples from the charged fluid and we recover the trend followed by the neutral gas in the uncoupled case. Now, at large values of the coupling parameter ($v_{\text{cn}} \gg \Omega_{\text{cut}}$), the two fluids are strongly coupled at all $k$ and follow the trend of the charged fluid in the uncoupled case, but with an effective gravity field $g_{\text{bi}}$ and an effective $k_{\text{bi,cut}}$:

$$ k_{\text{bi,cut}} \equiv k_{\text{cut}} \frac{g_{\text{bi}}}{g_c}. \qquad (79) $$

For the largest value plotted, we observe a dip given by:

$$ \Omega = \sqrt{\frac{\rho_c + \rho_n}{\rho_c} \frac{g_{\text{bi}}}{c_{\text{Ac}1}}} \times \frac{(1-d) \sin \theta}{d^{1/2} (1+d^{1/2}) (1 - \sin^2 \alpha \cos^2 \theta)}. \quad (80) $$

As $v_{\text{cn}}$ increases, the range of $k$ such that $\Omega$ follows this equation expands, so that when the coupling becomes infinite, we can obtain an asymptote for large values of $k$.

### 3.2.3. Discussion

This study bears similarities to the one conducted by Díaz et al. (2012), investigating MRTI with two fluids, along with the effect of magnetic field, but having it assigned to the interface. The authors also obtained two modes: a 'neutral' one and an 'ionic' one (see Fig.5 of their paper). As in the present study, the presence of coupling with the charged fluid does not completely suppress the instability in the neutral mode (which we refer here as to the
348 'HD' mode). Just as we have found, the presence of coupling with the charged fluid does not completely suppress the insta- bility in the neutral mode (i.e. our 'HD' mode). Moreover, as 349 the coupling increases, these authors also observed that the HD 350 mode tends to weaken, similarly to the behaviour shown by the 351 yellow curve in Fig. 2. It is less straightforward to compare the 352 case of charge modes in Diaz's work with the MHD mode here, 353 since the orientation of the magnetic field differs between the 354 two cases. However, it is worth noting that they also observe an 355 attenuation of the charge mode as the coupling increases. 356

358 3.3. Case of g<sub>n</sub> = 0

359 3.3.1. Dispersion relation

360 If we assume that the radiative forcing is only impacting the 361 charges, we can set g<sub>n</sub> = 0. With some algebraic manipulation, 362 the matrix of the boundary conditions can be expressed as:

$$
\begin{pmatrix}
1 & 1 & 0 & -1 & -1 & 0 \\
ik_x/k & 0 & 1 & ik_x/k & 0 & -1 \\
k & m_{1-} & 0 & k & -m_{2+} & 0 \\
k^2 & m_{1-}^2 & 0 & -k^2 & -m_{2+}^2 & 0 \\
ik_x & 0 & m_{1-} & -ik_x & 0 & -m_{2+} \\
\tilde{\alpha}_1 & \tilde{\beta}_1 & \tilde{\gamma}_1 & -\tilde{\alpha}_2 & -\tilde{\beta}_2 & -\tilde{\gamma}_2
\end{pmatrix}
\begin{pmatrix}
A_1 \\
B_1 \\
C_1 \\
A_2 \\
B_2 \\
C_2
\end{pmatrix}
= 0, \quad (81)
$$

363 where the greek tilde letters represents the same coefficients as 364 the one of the V20 matrix with the difference that V<sub>A</sub> is here 365 c<sub>Ac</sub>&phi;<sub>c</sub>, g is g<sub>c</sub>&phi;<sub>c</sub><sup>2</sup> and the sign in m<sub>2+</sub> and m<sub>1-</sub> are inverted, since 366 we choose to study &Omega; = -i&omega; rather than &omega;:

$$
\tilde{\alpha}_1 \equiv k_x \phi_c^2 g_c - \frac{k_x}{k} (\phi_c^2 c_{Ac1}^2 k^2 \sin^2 \theta - \Omega^2), \quad (82)
$$

$$
\tilde{\alpha}_2 \equiv d^{-1} \left( k_x \phi_c^2 g_c + \frac{k_x}{k} \left( \phi_c^2 d c_{Ac1}^2 k^2 \sin^2 \theta - \Omega^2 \right) \right), \quad (83)
$$

$$
\tilde{\beta}_1 \equiv k_x \phi_c^2 g_c, \tag{84}
$$

$$
\tilde{\beta}_2 \equiv d^{-1} k_x \phi_c^2 g_c, \qquad (85)
$$

$$
\tilde{\gamma}_1 \equiv i (\phi_c^2 c_{Ac1}^2 m_{1-}^2 - \sin^2 \theta - \Omega^2), \quad (86)
$$

$$
\tilde{\gamma}_2 \equiv id^{-1} (\phi_c^2 d c_{Ac1}^2 m_{2+}^2 \sin^2 \theta - \Omega^2). \quad (87)
$$

367 Therefore the dispersion relation is similar to the one derived 368 V20, but with the substitution mentioned above. Thus, we have 369 to solve:

$$
\Delta_{\text{MHD}} = 0, \text{ with } \tilde{g} \equiv g_c. \tag{88}
$$

370 This determinant is slightly more complex than the one ob-
371 tained by V20 since both $\tilde{g}$ and $\phi_c$ depend upon $\Omega$. In order to
372 centralise the discussion on the coupling effect, we decide to use
373 normalised quantity $k/k_{\text{cut}}$ and $\Omega/\Omega_{\text{cut}}$. With those quantities we
374 can rewrite Eq. (88) as:

$$
0 = \left(\frac{\tilde{\Omega}}{\tilde{\phi}_c}\right)^3 + \left(\frac{\tilde{\Omega}}{\tilde{\phi}_c}\right)^2 \tilde{k} \sqrt{2} a \sin\theta \\
+ \left(\frac{\tilde{\Omega}}{\tilde{\phi}_c}\right) \left(-\tilde{k} + \tilde{k}^2 (1 - \sin^2 \alpha \cos^2 \theta)\right) - \tilde{k}^2 \frac{\sqrt{2} \sin\theta}{a}, \quad (89)
$$

375 using the notation $a \equiv (1 + d^{1/2})/(1 + d)^{1/2}$.

376 We plot the solution of the dispersion relation normalised 377 by Ωcut in Fig. 3 for various values of the coupling parameter. 378 To see the impact of the friction of the neutral fluid, we 379 set ρ₀c/ρ₀n = 10⁻⁴. Without any coupling, we obtained the

<figure><img src="image_4.png" /><figcaption>Fig. 3: Solutions to the dispersion relation for various values of the coupling parameters. The solutions are normalised by Ω<sub>cut</sub>. We set g<sub>n</sub> = 0, θ = 1°, α = 0°, d = 0.25, and ρ<sub>0c</sub>/ρ<sub>0n</sub> = 10<sup>-4</sup>.</figcaption></figure>

same solution as before. The solution for a strong coupling (v<sub>nc</sub> ≫ Ω<sub>cut</sub>) is the same as the one without coupling but with a global attenuation.

The asymptotic values observed follow the same behaviour
as we see in Eq. 75. We can define k<sub>∞</sub> as an approximation of the
value of k from which this asymptote occurs. An approximation
of the value of k<sub>∞</sub> can be given by k<sup>2</sup> term dominating the k term
in the Ω term of the equation above. This yields

$$
k_{\infty} \equiv \frac{k_{\text{cut}}}{1 - \sin^2 \alpha \cos^2 \theta}. \qquad (90)
$$

Therefore, with α = 0° and θ = 10°, the plateau is reached
for k ≫ k∞. We can define

$$
\Omega_{\infty} \equiv \Omega(k \gg k_{\infty}). \qquad (91)
$$

Even though it is challenging to fully grasp the impact of the
coupling just by examining the dispersion relation, we can define
$\tilde{\Omega} = \Omega/\phi_c$, so that this variable satisfies V20 dispersion relation.
In particular, we can first look at the asymptotic expressions for
large and small values of the coupling parameters:

For small values of v<sub>nc</sub> ≪ Ω, the term φ<sub>c</sub> tends towards 1,
which means that the solution of our equation approaches Ω<sub>MHD</sub>;
if the coupling is much slower than the development of the in-
stability, the charged fluid undergoes the instability without be-
ing influenced by the neutrals and therefore follows the equa-
tions of V20. Conversely, when v<sub>nc</sub> ≫ Ω, the term φ<sub>c</sub> tends to-
wards √(ρ<sub>c</sub>/(ρ<sub>c</sub> + ρ<sub>n</sub>)). The equations are then identical to those
of the decoupled case but with an overall attenuation factor:
since the coupling is strong, the charges must first drag the neu-
trals, and thus the gravity g<sub>c</sub> is diluted. These two cases can be
decomposed into two regimes. Before k̃ = 1, the growth follows
k<sup>1/2</sup>, while beyond this point, it reaches the previously described
asymptotes. In a more general case, at large scales (i.e. small
values of k), ions are fully coupled to neutrals and, thus, we find
a growth as k<sup>1/2</sup>, following the perfect coupling model. As soon
as we enter the range where v<sub>nc</sub> ≪ Ω, friction loses its signif-
icance and the growth rate eventually converges to that of the
collision-less case.
### 3.3.2. Anisotropy of grown instabilities

The magnetic field naturally produces anisotropic filamentary structures (Goldreich & Sridhar 1995). The objective of this section is therefore to quantify this anisotropy in the magnetic RT instability and the impact of the ambipolar diffusion on the anisotropy.

In light of the discussions linked to Eq. (88), it is important to first recall the behaviour of the instability in the MHD framework. Indeed, the more complex case with coupling appears to be similar to this initial framework but with a few adjustments. In Fig. 4, we plot the map of $\Omega/\Omega_{\text{cut}}$ values without coupling in the $(k_x/k_{\text{cut}}, k_y/k_{\text{cut}})$ plane. From left to right, the inclination angle of the magnetic field gradually increases from 1° to 60°. A clear trend in growth rate anisotropy can be seen when comparing the first and last maps. Indeed, while the last map is almost uniform, the first map shows very large contrast differences. Between these two maps, the contrast decreases uniformly: the more constraining it becomes, then the larger the magnetic tension. The shape of the map also shows that the situation is not symmetric in $k_x$ and $k_y$: for $k_x/k_{\text{cut}} > 1$, we observe a strong slowdown in the growth rate for small $\theta$, which is not visible for $k_y/k_{\text{cut}} > 1$.

To quantify this anisotropy, we define a function $A_\alpha$:

$$A_{\alpha} \equiv \frac{\Omega(\alpha = 90^{\circ}) - \bar{\Omega}}{\Omega(\alpha = 90^{\circ}) + \bar{\Omega}}, \qquad (92)$$

$$\bar{\Omega} \equiv \frac{2}{\pi} \int_{0}^{\pi/2} \Omega(\alpha) \, d\alpha, \qquad (93)$$

where $\bar{\Omega}$ is the mean value of $\Omega$ over $\alpha$. Physically, this choice is guided by the understanding that in a given medium, what we observe is the variation or contrast in growth across different directions. By defining the anisotropy and averaging it over all values of $\alpha$, we can capture the contrast differences in each direction and better understand which modes can genuinely develop. When the anisotropy function is close to 1, it indicates strong directional dependence: only specific values of $\alpha$, primarily around $\pi/2$, can effectively develop. This highlights a pronounced anisotropic effect where growth is concentrated in certain directions. On the other hand, an anisotropy value close to 0 implies nearly uniform development in all directions, suggesting that the magnetic field exerts minimal influence on the directional growth patterns.

We plot the anisotropy $A_\alpha$ in Fig. 5 as function of $k/k_{\text{cut}}$ for various values of $\theta$. The initial observation is the distinct variation in behaviour across angles. For $\theta = \pi/2$, the anisotropy remains consistently at the level of zero, indicating that the magnetic field is unable to influence the RTI in any way other than isotropically. In other words, there is no directional dependence when $\theta$ is perpendicular. In contrast, as $\theta$ drops towards 0, we observe an increase in the anisotropy parameter, reaching its maximum at $\theta = 0$. This confirms that when the magnetic field is aligned with the instability (parallel configuration), it has a pronounced impact on the growth in different directions. The role of $k/k_{\text{cut}}$ is also noteworthy. For values below $k/k_{\text{cut}} = 1$, anisotropy remains almost negligible. However, as $k/k_{\text{cut}}$ increases, the anisotropy progressively rises, eventually approaching an asymptotic value at larger $k/k_{\text{cut}}$ values, indicating that the influence of the magnetic field on directional growth stabilises at high $k/k_{\text{cut}}$.

Next, we turn to the coupling effect and in Fig. 6, we plot the solutions to the dispersion relations for various values of the coupling parameter and various values of the angle $\alpha$. We plot

the two extreme cases, $v_{nc}/\Omega = 0$ (i.e. MHD mode in black) and $v_{nc}/\Omega = \infty$ (i.e. charged and neutral fluids totally coupled in green). Between these two limits, we can draw two values of the coupling parameter. As stated above, we can see the impact of the coupling parameter on the form of the curve. Following, for instance, the red curve: at small scales, the charges do not interact with the neutrals, since the growth time of the instability is much shorter than the ion-neutral collision time ($\Omega > v_{nc}$). As a result, the coupling cannot be established during the development of the instability. Then, when $\Omega < v_{nc}$, the charges start to feel the friction with the neutrals. Finally, at even larger scales, when $\Omega < v_{cn}$, the neutrals begin to respond to the charges, and the two fluids behave as a single coupled fluid. In addition to this behaviour, this figure provides a first insight into the difference in anisotropy: for large values of $v_{nc}$, the curves are identical up to an attenuation factor, meaning the anisotropy remains the same. Since decoupling occurs at increasingly larger scales as $v_{nc}$ decreases and given that the anisotropy caused by the magnetic field only appears at small scales, no impact on anisotropy will be observed for small $v_{nc}$. Between these scales, for the $v_{nc}$ values plotted, we observe that at $k/k_{\text{cut}}$, the gap between the curve for $\alpha = \pi/2$ and $\alpha = 0$ is larger than in the MHD case. Therefore, to address the anisotropy, we go on to consider $v_{nc}/\Omega_{\text{cut}} = 10$.

We plot in Fig. 7 (similarly to what is shown in Fig. 5), the anisotropy as a function of $k/k_{\text{cut}}$ for different values of $\theta$. We see that regardless of the angle $\theta$ or the value of $k/k_{\text{cut}}$, the anisotropy curve with coupling consistently lies above that of the MHD case. This confirms our initial observations: coupling tends to enhance the anisotropy induced by the magnetic field. The asymptotic behaviour of the anisotropy could be explained by Fig. 6: for larger $k/k_{\text{cut}}$, whatever the value of the angle $\alpha$ the growth rate reaches an asymptotic value given by Eq. (75), which leads to a plateau for the anisotropy.

To conclude this analysis of the anisotropy, we plotted the anisotropy as a function of the coupling parameter for various values of $k/k_{\text{cut}}$ and $\theta$ in Fig. 8. First, we note that the larger the $k/k_{\text{cut}}$, the stronger the anisotropy, which can be anticipated based on the previous figures. Specifically, for high of values $k/k_{\text{cut}}$ with $\theta \approx 0$, we find a constant anisotropy equal to 1. Conversely, for $k/k_{\text{cut}}$ values below 1, we observe minimal anisotropy. A noteworthy point across this plot is the particular dependence of anisotropy on the coupling parameter. Except for the previously described extreme cases, all curves exhibit a similar behaviour: at extreme values of $v_{nc}/\Omega_{\text{cut}}$, they reach flat asymptotes, while within the previously specified range where decoupling and thus ambipolar diffusion occur, anisotropy reaches a peak. This figure therefore clearly highlights the significant role that ambipolar diffusion can play in enhancing anisotropy.

If we shift from a theoretical framework to a more astrophysical context, we must consider that the magnetic field, the effective gravity, and the coupling strength between neutrals and charges are no longer variables but rather fixed parameters. The range in which anisotropy manifests significantly holds almost independently of the observational scale or the inclination angle of the magnetic field.

## 4. The RTI in the magnetised, partially ionised ISM

### 4.1. Context and physical justification of the model

The astrophysical context we aim to study here is specifically that of the ISM. We restricted our analysis to the multiphase, diffuse and very weakly ionised ISM. In this context, we con-
<figure><img src="image_5.png" /><figcaption>Fig. 4: Maps of the value of $\Omega/\Omega_{\text{cut}}$ without coupling ($v_{\text{nc}} = v_{\text{cn}} = 0$) in the $(k_x/k_{\text{cut}}, k_y/k_{\text{cut}})$ plane. The magnetic field inclination $\theta$ increases from left to right, taking the values $1^\circ, 10^\circ, 30^\circ$, and $60^\circ$.</figcaption></figure>

<figure><img src="image_6.png" /><figcaption>Fig. 5: Anisotropy $A_\alpha$ as a function of $k/k_{\text{cut}}$, for various values of the magnetic field inclination $\theta$ (different line style).</figcaption></figure>

<figure><img src="image_7.png" /><figcaption>Fig. 6: Solutions to the dispersion relation for various values of the coupling parameter $v_{\text{nc}}/\Omega$: 0 (black, MHD limit), 1 (red), $10^2$ (blue), and $\infty$ (green, perfect coupling), for two values of the angle $\alpha$: $90^\circ$ (solid) and $0^\circ$ (dashed). The solutions are normalised by $\Omega_{\text{cut}}$. We set $g_n = 0$, $\theta = 1^\circ$, $d = 0.25$, $\rho_{0c}/\rho_{0n} = 10^{-4}$ as in Fig. 3.</figcaption></figure>

sidered three species: neutrals (either H or H$_2$), ions (either H$^+$ or C$^+$), and dust grains. We therefore aimed to apply the model developed above to this type of medium. The incompressibility assumption makes sense while considering diffuse ISM phases. In atomic (warm and cold neutral) phases that are the subject of our interest here, turbulence properties is deduced from the 21 cm Hi line emission (Heiles & Troland 2005). In the cold neutral medium (CNM), it is found that the turbulent velocity dispersion has a median value of 2.8 km/s, which implies a moderately supersonic turbulence with sonic Mach numbers $M_s \sim 3$. In the warm neutral medium (WNM), higher turbulent velocities induce slightly subsonic to transonic flows because the sound velocity is also higher $c_s \sim 10$ km/s (Saury et al. 2014). The total magnetic field derived from Zeeman splitting effect in Hi emission (WNM) or absorption (CNM) lines have a median values of $B = 6$ µG in the CNM, while measurements provide less accurate estimates of $B_{\text{tot}} = 5 - 10$ µG in the WNM (Heiles & Crutcher 2005). We find typical Alfvénic Mach numbers are $M_a \sim 3$ in the CNM and $\lesssim 1$ in the WNM. We conclude due to these moderate Mach numbers that CNM (but also likely WNM) show a moderate compressibility. In any case, the compressibility applies mostly at scales close to the turbulence injection

scale. As scales of the ambipolar diffusion are concerned they are below the sonic scale where the gas turbulent velocity matches the ambient sound speed. Hence, we are dealing with a part of the turbulent cascade lying in the incompressible regime (see e.g. Beattie et al. (2024)).

## 4.2. Modelling ions and dust as a single charged fluid 560

In the diffuse ISM, we need to consider a priori a tri-fluid model for our system composed of neutrals, ions, and dust. If we want to use our bi-fluid model and the corresponding results presented in Section 3, we need to adopt a simpler approach by gathering all the charges into a single charged fluid. This approach is motivated by Zaqarashvili et al. (2011) who demonstrated how to unify a H$^+$/He$^+$ plasma in a single fluid. The condition is that $\Omega$, the gyrofrequencies of each charged species, is much larger than the collision frequency, $v$, between distinct charged species. In terms of physics, this means that the charged species are
<figure><img src="image_8.png" /><figcaption>Fig. 7: Variation of the anisotropy as a function of $k/k_{\text{cut}}$, for various values of $\theta$ (plain $\theta = 1^{\circ}$, dashed $\theta = 10^{\circ}$, dash-dotted $\theta = 30^{\circ}$), for $\nu_{nc}/\Omega_{\text{cut}} = 0$ (black) and $\nu_{nc}/\Omega_{\text{cut}} = 10$ (red).</figcaption></figure>

<figure><img src="image_9.png" /><figcaption>Fig. 8: Plot of the anisotropy for various values of $k/k_{\text{cut}}$ different colours (blue $k/k_{\text{cut}} = 5$, red $k/k_{\text{cut}} = 0.5$) and for various values of $\theta$ different line styles (plain $\theta = 1^{\circ}$, dashed $\theta = 10^{\circ}$, dash-dotted $\theta = 30^{\circ}$) as a function of the coupling parameters.</figcaption></figure>

coupled due to a strong magnetic field and not due to collisions. We describe below the different components of this charged fluid (i.e. ions, PAHs, and dust) and show that the conditions are satisfied to gather them into a single charged fluid, described with the index c.

Dust grains in the diffuse ISM cover a size distribution ranging from subnanometer particles (i.e. polycyclic aromatic hydrocarbons, PAHs) to submicrometer grains (Hensley & Draine 2023). Through the competition between electron collection, recombination of ions on the surface of grains, and photoemission of surface electrons, dust grains acquire an electric charge of $q = Ze$ (Weingartner & Draine 2001b). In the presence of a significant radiation field ($G_0 > 1$), the grain charge is on average positive and increases with the grain radius, a. Through the Lorentz force, dust grains of a mass, m, are coupled to the magnetic field of intensity, B, on a timescale of $\tau_{\text{Lar}} = 1/\Omega_{\text{Lar}}$, where

$\Omega_{\text{Lar}} = qB/mc$ is the Larmor frequency and c is the speed of light. Through collisions with gas particles, they are also coupled to the neutral gas on a timescale of $\tau_{\text{drag}} = m/(\rho_n\sigma c_s)$ where $\sigma = \pi a^2$ is the grain cross-section and $c_s$ is the sound velocity in the gas. The ratio of these two timescales, the Hall factor $\Gamma = \Omega_{\text{Lar}}\tau_{\text{drag}}$, expresses in the CNM:

$$ \Gamma^{\text{CNM}} \approx 20Z \left( \frac{n_n}{50 \text{ cm}^{-3}} \right)^{-1} \left( \frac{a}{0.1 \mu\text{m}} \right)^{-2} \left( \frac{B}{6 \mu\text{G}} \right) \left( \frac{T}{100 \text{ K}} \right)^{-1/2} \quad (94) $$

It is much larger than 1 for all sizes relevant to the diffuse ISM, implying that dust grains are all strongly coupled to the magnetic field and follow on average the dynamic of ions (Guillet et al. 2007). We could argue that neutral grains can not be coupled to the magnetic field. This remark is particularly relevant for PAHs because the smaller the grains, the larger the probability for the grain to be electrically neutral. Actually, grain charges tend to fluctuate. Their neutral state is only temporary, before they acquire a charge again. This stochastic process is characterised by a timescale, $\tau_Z$, that is much shorter than $\tau_{\text{drag}}$ (Yan et al. 2004, their Fig. 1). This implies that grains do not have the time to couple to the gas during their neutral phase and remain coupled to the magnetic field whatever the value of their instantaneous charge. It is therefore justified to assume that all dust grains (including PAHs) always remain coupled to the magnetic field whatever their charge state. We will therefore consider ions and dust grains as a single fluid that is 'frozen' in the magnetic field lines.

The properties of our charged fluid is the sum of the properties of ions, PAHs, and large (5 – 250 nm) dust grains. For the CNM and WNM, its total mass is dominated by dust grains: $\rho_d = 6.410^{-3}\rho_n$, $\rho_P = 7.10^{-4}\rho_n$ (Hensley & Draine 2023), and $\rho_i = n_i m_i$ (see Table 1 for the values of $n_i$ and $m_i$ in each phase). The collisional cross-section carried by PAHs is six times higher than the cross-section carried by grains (Hensley & Draine 2023). Moreover, the electric dipole induced by a charged PAH on a gas particle increases the collisional cross-section and makes it independent of the collisional velocity or PAH size (the Langevin cross-section $\langle\sigma v\rangle_L = 2.1 \times 10^{-9} \text{ cm}^3/\text{s}$, see Jean et al. (2009); Brahimi et al. (2020)). The same cross-section applies to ions. Next, we compute the charged-neutral collisional frequency, $\nu_{nc}$, as

$$ \rho_n \nu_{cn} = \rho_c \nu_{nc} = \rho_i \nu_{ni} + \rho_P \nu_{nP} + \rho_d \nu_{nd}, \quad (95) $$

$$ \rho_i \nu_{ni} = n_n \frac{m_i m_n}{m_i + m_n} n_i \langle \sigma v \rangle_L, \quad (96) $$

$$ \rho_P \nu_{nP} \approx \rho_n n_P \langle \sigma v \rangle_L, \quad (97) $$

$$ \rho_d \nu_{nd} \approx \rho_n \langle n\sigma \rangle_d c_s, \quad (98) $$

where $n_n$ (resp. $n_i$, $n_P$ and $n_d$) is the density of neutral particles (resp. ions, PAHs, and dust grains), $m_n$ (resp. $m_i$) are the mass of neutral particles (resp. ions), and $\langle n\sigma \rangle_d = 3.10^{-21} n_H$ (Weingartner & Draine 2001a). We find the contribution of PAHs and dust grains to be negligible compared to the ion term both for the CNM and WNM and, thus, we can neglect them. We then obtain:

$$ \nu_{cn} \approx \frac{\rho_i}{\rho_n} \nu_{ni} \approx \frac{\rho_i}{m_i + m_n} \langle \sigma_{cn} \nu \rangle, \quad (99) $$

$$ \nu_{nc} = \frac{\rho_n}{\rho_i + \rho_P + \rho_d} \nu_{cn}. \quad (100) $$
<table><caption>Table 1: Values of the parameters in the warm neutral medium (WNM) and cold neutral medium (CNM)</caption><thead><tr><th>Phase</th><th>WNM</th><th>CNM</th></tr></thead><tbody><tr><td><i>B</i><sub>0</sub> [G]</td><td>5 · 10<sup>-6</sup></td><td>6 · 10<sup>-6</sup></td></tr><tr><td><i>n</i><sub>n</sub> [cm<sup>-3</sup>]</td><td>0.2</td><td>50</td></tr><tr><td><i>n</i><sub>i</sub> [cm<sup>-3</sup>]</td><td>7 · 10<sup>-3</sup></td><td>2.4 · 10<sup>-2</sup></td></tr><tr><td>Ions</td><td>H<sup>+</sup></td><td>C<sup>+</sup></td></tr><tr><td>Neutral</td><td>H+He</td><td>H+He</td></tr><tr><td><i>m</i><sub>i</sub> [m<sub>H</sub>]</td><td>1</td><td>12</td></tr><tr><td><i>m</i><sub>n</sub> [m<sub>H</sub>]</td><td>1.21</td><td>1.21</td></tr><tr><td><i>T</i>(K)</td><td>8000 K</td><td>50 K</td></tr><tr><td><i>c</i><sub>Ac</sub> [cm·s<sup>-1</sup>]</td><td>1.2 · 10<sup>7</sup></td><td>1.55 · 10<sup>6</sup></td></tr></tbody></table>

Notes. $B_0$ is the magnetic field strength, $n_n$ (resp. $n_c$) is the density in the neutral (resp. charged) fluid. The temperatures are mean values. (Brahimi et al. 2020).

### 4.3. Effective gravity due to radiation pressure

In the diffuse ISM, we can a priori assume that the gravity of the surrounding stars, as well as self-gravity, can be ignored. Several mechanisms could justify the appearance of a RTI by creating an effective gravity. The one we aim to explore relies on radiation pressure. We consider the following framework: a star irradiates an ISM cloud. In doing so, the cloud will be subject to radiation. Radiation can first push the neutrals and simply ionise them, leading to them being absorbed. Beyond a distance where all this radiation is absorbed, we can assume that the neutrals will no longer experience any radiation pressure (see chapters 6-7 in Osterbrock & Ferland (2006)). The remaining radiation pressure will therefore primarily affect dust particles. Jacquet & Krumholz (2011) demonstrated that this pressure generates a steady acceleration of the medium a that, through the equivalence principle, is equivalent to a gravity field of $g = -a$ exerting its effects in an inertial frame, which can trigger an RT instability. Porth et al. (2014) also explored this in the case of the Crab Nebula, where supernova ejecta formed a shell subjected to a matter wind from the pulsar, enabling the appearance of an RTI.

From the energy conservation, the total energy emitted by dust is equal to the total starlight energy absorbed by dust. In the neighbourhood of a bright star, the interstellar radiation field is highly anisotropic. The transfer rate of momentum by radiation pressure acting on dust is therefore proportional to the total power, $L_d$, emitted by dust per H (in erg/s/H), here computed for the interstellar radiation field (ISRF, $G_0 = 1$). The same is true for PAHs ($R_P$), but not for ions which do not absorb starlight. The Hensley & Draine (2023) diffuse dust model gives $L_d = 3.4 \times 10^{-24}$ erg/s/H and $L_P = 1.6 \times 10^{-24}$ erg/s/H. The effective gravitational field $g_c$ for the charged fluid is then computed from the conservation of its momentum:

$$ \rho_c g_c = G_0 \frac{(L_P + L_d) n_H}{c}. \qquad (101) $$

We find

$$ g_c \approx 5 \times 10^{-9} G_0. \qquad (102) $$

### 4.4. Application to the CNM and WNM

Using Table 1 and Eq. 73, we can estimate the cutoff wavelength in each phase:

$$ L_{\text{cut}}^{\text{CNM}} \equiv 2\pi/k_{\text{cut}}^{\text{CNM}} \approx \frac{0.5}{G_0} \text{kpc}, \qquad (103) $$

$$ L_{\text{cut}}^{\text{WNM}} \equiv 2\pi/k_{\text{cut}}^{\text{WNM}} \approx \frac{100}{G_0} \text{kpc}. \qquad (104) $$

It then appears that this length for the CNM is typically on 670 the order of tens of parsecs for a $G_0$ value of 10. With these 671 considerations, the only instabilities that can grow are those 672 described by the asymptotic behaviour (i.e for $k \gg k_\infty$).

Similarly, we can estimate the ratio of the collisional frequency between the neutral and charged fluids to the cut-off frequency of RT instability (Eq. 74) as:

$$ (v_{\text{nc}}/\Omega_{\text{cut}})^{\text{CNM}} = \frac{1.1 \cdot 10^6}{G_0}, \qquad (105) $$

$$ (v_{\text{nc}}/\Omega_{\text{cut}})^{\text{WNM}} = \frac{6.0 \cdot 10^4}{G_0}. \qquad (106) $$

It thus becomes clear that unless there is an extremely high 678 radiation pressure, coupling will always occur on timescales 679 much shorter than those over which the RT instability could 680 develop. As a result, with respect to this instability, the two 681 fluids will behave as if they are perfectly coupled. Therefore, we 682 cannot expect any observable signature of ambipolar diffusion 683 in the development of this instability: what could be observed is 684 merely strong scale-independent anisotropy.

We can go on to estimate the growth time of the instability:

$$ t^{\text{CNM}} \equiv \frac{1}{\Omega_{\infty}^{\text{CNM}}(\alpha = \frac{\pi}{2})} \approx \frac{1.4}{G_0} \text{Myr}, \qquad (107) $$

$$ t^{\text{WNM}} \equiv \frac{1}{\Omega_{\infty}^{\text{WNM}}(\alpha = \frac{\pi}{2})} \approx \frac{18}{G_0} \text{Myr}. \qquad (108) $$

With realistic values (for $G_0$ (1-100)), we expect growth times 688 on the order of, at best, 0.1 Myr. To better grasp the order 689 of magnitude of this growth time, we applied it to the case 690 of the environment of the star Merope, a star in the Pleiades, 691 which is surrounded by an interstellar medium cloud from the 692 Taurus region. The interaction between these two systems has 693 been documented by (Gibson & Nordsieck 2003). We can thus 694 recall that the Pleiades came into contact with this cloud at 695 $\tau \approx 0.1$ Myr. Moreover, the radiation field from Merope can be 696 estimated as $G_0 = 20$. With this radiation field, we obtain a 697 growth time for the instability comparable to the characteristic 698 interaction time with the star. Therefore, the instability could, in 699 principle, develop. Two limitations must be taken into account: 700 on the one hand, the calculated time corresponds to the fastest- 701 growing mode ($\alpha = 90^\circ$), meaning that few other modes (with 702 different angle values $\alpha$) have sufficient time to fully develop; 703 on the other hand, since the two timescales are comparable, the 704 actual development of this mode remains rather marginal.

The conclusion drawn here is, at first glance, rather negative. 705 It is therefore important to recall that the analysis presented here 706 is a linear analysis. However, as pointed out by Soler & Ballester 707 (2022), ambipolar diffusion can play a significant role when the 708 instability enters a nonlinear regime with large-scale effects. In- 709 deed, Khomenko, E. et al. (2014) found, using a 2.5D numeri- 710 cal simulation, that the growth rate of small-scale modes during 711
714 the nonlinear phase can be significantly higher with ambipolar diffusion than without. Other, more recent works also support this view (Popescu Braileanu, B. et al. 2021a,b). These results highlight the non-negligible impact of ambipolar diffusion in the structuring of the instability during the nonlinear phase.

## 5. Conclusion

720 We analytically studied the effect of ambipolar diffusion on a MRTI for an oblique magnetic field. We used a two-fluid incompressible formalism. To study an instability triggered by radiation pressure rather than gravity, we chose to differentiate the acceleration imparted on each species. This approach leads to the following conclusions:

726 1. In the presence of acceleration on both species: in this case, the dispersion relation yields two modes as solutions: a modified HD mode and a modified MHD mode the behaviours of the two fluids will be linear combinations of these two modes. These two modes are defined such that in the asymptotic case, where the coupling is zero, they correspond to the HD and MHD growth rates, respectively.

732 2. When gravity acts only on the charged species, the system presents a single solution for the growth rate of the instability. This solution transitions from an uncoupled solution to a fully coupled one, which is simply the uncoupled solution attenuated by a global factor. Between these two asymptotic solutions, the signature of the ambipolar diffusion is a range in length within which the growth rate changes slow.

740 Applying the results of the dispersion relation to the study of the interstellar medium, we showed that the magnetic field tends to enhance the anisotropy in the structuring of the RTI, especially when the field is tangent to the interface. We also demonstrated that this anisotropy becomes stronger with coupling, but only for intermediate coupling values.

747 We developed a representative framework for the interstellar medium by calculating an effective gravity that could be provided by the radiation pressure of a star located behind an interstellar cloud on our line of sight. This calculation revealed that the cutoff length below which the instability reaches its asymptote for all wavelengths is very large. Consequently, the anisotropy developed by the instability is then very strong but independent of the scale, as it saturates below a certain length scale. Under these parameters, the RTI develops over a timescale significantly longer than the characteristic time of complete coupling between neutrals and ions. Therefore, there is no clear signature of ambipolar diffusion in the growth of the RTI for the estimated values. However, it is conceivable that stars with stronger radiation or another force acting as a constant anisotropic acceleration could allow for the development of the RTI within the framework we describe in this article.

763 This work is a linear analysis of the RTI that does not capture the full complexity of the actual behaviour behind the instability. In effect, a nonlinear study could be of great interest, as partial ionisation may have a much more significant impact when nonlinear effects come into play (Soler & Ballester 2022). This remark will be the basis for a forthcoming numerical study.

770 Acknowledgements. The authors are grateful to O. Berné, F. Boulanger, K. Ferrière for fruitful discussions. We also thank A. Cogez for his contribution during his internship. This work was supported by the Thematic Action “Programme National Physique Stellaire” (PNPS) of INSU Programme National “Astro”, with contributions from CNRS Physique & CNRS Chimie, CEA, and CNES.

## References

775

776 Beattie, J. R., Federrath, C., Klessen, R. S., Cielo, S., & Bhattacharjee, A. 2024, arXiv e-prints, arXiv:2405.16626

777 Brahimi, L., Marcowith, A., & Ptuskin, V. S. 2020, Astronomy and ; Astrophysics, 633, A72

778 Bucciantini, N., Amato, E., Bandiera, R., Blondin, J. M., & Del Zanna, L. 2004, A&A, 423, 253

780 Chandrasekhar, S. 1961, Hydrodynamic and hydromagnetic stability

781 Díaz, A. J., Khomenko, E., & Collados, M. 2014, A&A, 564, A97

782 Díaz, A. J., Soler, R., & Ballester, J. L. 2012, ApJ, 754, 41

783 Gibson, S. J. & Nordsieck, K. H. 2003, The Astrophysical Journal, 589, 362

784 Goldreich, P. & Sridhar, S. 1995, ApJ, 438, 763

785 Guillet, V., Pineau Des Forêts, G., & Jones, A. P. 2007, A&A, 476, 263

786 Gull, S. F. 1975, MNRAS, 171, 263

787 Heiles, C. & Crutcher, R. 2005, in Lecture Notes in Physics, Berlin Springer Verlag, Vol. 664, Cosmic Magnetic Fields, ed. R. Wielebinski & R. Beck, 137

790 Heiles, C. & Troland, T. H. 2005, ApJ, 624, 773

791 Hensley, B. S. & Draine, B. T. 2023, ApJ, 948, 55

792 Hester, J. J., Stone, J. M., Scowen, P. A., et al. 1996, ApJ, 456, 225

793 Hillier, A. S. 2016, MNRAS, 462, 2256

794 Jacquet, E. & Krumholz, M. R. 2011, The Astrophysical Journal, 730, 116

795 Jean, P., Gillard, W., Marcowith, A., & Ferrière, K. 2009, Astronomy andamp; Astrophysics, 508, 1099–1116

796 Jun, B.-I., Jones, T. W., & Norman, M. L. 1996, ApJ, 468, L59

797 Jun, B.-I., Norman, M. L., & Stone, J. M. 1995, ApJ, 453, 332

800 Khomenko, E., Díaz, A., de Vicente, A., Collados, M., & Luna, M. 2014, A and A, 565, A45

801 Krolik, J. H. 1977, The Physics of Fluids, 20, 364

802 Mouschovias, T. C., Shu, F. H., & Woodward, P. R. 1974, A&A, 33, 73

803 Osterbrock, D. E. & Ferland, G. J. 2006, Astrophysics of gaseous nebulae and active galactic nuclei

804 Popescu Braileanu, B., Lukin, V. S., Khomenko, E., & de Vicente, Á. 2021a, A and A, 646, A93

805 Popescu Braileanu, B., Lukin, V. S., Khomenko, E., & de Vicente, Á. 2021b, A and A, 650, A181

810 Porth, O., Komissarov, S. S., & Keppens, R. 2014, MNRAS, 443, 547

811 Priest, E. R. 1987, Solar magneto-hydrodynamics.

812 Rosen, A. & Bregman, J. N. 1995, ApJ, 440, 634

813 Saury, E., Miville-Deschênes, M.-A., Hennebelle, P., Audit, E., & Schmidt, W. 2014, A&A, 567, A16

814 Soler, R. & Ballester, J. L. 2022, Frontiers in Astronomy and Space Sciences, 9, 789083

815 Soler, R. & Ballester, J. L. 2022, Frontiers in Astronomy and Space Sciences, 9

816 Vickers, E., Ballai, I., & Erdélyi, R. 2020, A and A, 634, A96

817 Weingartner, J. C. & Draine, B. T. 2001a, ApJ, 553, 581

820 Weingartner, J. C. & Draine, B. T. 2001b, ApJS, 134, 263

821 Yan, H., Lazarian, A., & Draine, B. T. 2004, ApJ, 616, 895

822 Zaqarashvili, T. V., Khodachenko, M. L., & Rucker, H. O. 2011, A and A, 534, A93