# Linear analysis of the Kelvin-Helmholtz instability in relativistic
magnetized symmetric flows

Anthony Chow1, Michael E. Rowan2, Lorenzo Sironi1, Jordy Davelaar3,1,
Gianluigi Bodo4, Ramesh Narayan5
1Department of Astronomy and Columbia Astrophysics Laboratory, Columbia
University, New York, NY 10027, USA
2Advanced Micro Devices, Inc., Santa Clara, CA, USA
3Center for Computational Astrophysics, Flatiron Institute, 162 Fifth Avenue,
New York, NY 10010, USA
4INAF Osservatorio Astrofisico di Torino, Via Osservatorio, 20, 10025 Pino
Torinese TO, Italy
5Harvard-Smithsonian Center for Astrophysics, 60 Garden Street, Cambridge, MA
02138, USA E-mail<EMAIL_ADDRESS>
(Accepted XXX. Received YYY; in original form ZZZ)

###### Abstract

We study the linear stability of a planar interface separating two fluids in
relative motion, focusing on the symmetric configuration where the two fluids
have the same properties (density, temperature, magnetic field strength, and
direction). We consider the most general case with arbitrary sound speed
$c_{\rm s}$, Alfvén speed $v_{\rm A}$, and magnetic field orientation. For the
instability associated with the fast mode, we find that the lower bound of
unstable shear velocities is set by the requirement that the projection of the
velocity onto the fluid-frame wavevector is larger than the projection of the
Alfvén speed onto the same direction, i.e., shear should overcome the effect
of magnetic tension. In the frame where the two fluids move in opposite
directions with equal speed $v$, the upper bound of unstable velocities
corresponds to an effective relativistic Mach number $M_{\rm re}\equiv
v/v_{\rm{f}\perp}\sqrt{(1-v_{\rm{f}\perp}^{2})/(1-v^{2})}\cos\theta=\sqrt{2}$,
where $v_{\rm{f}\perp}=[v_{\rm{A}}^{2}+c_{\rm s}^{2}(1-v_{\rm{A}}^{2})]^{1/2}$
is the fast speed assuming a magnetic field perpendicular to the wavevector
(here, all velocities are in units of the speed of light), and $\theta$ is the
laboratory-frame angle between the flow velocity and the wavevector projection
onto the shear interface. Our results have implications for shear flows in the
magnetospheres of neutron stars and black holes — both for single objects and
for merging binaries — where the Alfvén speed may approach the speed of light.

###### keywords:

instability – relativistic processes – (magnetohydrodynamics) MHD – plasmas –
waves – methods: analytical

††pubyear: 2015††pagerange: Linear analysis of the Kelvin-Helmholtz
instability in relativistic magnetized symmetric flows–B

##  1 Introduction

The Kelvin-Helmholtz instability (KHI) (Von Helmholtz & Monats, 1868; Lord
Kelvin, 1871) — growing in between two fluids in relative motion — is one of
the most common instabilities in space and astrophysical flows. Since the
pioneering works of Chandrasekhar (1961), the linear theory of the KHI has
been investigated under a variety of conditions (Turland & Scheuer, 1976;
Blandford & Pringle, 1976; Ferrari et al., 1980; Pu & Kivelson, 1983; Kivelson
& Zu-Yin, 1984; Bodo et al., 2004; Osmanov et al., 2008; Blumen et al., 1975;
Ferrari et al., 1978; Sharma & Chhajlani, 1998; Prajapati & Chhajlani, 2010;
Sobacchi & Lyubarsky, 2018; Berlok & Pfrommer, 2019; Rowan, 2019; Hamlin &
Newman, 2013; Bodo et al., 2013, 2016, 2019; Pimentel & Lora-Clavijo, 2019),
depending on whether the two fluids have comparable or different properties
(respectively, “symmetric” or “asymmetric” configuration), whether the flow is
incompressible or compressible, and whether or not the fluids are magnetized.
Most earlier works assumed that the relative velocity between the two fluids
is non-relativistic.

In Chow et al. (2023), we studied the linear properties of the KHI for the
asymmetric conditions expected at the boundaries of relativistic astrophysical
jets. We assumed a magnetically-dominated jet and a gas-pressure-dominated
wind, derived the most general form of the dispersion relation, and provided
an analytical approximation of its solution for wind sound speeds much smaller
than the jet Alfvén speed, as appropriate for realistic astrophysical jets.

Although asymmetric conditions are more common in astrophysical systems, it is
still interesting to consider the (potentially simpler) symmetric case in
which the two fluids have the same properties (sound speed, magnetic field
strength, and orientation). Most previous works on the KHI in the relativistic
regime focused on the symmetric case. For relativistic hydrodynamic systems
(i.e., with no magnetic fields), Bodo et al. (2004) obtained an analytical
expression for the instability growth rate and derived an absolute upper bound
on the range of unstable shear velocities. For symmetric relativistic
magnetohydrodynamical flows, Osmanov et al. (2008) considered the case of
magnetic field aligned with the direction of relative motion.

In this paper, we study the most general symmetric configuration with
arbitrary sound speed, Alfvén speed, and magnetic field orientation. We relax
the restriction of a flow-aligned magnetic field adopted by Osmanov et al.
(2008) and consider the most general scenario where both the unstable
wavevector and the magnetic field can have arbitrary directions. We derive the
most general form of the dispersion relation and present several analytical
results in the astrophysically-relevant regime where magnetic pressure
dominates over gas pressure.

The rest of the paper is organized as follows: section 2 describes the setup
of our study; section 3 derives the dispersion relation for the most general
case, having arbitrary Alfvén speed and sound speed; section 4 anticipates the
range of unstable velocities, and provides a physical interpretation for the
lower bound of the unstable region; section 5 focuses on magnetically-
dominated flows (i.e., with magnetic pressure much larger than gas pressure)
and derives analytical results for the growth rate and the range of unstable
velocities; section 6 discusses results obtained for the general case in which
the gas pressure cannot be neglected. We conclude in section 7, where we also
discuss some astrophysical implications.

##  2 Setup

In this work, we consider a magnetized planar vortex sheet in the $x$–$z$
plane at $y=0$, as shown in Figure 1. The magnetic field in fluid I ($y>0$)
and II ($y<0$) is given by $\mathbf{B}_{0}=(B_{{\rm 0x}},0,B_{{\rm 0z}})$,
where the $x$ component of the magnetic field, $B_{{\rm 0x}}$, is along the
flow direction and the $z$ component, $B_{{\rm 0z}}$, is perpendicular to the
flow direction. The magnetic field $\mathbf{B}_{0}$ is defined in the fluid
rest frame. We solve the system in the laboratory frame where fluid I moves
with velocity $\mathbf{v}=+v\hat{x}$ and fluid II moves with velocity
$\mathbf{v}=-v\hat{x}$. We use subscripts ‘+’ and ‘-’ to denote fluid I or II.
We adopt Gaussian units such that $c=4\pi=1$ and define all velocities (i.e.,
the shear velocity, the sound speed and the Alfvén speed) in units of the
speed of light $c$.

Figure 1: Schematic diagram of the interface separating two fluids in relative
motion. The interface (grey color) is located in the $x-z$ plane. Above and
below the interface are fluids I and II with the same Alfvén speed
$v_{\rm{A}}$ and sound speed $c_{\rm s}$. $\mathbf{q}_{\parallel}$ is the
projection of wavevector $\mathbf{q}$ onto the interface. We consider the
system in the laboratory frame where fluids I and II have shear velocities
$+v\hat{x}$ and $-v\hat{x}$, respectively. $\mathbf{B}_{0}$ is the magnetic
field vector in the fluid rest frame. $\theta$ is the angle between
$\mathbf{q}_{\parallel}$ and $\hat{x}$ (in the laboratory frame) while
$\Omega$ is the angle between $\mathbf{B}_{0}$ and $\hat{x}$ (in the fluid
frame).

We describe the flow with the equations of ideal relativistic
magnetohydrodynamics (RMHD) (e.g., Anile, 1990; Mignone et al., 2018; Rowan,
2019):

$\displaystyle\frac{\partial(\rho\gamma)}{\partial
t}+\nabla\cdot(\rho\gamma\mathbf{v})=0$ (1a)
$\displaystyle\frac{\partial}{\partial
t}(w\gamma^{2}\mathbf{v})+\nabla\cdot(w\gamma^{2}\mathbf{vv})+\nabla
p=\rho_{\rm e}\mathbf{E}+\mathbf{J}\times\mathbf{B}$ (1b)
$\displaystyle\frac{\partial\mathbf{B}}{\partial t}+\nabla\times\mathbf{E}=0$
(1c) $\displaystyle\frac{\partial\mathbf{E}}{\partial
t}-\nabla\times\mathbf{B}=-\mathbf{J}$ (1d)
$\displaystyle\frac{\partial}{\partial
t}(w\gamma^{2}-p)+\nabla\cdot(w\gamma^{2}\mathbf{v})=\mathbf{J}\cdot\mathbf{E}$
(1e)

supplemented with the divergence constraints

$\displaystyle\nabla\cdot\mathbf{E}=\rho_{\rm
e},\quad\nabla\cdot\mathbf{B}=0~{}.$ (2)

Here, $\rho$ is the rest mass density, $\rho_{\rm e}$ the charge density,
$\mathbf{J}$ the current density vector, $\mathbf{v}$ the fluid velocity
vector, $\gamma$ the Lorentz factor, $\mathbf{B}$ the magnetic field vector,
$\mathbf{E}$ the electric field vector, $w$ the enthalpy density and $p$ the
gas pressure. For an ideal gas with constant adiabatic index $\Gamma$, the
enthalpy can be written as $w=\rho+\Gamma p/(\Gamma-1)$.

The sound speed (see e.g. Mignone et al., 2018) and the Alfvén speed, both
defined in the fluid rest frame, can be written as

$\displaystyle c_{s}=\sqrt{\frac{w_{0}-\rho_{0}(\partial
w_{0}/\partial\rho_{0})}{(\partial w_{0}/\partial
p_{0})-1}\frac{1}{w_{0}}}=\sqrt{\Gamma\frac{p_{0}}{w_{0}}}~{}$ (3)

and

$\displaystyle{\bf{v}}_{\rm A}=\frac{{\bf B}_{0}}{\sqrt{w_{0}+B_{{\rm
0x}}^{2}+B_{{\rm 0z}}^{2}}}~{},$ (4)

where “0” subscripts indicate variables in the equilibrium state. We define
$v_{\rm A}$ as the magnitude of the Alfvén speed.

##  3 General Dispersion relation

The dispersion relation of the surface wave at the interface can be obtained
from the dispersion relation of the body waves of each fluid (Alfvén, slow,
fast modes) together with the pressure balance and displacement matching
across the interface (e.g., Bodo et al., 2004). The dispersion relation of the
body wave can be found by linearization of Equation 1 in the fluid rest frame,
by substituting the following perturbed variables

$\displaystyle{\rho}$ $\displaystyle\approx{\rho_{0}}+{\rho_{1}}~{},$ (5a)
$\displaystyle{p}$ $\displaystyle\approx{p}_{0}+{p}_{1}~{},$ (5b)
$\displaystyle{\mathbf{v}}$
$\displaystyle\approx\mathbf{0}+{\mathbf{v}}_{1}~{},$ (5c)
$\displaystyle{\mathbf{E}}$
$\displaystyle\approx\mathbf{0}+{\mathbf{E}}_{1}~{},$ (5d)
$\displaystyle\mathbf{B}$
$\displaystyle\approx\mathbf{B}_{0}+\mathbf{B}_{1}~{},$ (5e)
$\displaystyle{\rho}_{\rm e}$ $\displaystyle\approx 0+{\rho}_{\rm e1}~{},$
(5f)

where “$1$” subscripts indicate first-order perturbed variables, all defined
in the fluid rest frame. The perturbed electric field in the fluid is given by
${\mathbf{E}}_{1}=-{\mathbf{v}}_{1}\times\mathbf{B}_{0}$ in the ideal MHD
limit, and the current density is eliminated from Equation 1 by using Equation
1d. We consider perturbed variables ${X}_{1}$ to be plane waves, of the form
${X}_{1}\propto e^{i(\tilde{\mathbf{q}}\cdot\mathbf{x}-\tilde{\omega}t)}$
where $\tilde{\mathbf{q}}=(\tilde{k},\tilde{l},\tilde{m})$ is the complex
wavevector and $\tilde{\omega}$ is the complex angular frequency, both defined
in the rest frame of each fluid (we use overtilde to denote rest-frame
frequencies and wavevectors). By solving the linearized system of RMHD
equations, the dispersion relation of the body waves can then be written as

$\displaystyle\tilde{\omega}_{\pm}(\tilde{\omega}_{\pm}^{2}$
$\displaystyle-\tilde{q}_{\pm}^{2}v_{\rm{A}}^{2}\cos^{2}\tilde{\varphi}_{\pm})[\tilde{\omega}_{\pm}^{4}-(v_{\rm{A}}^{2}+c_{\rm
s}^{2}(1-v_{\rm{A}}^{2}))\tilde{q}_{\pm}^{2}\tilde{\omega}_{\pm}^{2}$
$\displaystyle+(\tilde{q}_{\pm}^{2}-\tilde{\omega}_{\pm}^{2})c_{\rm
s}^{2}v_{\rm{A}}^{2}\tilde{q}_{\pm}^{2}\cos^{2}\tilde{\varphi}_{\pm}]=0~{},$
(6)

where

$\displaystyle\cos\tilde{\varphi}_{\pm}=\frac{\tilde{k}_{\pm}B_{{\rm
0x}}+\tilde{m}_{\pm}B_{{\rm
0z}}}{\sqrt{\tilde{k}_{\pm}^{2}+\tilde{l}_{\pm}^{2}+\tilde{m}_{\pm}^{2}}\sqrt{B_{{\rm
0x}}^{2}+B_{{\rm 0z}}^{2}}}~{},$ (7)

is the cosine of the angle between the magnetic field and the fluid-frame
wavevector $\tilde{\mathbf{q}}$. The second term of section 3 (in round
brackets) corresponds to Alfvén waves, which lead to stable solutions. The
third term in section 3 (in square brackets) has solutions

$\displaystyle\frac{\tilde{\omega}_{\pm}^{2}}{\tilde{q}_{\pm}^{2}}$
$\displaystyle=\frac{1}{2}\left(v_{\rm{A}}^{2}+c_{\rm
s}^{2}(1-v_{\rm{A}}^{2})+v_{\rm{A}}^{2}c_{\rm
s}^{2}\cos^{2}\tilde{\varphi}_{\pm}\right.$
$\displaystyle\quad\left.\pm\sqrt{[v_{\rm{A}}^{2}+c_{\rm
s}^{2}(1-v_{\rm{A}}^{2}\sin^{2}\tilde{\varphi}_{\pm})]^{2}-4v_{\rm{A}}^{2}c_{\rm
s}^{2}\cos^{2}\tilde{\varphi}_{\pm}}\right)~{},$ (8)

which correspond to fast and slow magnetosonic modes (plus and minus sign,
respectively), which are potentially unstable (Osmanov et al., 2008). By
setting $\tilde{\varphi}_{\pm}=\pi/2$ in the fast mode, we can define a frame-
independent quantity

$\displaystyle v_{\rm{f}\perp}=\sqrt{v_{\rm{A}}^{2}+c_{\rm
s}^{2}(1-v_{\rm{A}}^{2})}~{},$ (9)

which is the fast magnetosonic speed in the case that the fluid-frame
wavevector and the magnetic field are orthogonal. Hereafter, we will use
$v_{\rm{f}\perp}$ as defined in Equation 9 to parameterize our results, even
in the general case when the wavevector and the magnetic field are not
orthogonal.

Notice that the transverse components of the wavenumber, $l$ and $m$, are
Lorentz invariant. By applying the following transformations

$\displaystyle\tilde{\omega}_{\pm}=\gamma(\omega\mp
kv),\quad\tilde{k}_{\pm}=\gamma(k\mp\omega
v),\quad\tilde{l}_{\pm}=l_{\pm},\quad\tilde{m}_{\pm}=m,$ (10)

to section 3, we can construct the Lorentz invariant ratio
$l_{+}^{2}/l_{-}^{2}$:

$\displaystyle\frac{l_{+}^{2}}{l_{-}^{2}}=f(B_{{\rm 0x}}/B_{{\rm
0z}},v,c_{s},v_{\rm{A}},k,m,\omega)~{},$ (11)

where the function $f$ can be written explicitly as a function of the
parameters indicated in parenthesis. An independent way of obtaining
$l_{+}/l_{-}$ is to solve the linearized RMHD equations, Eqs. (1), together
with two additional constraints. First, we enforce the first-order pressure
balance between the two sides of the interface:

$\displaystyle p_{1+}\\!+\\!B_{{\rm 0x}}B_{\rm 1x,+}\\!+\\!B_{{\rm 0z}}B_{\rm
1z,+}\\!=\\!p_{1-}\\!+\\!B_{{\rm 0x}}B_{\rm 1x,-}\\!+\\!B_{{\rm 0z}}B_{\rm
1z,-}~{}.$ (12)

Second, the fluid displacements in the $y$ direction must match at the
interface. Since the Lagrangian derivative of the displacement is equal to the
transverse velocity $v_{1\rm y}$ of the fluid, matching the displacements is
equivalent to

$\displaystyle\frac{v_{\rm 1y,-}}{\gamma(\omega+kv)}=\frac{v_{\rm
1y,+}}{\gamma(\omega-kv)}~{}.$ (13)

This results in an independent solution for the ratio $l_{+}/l_{-}$:

$\displaystyle\frac{l_{+}}{l_{-}}=g(B_{{\rm 0x}}/B_{{\rm 0z}},v,c_{\rm
s},v_{\rm{A}},k,m,\omega)~{}.$ (14)

where the function $g$ can be written as a function of the parameters
indicated in parenthesis.

We can finally equate Equation 11 and the square of Equation 14 to obtain the
dispersion relation in an implicit form:

$\displaystyle f=g^{2}~{}.$ (15)

We can then rearrange Equation 15 by introducing the dimensionless ratios

$\displaystyle\phi=\frac{\omega}{v_{\rm{f}\perp}\sqrt{k^{2}+m^{2}}},\quad
M=\frac{v}{v_{\rm{f}\perp}},$ (16)

where $M$ is the flow fast Mach number, i.e., the flow speed in units of the
fast magnetosonic speed defined in Equation 9. With this parameterization, the
dispersion relation in Equation 15 takes the form

$\displaystyle\left(\frac{\phi}{M}\right)\left[c_{0}+c_{1}\left(\frac{\phi}{M}\right)^{2}+c_{2}\left(\frac{\phi}{M}\right)^{4}+c_{3}\left(\frac{\phi}{M}\right)^{6}+c_{4}\left(\frac{\phi}{M}\right)^{8}\right]=0,$
(17)

where $c_{0},c_{1},c_{2},c_{3}$ and $c_{4}$ are constants with respect to
$\phi$. Note that Im$(\phi)>0$ implies that the amplitude of the perturbation
grows exponentially, i.e., the system becomes unstable. Since the factor
$\phi/M$ in front of Equation 17 corresponds to a stable mode, Equation 17 is
effectively a quartic equation in $(\phi/M)^{2}$ and can be solved
analytically (Weisstein, 2020), giving a total of eight (generally, complex)
roots. However, not all of them may be acceptable. First, since we have
introduced spurious roots when squaring it, i.e., not all solutions will
satisfy Equation 14. Also, only outgoing waves should be retained by the
Sommerfeld radiation condition (Sommerfeld, 1912). This requires accepting
only those roots that satisfy the conditions:

$\displaystyle\text{Im}(l_{+})>0$ (18) $\displaystyle\text{Im}(l_{-})<0.$ (19)

In the rest of the paper, we present solutions to the dispersion relation. We
first focus on the special case of magnetically-dominated flows, i.e., $c_{\rm
s}\ll v_{\rm A}$, and then we address the most general case. We employ the
following parameters:

$\displaystyle M_{\rm
r}=M\frac{\sqrt{1-v_{\rm{f}\perp}^{2}}}{\sqrt{1-v^{2}}},\quad\cos\Omega=\frac{B_{{\rm
0x}}}{\sqrt{B_{{\rm 0x}}^{2}+B_{{\rm
0z}}^{2}}},\quad\cos\theta=\frac{k}{\sqrt{k^{2}+m^{2}}},$ (20)

where $M_{\rm r}$ is a “relativistic” Mach number defined with the spatial
part of the 4-velocity, as in Osmanov et al. (2008), $\Omega$ is the angle
between the $+\hat{x}$ direction and the magnetic field (in the fluid rest
frame), and $\theta$ is the angle (measured in the laboratory frame) between
$+\hat{x}$ and the wavevector projection onto the $y=0$ interface. It is also
convenient to absorb the dependence on $\cos\theta$ by defining an effective
Mach number $M_{\rm re}$ and an effective flow velocity $v_{\rm e}$ similarly
to the purely hydrodynamic case discussed by Bodo et al. (2004):

$\displaystyle M_{\rm re}=M_{\rm
r}\cos\theta=\frac{v}{v_{\rm{f}\perp}}{}\frac{\sqrt{1-v_{\rm{f}\perp}^{2}}}{\sqrt{1-v^{2}}}\cos\theta,\quad
v_{\rm e}=v\cos\theta.$ (21)

In the following, we usually write our solution for $\phi$ in terms of $M_{\rm
re}$ and $v_{\rm e}$ (clearly, when $\cos\theta=1$, we have $M_{\rm
re}=M_{\rm{r}}$ and $v_{\rm e}=v$). For $c_{\rm s}=0$, in some cases, it will
be convenient to use $v_{\rm{A}}$ instead of $v_{\rm e}$. Assuming $c_{\rm
s}=0$, the two are related by the following identity

$v_{\rm e}=\frac{v_{\rm{A}}M_{\rm re}\cos\theta}{\sqrt{\cos^{2}\theta-
v_{\rm{A}}^{2}\cos^{2}\theta+M_{\rm re}^{2}v_{\rm{A}}^{2}}}~{}.$ (22)

##  4 Boundaries of the unstable region

We anticipate here some general remarks on the range of unstable velocities
$v$ (or equivalently, of $M_{\rm re}$), which will be confirmed below in
several specific cases. For the upper bound, we expect the system to be
unconditionally stable if $M_{\rm re}\geq\sqrt{2}$. This bound has already
been discussed in the purely hydrodynamic case $v_{\rm{A}}=0$ by Bodo et al.
(2004), where $M_{\rm re}$ was defined with the sound speed. Here, we find
that the same bound holds, as long as $M_{\rm re}$ is defined with the fast
magnetosonic speed as in Equation 21.

A physically-insightful lower bound can be found in the case of instability
associated with the fast mode. As shown below, this case (as opposed to the
case of instability associated with the slow mode) is the most important one
since it leads to faster growth rates, and it is the only one that persists
for magnetically-dominated flows with $c_{\rm s}=0$. As we show in Appendix B,
unstable modes associated with the fast magnetosonic wave have vanishing phase
speeds (so, zero real part of $\omega$). The boundaries of the unstable region
will also have, by construction, a vanishing growth rate (so, zero imaginary
part of $\omega$). This implies that $\tilde{k}_{\pm}=\gamma k$ (see Equation
10) so that, in the fluid comoving frame, the angle $\tilde{\theta}$ between
the “projected wavevector” $\tilde{\bf q}_{\parallel}$ (i.e., the projection
of $\tilde{\bf q}$ onto the $x-z$ interface) and the $\hat{x}$ direction of
the shear velocity satisfies

$\displaystyle\gamma\tan\tilde{\theta}=\tan{\theta}~{}.$ (23)

We expect the system to be unstable only if the projection of the shear
velocity onto the direction of $\tilde{\bf{q}}_{\parallel}$ is larger than the
projection of the Alfvén speed onto the same direction. Equivalently,

$v\cos\tilde{\theta}>v_{\rm A}\cos(\Omega-\tilde{\theta})~{}.$ (24)

In other words, instability happens when the shear can overcome the effect of
magnetic tension (for particular cases of this general criterion, see Pu &
Kivelson (1983); Roychoudhury & Lovelace (1986)). Using Equation 23, this can
be written as a function of lab-frame quantities as

$\frac{v}{v_{\rm{A}}}>\cos\Omega+\frac{\sin\Omega\,\tan\theta}{\gamma}~{}.$
(25)

We demonstrate below that this successfully describes the lower bound of
unstable $v$ in the case of fast-mode instability for all $c_{\rm s}$,
$v_{\rm{A}}$, and field orientations.

##  5 Magnetically-dominated flows

Although the dispersion relation in Equation 17 for the most general symmetric
shear flows is analytically solvable, the solution is too complex for further
analysis of the instability region. However, in certain regimes of our
parameter space, the solution simplifies and is physically interpretable. This
section considers the special case where the flow is magnetically dominated
($c_{\rm s}\ll v_{\rm{A}}$). We then assume the cold plasma limit ($c_{\rm
s}=0$) and the solution is analytically tractable and allows for physical
insight.

For magnetically-dominated flows, we assume $c_{\rm s}=0$ so that
$v_{\rm{f}\perp}=v_{\rm{A}}$ and the dispersion relation (Equation 17) is
reduced to a quadratic equation in $(\phi/M)^{2}$ with solutions:

$\displaystyle\frac{\phi^{2}}{M^{2}}$
$\displaystyle=\frac{\mu_{1}\pm\sqrt{\mu_{2}}}{\mu_{3}}$ (26)

where $\mu_{0},\mu_{1}$ and $\mu_{2}$ are given in Appendix A for
completeness.

We now consider a few special cases. First, we consider the case where the
projected wavevector is parallel to the shear flow ($\cos\theta=1$) for
arbitrary $v_{\rm{A}}$ and $\cos\Omega$. Then, we retain the dependence on
$\cos\theta$ and $v_{\rm{A}}$, but consider the specific cases of in-plane or
out-of-plane fields ($B_{\rm 0z}=0$ or $B_{\rm 0x}=0$, respectively;
equivalently, $\cos\Omega=1$ or $\cos\Omega=0$, respectively). Finally, we
retain the angular dependence on $\cos\theta$ and $\cos\Omega$, but we
consider the extreme limit $v_{\rm{A}}\rightarrow 1$. These cases are
summarized in Table 1.

Table 1: Magnetically-dominated cases. $v_{\rm{A}}$ | $\cos\theta$ | $\cos\Omega$ | Growth Rate | Unstable Range
---|---|---|---|---
Arbitrary | 1 | Arbitrary | Equation 27 | Equation 28
Arbitrary | Arbitrary | 1 | Equation 30 | Equation 31
Arbitrary | Arbitrary | 0 | Equation 33 | Equation 34
1 | Arbitrary | Arbitrary | Equation 36 | Equation 37

###  5.1 Magnetically-dominated flows with projected wavevector along the
shear velocity

We now consider the case where the projected wavevector is parallel to the
shear flow ($\cos\theta=1$, in which case $M_{\rm re}=M_{\rm{r}}$ and $v_{\rm
e}=v$), where the solution of Equation 26 can be written as

$\displaystyle\frac{\phi^{2}}{M^{2}}$
$\displaystyle=\frac{\kappa_{1}\pm\sqrt{\kappa_{2}}}{\kappa_{3}},$ (27)

where

$\displaystyle\kappa_{1}$ $\displaystyle=M_{\rm
r}^{4}v_{\rm{A}}^{2}(1-v_{\rm{A}}^{2}\cos^{2}\Omega)^{2}+M_{\rm
r}^{2}(1-v_{\rm{A}}^{2})(1-v_{\rm{A}}^{2}\cos^{2}\Omega)^{2}$
$\displaystyle+(1-v_{\rm{A}}^{2})(1-v_{\rm{A}}^{2}\cos^{4}\Omega)$
$\displaystyle\kappa_{2}$
$\displaystyle=(1-v_{\rm{A}}^{2}\cos^{2}\Omega)^{2}[(1-\cos^{2}\Omega)^{2}(1-v_{\rm{A}}^{2})^{2}$
$\displaystyle+4M_{\rm r}^{2}(1-\cos^{2}\Omega
v_{\rm{A}}^{2})^{2}(1-v_{\rm{A}}^{2}+M_{\rm r}^{2}v_{\rm{A}}^{2})]$
$\displaystyle\kappa_{3}$ $\displaystyle=M_{\rm r}^{2}[1-v_{\rm{A}}^{2}+M_{\rm
r}^{2}v_{\rm{A}}^{2}(1-v_{\rm{A}}^{2}\cos^{2}\Omega)][1+(1-2\cos^{2}\Omega)v_{\rm{A}}^{2}$
$\displaystyle+M_{\rm
r}^{2}v_{\rm{A}}^{2}(1-v_{\rm{A}}^{2}\cos^{2}\Omega)]/(1-v_{\rm{A}}^{2}+M_{\rm
r}^{2}v_{\rm{A}}^{2}).$

Since $\kappa_{2}$ and $\kappa_{3}$ are always positive, the condition for
instability can be found from the inequality $\kappa_{1}^{2}<\kappa_{2}$. In
this case, the solution is purely imaginary, and the instability range in
$M_{\rm r}$ reads

$\displaystyle\frac{\cos\Omega\sqrt{1-v_{\rm{A}}^{2}}}{\sqrt{1-v_{\rm{A}}^{2}\cos^{2}\Omega}}<M_{\rm
r}<\frac{\sqrt{2-\cos^{2}\Omega-
v_{\rm{A}}^{2}\cos^{2}\Omega}}{\sqrt{1-v_{\rm{A}}^{2}\cos^{2}\Omega}}~{}.$
(28)

Equation 28 completely characterizes the boundary of the instability region
for cold flows when $\cos\theta=1$. Note that the left inequality in Equation
28 can be cast as $v_{\rm{A}}\cos\Omega<v$, which implies that the interface
is unstable only when the shear speed is greater than the projection of the
Alfvén speed onto the flow direction (in this case, the latter coincides with
the direction of the projected wavevector), in agreement with Equation 25 for
${\theta}=0$.

By setting $\cos\Omega=0$ in Equation 28 (i.e., magnetic field perpendicular
to the shear flow), the instability region attains its maximum range.
Similarly to the purely hydrodynamic case (but with the fast speed now
replacing the sound speed in $M_{\rm re}$), we find

$\displaystyle 0<M_{\rm r}<\sqrt{2}~{}.$ (29)

###  5.2 Magnetically-dominated flows with in-plane magnetic field

When the magnetic field is parallel to the flow velocity, i.e.,
$\cos\Omega=1$, Equation 26 can be rewritten as:

$\displaystyle\frac{\phi^{2}}{M^{2}}$
$\displaystyle=\frac{\nu_{1}\pm\sqrt{\nu_{2}}}{\nu_{3}}\cos^{2}\theta$ (30)

where

$\displaystyle\nu_{1}=M_{\rm re}^{4}+v_{\rm e}^{2}+M_{\rm re}^{2}(1-v_{\rm
e}^{2}),$ $\displaystyle\nu_{2}=\sin^{4}\theta(M_{\rm re}^{2}+v_{\rm
e}^{2})^{2}+4M_{\rm re}^{4}[v_{\rm e}^{2}+M_{\rm re}^{2}(1-v_{\rm e}^{2})],$
$\displaystyle\nu_{3}=(M_{\rm re}^{2}+v_{\rm e}^{2})^{2}.$

Since $\nu_{2}$ and $\nu_{3}$ are always positive, the condition for
instability can be found from the inequality $\nu_{1}^{2}<\nu_{2}$. In this
case, the solution is purely imaginary, and the instability range in $M_{\rm
re}$ is

$\displaystyle\cos\theta<M_{\rm re}<\sqrt{2-\cos^{2}\theta}~{},$ (31)

where the case $\cos\theta=0$ should be excluded since it yields a vanishing
growth rate. Equation 31 completely characterizes the boundary of the
instability region for cold flows when $\cos\Omega=1$. Note that the left
inequality in Equation 31 can be cast as $v_{\rm{A}}<v$, implying that the
interface is unstable only when the flow is super-Alfvénic. This is consistent
with the general expression in Equation 25 for the case $\Omega=0$ considered
here.

By setting $\cos\theta\rightarrow 0$ in Equation 31, the instability region
attains its maximum range

$\displaystyle 0<M_{\rm re}<\sqrt{2}~{}.$ (32)

###  5.3 Magnetically-dominated flows with out-of-plane magnetic field

We now consider the case of the magnetic field perpendicular to the shearing
flow, i.e., $\cos\Omega=0$. Since the instability boundaries depend explicitly
on $v_{\rm{A}}$, we choose to express the solution in terms of $v_{\rm{A}}$
instead of $v_{\rm e}$. By setting $\Omega=\pi/2$ in Equation 26, we obtain
the solution for $\phi$ in the form

$\displaystyle\frac{\phi^{2}}{M^{2}}$
$\displaystyle=\frac{\lambda_{1}\pm\sqrt{\lambda_{2}}}{\lambda_{3}}\cos^{2}\theta~{},$
(33)

where

$\displaystyle\lambda_{1}$ $\displaystyle=M_{\rm re}^{2}(1+M_{\rm
re}^{2})v_{\rm{A}}^{2}-\cos^{4}\theta v_{\rm{A}}^{2}(1-v_{\rm{A}}^{2})$
$\displaystyle+\cos^{2}\theta[1-v_{\rm{A}}^{4}+M_{\rm
re}^{2}(1-2v_{\rm{A}}^{2})]$ $\displaystyle\lambda_{2}$
$\displaystyle=\cos^{8}\theta(1\\!-\\!v_{\rm{A}}^{2})^{2}+M_{\rm
re}^{4}v_{\rm{A}}^{2}(4\cos^{2}\theta\\!+\\!4M_{\rm re}^{2}\sin^{2}\theta
v_{\rm{A}}^{2}+\sin^{4}\theta v_{\rm{A}}^{6})$ $\displaystyle+2\cos^{2}\theta
M_{\rm re}^{2}(1\\!-\\!v_{\rm{A}}^{2})[2M_{\rm re}^{2}\sin^{2}\theta
v_{\rm{A}}^{2}\\!+\\!\cos^{4}\theta
v_{\rm{A}}^{4}\\!+\\!\cos^{2}\theta(2\\!-\\!v_{\rm{A}}^{4})]$
$\displaystyle\lambda_{3}$ $\displaystyle=M_{\rm
re}^{4}v_{\rm{A}}^{2}\\!+\\!\cos^{2}\theta M_{\rm
re}^{2}(1\\!+\\!v_{\rm{A}}^{2})~{}.$

Since $\lambda_{2}$ and $\lambda_{3}$ are always positive, the criterion for
unstable modes can be derived from the inequality
$\lambda_{1}^{2}<\lambda_{2}$. The solution is purely imaginary, and the
instability range in $M_{\rm re}$ is

$\displaystyle\sqrt{(1-\cos^{2}\theta)(1-v_{\rm{A}}^{2})}<M_{\rm
re}<\sqrt{1+\cos^{2}\theta+v_{\rm{A}}^{2}(1-\cos^{2}\theta)}~{},$ (34)

where the case $\cos\theta=0$ should be excluded since it yields a vanishing
growth rate. Equation 34 completely characterizes the boundary of the
instability region for cold flows when $\cos\Omega=0$. Note that the left
inequality of Equation 34 can be cast as $v_{\rm{A}}\tan\theta<\gamma v$. Once
again, this is consistent with Equation 25 for the case $\Omega=\pi/2$
considered here.

By setting $\cos\theta=1$, the instability region attains its maximum range

$\displaystyle 0<M_{\rm re}<\sqrt{2}~{},$ (35)

as we had already derived in Equation 29 under similar assumptions
($\cos\theta=1$ and $\cos\Omega=0$).

###  5.4 Flows with Alfvén speed approaching the speed of light

The system remains unstable in the extreme limit $v_{\rm{A}}\rightarrow 1$. In
this case, Equation 26 reduces to

$\displaystyle\frac{\phi^{2}}{M^{2}}=\frac{\eta_{1}^{2}-\eta_{2}^{2}}{\eta_{3}}~{},$
(36)

where

$\displaystyle\eta_{1}$
$\displaystyle=\gamma\cos\theta\cos\Omega+\sin\theta\sin\Omega~{},$
$\displaystyle\eta_{2}$ $\displaystyle=\gamma v\cos\theta~{},$
$\displaystyle\eta_{3}$
$\displaystyle=\gamma^{2}v^{2}(1-v^{2}\cos^{2}\Omega)~{}.$

Since $\eta_{3}$ is always positive, the criterion for unstable modes can be
derived from the inequality $\eta_{1}^{2}<\eta_{2}^{2}$. In this case, the
solution is purely imaginary, and the instability range is

$\displaystyle\cos\Omega+\frac{\sin\Omega\tan\theta}{\gamma}<v<1~{}.$ (37)

The lower limit coincides with Equation 25 for $v_{\rm{A}}=1$.

In the special case of out-of-plane fields, i.e., $\cos\Omega=0$, the
instability growth rate is

$\displaystyle\frac{\phi^{2}}{M^{2}}=\frac{\sin^{2}\theta-v^{2}}{v^{2}}~{},$
(38)

whereas the instability range can be cast as

$\displaystyle{\sin\theta<v<1}~{}.$ (39)

##  6 General Case

In this section, we discuss the dependence of the growth rate on $c_{\rm s}$
and $v_{\rm{A}}$. The general solution for arbitrary $c_{\rm s}$ and
$v_{\rm{A}}$ is analytically available but too lengthy to be presented, so we
focus on numerical results.

###  6.1 Comparison to the hydrodynamic case

When the unstable modes propagate perpendicular to the magnetic field, we
expect magnetic tension to have no effect, and the solution should resemble
the hydrodynamic symmetric case discussed by Bodo et al. (2004), but now with
the fast magnetosonic speed instead of the sound speed in the definition of
$M_{\rm re}$. We demonstrate this by setting $\cos\Omega=0$ and $\cos\theta=1$
in the general dispersion relation (Equation 17) and the solution can be
written as

$\displaystyle\frac{\phi^{2}}{M^{2}}=\frac{1+M_{\rm
r}^{2}-v^{2}\pm\sqrt{4M_{\rm r}^{2}(1-v^{2})+(1+v^{2})^{2}}}{M_{\rm
r}^{2}+2v^{2}}$ (40)

which is essentially the same as Equation 15 in Bodo et al. (2004) for
relativistic hydrodynamic flows when the projected wavevector is parallel to
the shear flow. We conclude that even though the system is magnetized, in the
case when $\cos\Omega=0$ and $\cos\theta=1$, the instability behaves similarly
to the hydrodynamic case. Here, the magnetic field provides pressure but not
tension. We remind that the Mach numbers in Equation 40 are defined using
$v_{\rm{f}\perp}$, and so we allow for an arbitrary degree of gas versus
magnetic contributions to the pressure (in contrast to Bodo et al. (2004),
which assumed $v_{\rm{A}}=0$). In analogy to what found by Bodo et al. (2004),
the instability range is

$0<M_{\rm r}<\sqrt{2}~{}.$ (41)

Figure 2: Dependence of the instability growth rate Im($\phi$) on $M_{r}$ and
$\Omega$ when $\theta=0$ (i.e., wavevector along the shear flow), for three
choices of $v_{\rm{A}}$ and $c_{\rm s}$. Panels in the first, second, and
third columns represent $v_{\rm{A}}=0.1,0.5$ and $0.9$. Panels in the first,
second, and third rows represent $c_{\rm s}=0,1/2\sqrt{3},1/\sqrt{3}$. In all
the panels, Im($\phi$) is then normalized to its maximum value, which is
quoted in the panels themselves. The dashed cyan vertical lines represent the
common upper bound $M_{\rm re}=\sqrt{2}$ across all panels. The dashed white
lines in the top three panels represent the exact boundaries of the
instability region when $c_{\rm s}=0$ and $\cos\theta=1$, see Equation 28. The
dotted white lines in the middle and bottom rows represent the lower bound of
the fast-mode unstable region imposed by magnetic tension, as in Equation 25.
Figure 3: Dependence of the instability growth rate Im($\phi$) on $M_{\rm re}$
and $\theta$ when $\Omega=0$ (i.e., magnetic field along the flow direction),
for three choices of $v_{\rm{A}}$ and $c_{\rm s}$. The dashed cyan vertical
lines represent the common upper bound $M_{\rm re}=\sqrt{2}$ across all
panels. The dashed white lines in the top three panels represent the exact
boundaries of the instability region when $c_{\rm s}=0$ and $\cos\Omega=1$,
see Equation 31. See the caption of Figure 3 for further details. Figure 4:
Dependence of the instability growth rate Im($\phi$) on $M_{\rm re}$ and
$\theta$ when $\Omega=\pi/2$ (i.e., magnetic field perpendicular to the flow
direction), for three choices of $v_{\rm{A}}$ and $c_{\rm s}$. The dashed cyan
vertical lines represent the common upper bound $M_{\rm re}=\sqrt{2}$ across
all panels. The dashed white lines in the top three panels represent the exact
boundaries of the instability region when $c_{\rm s}=0$ and $\cos\Omega=0$,
see Equation 34. See the caption of Figure 3 for further details.

###  6.2 General case

We now consider the problem in its full generality, allowing for arbitrary
$c_{\rm s}$, $v_{\rm{A}}$, $\theta$, and $\Omega$. Figure 2 shows the
dependence of the instability growth rate Im($\phi$) on $M_{\rm re}$ and
$\cos\Omega$ when the wavevector is parallel to the shear flow
($\cos\theta=1$, in which case $M_{\rm re}=M_{r}$). Figure 3 and Figure 4 show
the dependence of the instability growth rate on $M_{\rm re}$ and $\cos\theta$
when the orientation of the magnetic field is parallel ($\cos\Omega=1$) and
perpendicular ($\cos\Omega=0$) to the shear flow, respectively. The case in
Figure 3 is the same as in Osmanov et al. (2008). In all three figures, we
consider different combinations of $v_{\rm{A}}$ and $c_{\rm s}$. The panels in
the left, middle, and right columns correspond to $v_{\rm{A}}=0.1,0.5$ and
$0.9$, respectively; the panels in the top, middle and bottom rows correspond
to $c_{\rm s}=0,1/\sqrt{12}$ and $1/\sqrt{3}$ respectively. The figures plot
the imaginary part of $\phi$ (i.e., the growth rate). We present the real part
of the solution (i.e., the phase speed) in Appendix B.

In all figures, we identify two separate instability regions in the middle-
middle, middle-right, and bottom-right panels. The narrow purple instability
stripes are associated with the slow-mode magnetosonic wave, while the large
instability regions correspond to the fast-mode wave,111More precisely, the
narrow purple stripes correspond to slow waves in one fluid and fast waves in
the other fluid; the large instability region corresponds to fast waves in
both fluids. see Osmanov et al. (2008). For $c_{\rm s}=0$ (top three panels in
all figures), only the fast-mode solution survives. The middle-left, bottom-
left, and bottom-middle panels do not have visible slow-mode instability
bands, but further analysis reveals that they exist near the left boundaries
of the main instability zones whenever $c_{\rm s}\neq 0$.

In Appendix B, we show that unstable modes associated with the slow wave have
nonzero phase speed. In contrast, unstable modes associated with the fast
magnetosonic wave have vanishing phase speed. As discussed in section 4, the
lower bound of unstable velocities — regarding the fast-mode instability —
should obey Equation 25. This is depicted as a dotted white line in the
figures (for the top row, Equation 25 corresponds to the left branch of the
dashed lines), confirming that Equation 25 successfully describes the lower
bound of unstable $M_{\rm re}$ of the fast-mode instability.

For the gas pressure dominated case with $c_{\rm s}=1/\sqrt{3}$ and
$v_{\rm{A}}=0.1$ (bottom-left panels in Figures 2-4), we expect the
orientation of the magnetic field (i.e., $\cos\Omega$) to have little effect
on the instability growth rate, as indeed confirmed by the bottom-left panel
of Figure 2. In the hydrodynamic case, Bodo et al. (2004) found that the
growth rate peaks at $\cos\theta=1$, it vanishes for $\cos\theta=0$, and that
the range of unstable Mach numbers is $0<M_{\rm re}<\sqrt{2}$ for all
$\cos\theta\neq 0$. This is in agreement with the near-hydro case displayed in
the bottom-left panels of Figure 3 and Figure 4.

For magnetically-dominated flows, we first consider the top row of Figure 3
where the plasma is cold ($c_{\rm s}=0$), and the magnetic field is along the
shear velocity ($\cos\Omega=1$). When the wavevector is perpendicular to the
flow ($\cos\theta=0$), the flow is stable, regardless of $v_{\rm{A}}$. When
the wavevector is parallel to the magnetic field ($\cos\theta=1$ in this
case), the effect of magnetic tension in suppressing the instability is the
greatest, and the flow is also stable. Since the system is stable when
$\cos\theta$ is either 0 or 1, the maximum growth rate is attained at
intermediate values of $\cos\theta$. The boundary of the instability region is
fully characterized by Equation 31, which is depicted as white dashed lines in
Figure 3. The left boundary of the unstable region corresponds to
$v_{\rm{A}}<v$, indicating that the system is unstable only when the shear can
overcome magnetic tension. The upper bound of the instability region is
determined by the right inequality in Equation 31. For $\cos\theta\rightarrow
0$, magnetic pressure dominates over tension, and the range of unstable Mach
numbers approaches the one of the hydrodynamic case (compare with the bottom-
left panel in Figure 3).

The magnetically-dominated flows with $\cos\Omega=0$, i.e., field orthogonal
to the shear velocity (top row in Figure 4) show instability profiles which
are vastly different than the case in Figure 3. At $\cos\theta=1$, the
magnetic field is perpendicular to the wavevector, parallel to the shear flow.
Since magnetic tension plays no role, both the instability growth rate and the
unstable Mach number range at $\cos\theta=1$ resemble those of the gas-
pressure dominated case (bottom left panel of Figure 4; see also our
discussion in subsection 6.1). When $\cos\theta\rightarrow 0$, the wavevector
aligns more and more closely to the magnetic field, and magnetic tension gets
stronger. Therefore, the instability growth rate gets suppressed when
$\cos\theta\rightarrow 0$, and the instability range shrinks. It is
interesting to note that at intermediate $\cos\theta$ between 0 and 1, the
range of unstable Mach numbers widens for higher $v_{\rm{A}}$. This has a
simple explanation: for fixed $\cos\theta$ and fixed $M_{\rm re}$, increasing
$v_{\rm{A}}$ implies an increase in $v$. In the fluid frame, the projected
wavevector of the fast-mode instability forms an angle $\tilde{\theta}$ with
respect to the shear velocity such that
$\gamma\tan\tilde{\theta}=\tan{\theta}$, see Equation 23. Therefore, with
increasing $v$, the fluid-frame wavevector becomes more aligned with the flow
and nearly orthogonal to the magnetic field. This implies that at fixed
$\cos\theta$ and increasing $v_{\rm{A}}$, the system can more easily overcome
the effect of magnetic tension.

For magnetically-dominated flows with arbitrary orientation of the magnetic
field and fixed $\cos\theta=1$ (top row in Figure 2), the instability profiles
at $\cos\Omega=1$ and $\cos\Omega=0$ obviously coincide with those of Figure 3
and Figure 4 for $\cos\theta=1$. The white dashed boundaries of the
instability region are described by Equation 28, where the left inequality can
be cast as $v_{\rm{A}}\cos\Omega<v$, indicating that the interface is unstable
only when the shear speed is greater than the projection of the Alfvén speed
onto the flow direction (which, in this case, coincides with the direction of
the projected wavevector).

For the general case where both $v_{\rm{A}}$ and $c_{\rm s}$ are nonzero
(middle row of Figures 2-4), the instability behavior can be understood as an
intermediate case between the gas-pressure-dominated case and the
magnetically-dominated case.

Finally, we consider the global boundary of the instability region in the
general case. Although the general dispersion relation in Equation 17 is
analytically solvable, the explicit expression of the unstable growth rate is
too complicated to derive an instability boundary in terms of $M_{\rm re}$
analytically. However, the instability criteria that we derived for
magnetically-dominated flows in section 5 allow us to argue for a global
instability boundary in the general case. In the top three panels of Figures
2-4, the white dashed lines correspond to the instability boundaries described
by Equation 28, 31 and 34 respectively. In each case, as well as in the purely
hydrodynamic case, the maximum instability range is bounded below by 0 and
above by $\sqrt{2}$, as shown in Equation 29, 32 and 35. The lines $M_{\rm
re}=\sqrt{2}$ (the cyan dashed lines in all panels) appear to provide an upper
limit to the unstable region in all cases. Therefore, we argue that the
maximum instability range in $M_{\rm re}$ for general symmetric flows is

$\displaystyle 0<M_{\rm re}<\sqrt{2}~{},$ (42)

i.e., that flows with $M_{\rm re}\geq\sqrt{2}$ are unconditionally stable.

##  7 Discussion and Conclusions

We have studied the linear stability properties of the KHI for relativistic,
symmetric, magnetized flows. We considered arbitrary sound speeds, Alfvén
speeds, and magnetic field orientations and derived the most general form of
the dispersion relation.

Our results show that, for $c_{\rm s}\neq 0$, there are two distinct unstable
regions, corresponding to the slow and fast modes of the magnetosonic wave,
consistent with the findings in Osmanov et al. (2008). For the instability
associated with the fast mode, which leads to greater growth rates, we
obtained a lower bound on the range of unstable shear velocities, see Equation
25. This lower bound reveals a necessary condition for the KHI growth, namely
that the projection of the shear velocity $v$ onto the direction of
$\tilde{\bf q}_{\parallel}$ (which, itself, is the fluid-frame wavevector
projected onto the interface) should be larger than the projection of the
Alfvén speed onto the same direction. In other words, for the fast mode to
become unstable, the free energy of the shear must overcome magnetic tension.
Concerning the upper bound of unstable shear velocities, we found that the
system is unconditionally stable if $M_{\rm re}\geq\sqrt{2}$, regardless of
$c_{\rm s}$, $v_{\rm{A}}$ or magnetic field direction. Here, $M_{\rm re}$ is
an effective Mach number defined with the fast magnetosonic speed, see
Equation 21.

By considering the cold plasma limit $c_{\rm s}=0$ (equivalently, the limit of
magnetically-dominated plasmas), we derived an analytically tractable
dispersion relation and computed the instability growth rate and the range of
unstable shear velocities for several special cases (see Table 1). In the
limiting case of Alfvén speed approaching the speed of light, we found that
the system is still unstable as long as the shear speed is greater than the
lower bound in Equation 37.

Although most shear flows in the Universe are expected to have an asymmetric
configuration — e.g., the boundaries of AGN jets (Chow et al., 2023) — the
results obtained here for a symmetric setup in the magnetically-dominated
regime $v_{\rm{A}}\gg c_{\rm s}$ may have important implications for the
magnetospheres of neutron stars and black holes — both for single objects and
for merging binaries.

Binary neutron stars in the latest stages of inspiral are unlikely to be
tidally locked (Bildsten & Cutler, 1992). Due to the orbital motion, the
interface between the two magnetospheres can be modeled as a magnetically-
dominated shear layer. This configuration is analogous to the symmetric setup
we adopted in this work. In the final stages of the inspiral, when the
surfaces of the two neutron stars come into contact, the KHI and ensuing
turbulence have also been invoked to govern rapid magnetic field amplification
via dynamo processes (e.g., Price & Rosswog (2006); Zrake & MacFadyen (2013);
Kiuchi et al. (2015)).

For single objects, shear motions in the azimuthal direction (equivalently,
discontinuities in the profile of the angular velocity) have been revealed by
fully kinetic and fluid simulations, both for black hole magnetospheres (e.g.,
Figure 3 in El Mellah et al., 2022) and for pulsar magnetospheres (e.g.,
Figure 4 in Timokhin, 2006). In addition, in active pulsars, a discontinuity
in the poloidal flow velocity is expected to occur near the star, between the
open field lines (loaded with radially streaming pairs) and the closed
magnetosphere.

We conclude with a few caveats. The plane-parallel approach we employed is
applicable only if the width of the shear layer is small compared to other
length scales of the system (e.g., for binary neutron stars, their
separation). Also, we have assumed that the fluids on the two sides of the
interface have the same sound speed, magnetic field strength, and orientation.
In the magnetically-dominated regime, pressure balance requires the magnetic
field strength on both sides to be the same. However, the field direction will
generally be different. In fact, in nearly all of the astrophysical
applications discussed above, the surface of velocity shear also presents a
significant magnetic shear. Our simplifying assumptions will be relaxed in a
future work.

## Acknowledgements

We are grateful to H. Hakobyan, J. Nättilä, and especially J. Mahlmann for
useful discussions. L.S. acknowledges support from the Cottrell Scholars
Award, and DoE Early Career Award DE-SC0023015. L.S. and J.D. acknowledge
support from NSF AST-2108201. This research was facilitated by Multimessenger
Plasma Physics Center (MPPC), NSF grant PHY-2206609. J.D. is supported by a
Joint Columbia University/Flatiron Research Fellowship, research at the
Flatiron Institute is supported by the Simons Foundation.

## Data Availability

No new data were generated or analysed in support of this research.

## References

  * Anile (1990) Anile A. M., 1990, Relativistic fluids and magneto-fluids. Cambridge University Press
  * Berlok & Pfrommer (2019) Berlok T., Pfrommer C., 2019, MNRAS, 485, 908
  * Bildsten & Cutler (1992) Bildsten L., Cutler C., 1992, ApJ, 400, 175
  * Blandford & Pringle (1976) Blandford R. D., Pringle J. E., 1976, Monthly Notices of the Royal Astronomical Society, 176, 443–454
  * Blumen et al. (1975) Blumen W., Drazin P. G., Billings D. F., 1975, Journal of Fluid Mechanics, 71, 305
  * Bodo et al. (2004) Bodo G., Mignone A., Rosner R., 2004, Phys. Rev. E, 70, 036304
  * Bodo et al. (2013) Bodo G., Mamatsashvili G., Rossi P., Mignone A., 2013, MNRAS, 434, 3030
  * Bodo et al. (2016) Bodo G., Mamatsashvili G., Rossi P., Mignone A., 2016, MNRAS, 462, 3031
  * Bodo et al. (2019) Bodo G., Mamatsashvili G., Rossi P., Mignone A., 2019, MNRAS, 485, 2909
  * Chandrasekhar (1961) Chandrasekhar S., 1961, Hydrodynamic and Hydromagnetic stability. Clarendon Press
  * Chow et al. (2023) Chow A., Rowan M., Davelaar J., Sironi L., 2023, The Kelvin-Helmholtz instability at the boundary of relativistic magnetized jets (arXiv:2209.13699)
  * El Mellah et al. (2022) El Mellah I., Cerutti B., Crinquand B., Parfrey K., 2022, A&A, 663, A169
  * Ferrari et al. (1978) Ferrari A., Trussoni E., Zaninetti L., 1978, A&A, 64, 43
  * Ferrari et al. (1980) Ferrari A., Trussoni E., Zaninetti L., 1980, Monthly Notices of the Royal Astronomical Society, 193, 469–486
  * Hamlin & Newman (2013) Hamlin N. D., Newman W. I., 2013, Phys. Rev. E, 87, 043101
  * Kiuchi et al. (2015) Kiuchi K., Cerdá-Durán P., Kyutoku K., Sekiguchi Y., Shibata M., 2015, Phys. Rev. D, 92, 124034
  * Kivelson & Zu-Yin (1984) Kivelson M. G., Zu-Yin P., 1984, Planetary and Space Science, 32, 1335–1341
  * Lord Kelvin (1871) Lord Kelvin 1871, Philos. Mag., 42, 362
  * Mignone et al. (2018) Mignone A., Mattia G., Bodo G., 2018, Physics of Plasmas, 25, 092114
  * Osmanov et al. (2008) Osmanov Z., Mignone A., Massaglia S., Bodo G., Ferrari A., 2008, A&A, 490, 493
  * Pimentel & Lora-Clavijo (2019) Pimentel O. M., Lora-Clavijo F. D., 2019, MNRAS, 490, 4183
  * Prajapati & Chhajlani (2010) Prajapati R. P., Chhajlani R. K., 2010, Physics of Plasmas, 17, 112108
  * Price & Rosswog (2006) Price D. J., Rosswog S., 2006, Science, 312, 719
  * Pu & Kivelson (1983) Pu Z.-Y., Kivelson M. G., 1983, Journal of Geophysical Research, 88, 841
  * Rowan (2019) Rowan M. E., 2019, PhD thesis, Harvard University, Massachusetts
  * Roychoudhury & Lovelace (1986) Roychoudhury S., Lovelace R. V. E., 1986, ApJ, 302, 188
  * Sharma & Chhajlani (1998) Sharma P. K., Chhajlani R. K., 1998, Physics of Plasmas, 5, 625
  * Sobacchi & Lyubarsky (2018) Sobacchi E., Lyubarsky Y. E., 2018, MNRAS, 473, 2813
  * Sommerfeld (1912) Sommerfeld A., 1912, Jhrber. German. Math.-Verein, 12, 309
  * Timokhin (2006) Timokhin A. N., 2006, MNRAS, 368, 1055
  * Turland & Scheuer (1976) Turland B. D., Scheuer P. A., 1976, Monthly Notices of the Royal Astronomical Society, 176, 421–441
  * Von Helmholtz & Monats (1868) Von Helmholtz H., Monats K., 1868, Preuss. Akad. Wiss. Berlin, 23, 215
  * Weisstein (2020) Weisstein E. W., 2020, Quartic equation, {https://mathworld.wolfram.com/QuarticEquation.html}
  * Zrake & MacFadyen (2013) Zrake J., MacFadyen A. I., 2013, ApJ, 769, L29

##  Appendix A General solution for magnetically-dominated flows

For magnetically-dominated flows, we assume $c_{\rm s}=0$ so that
$v_{\rm{f}\perp}=v_{\rm{A}}$ and the dispersion relation (Equation 17) is
reduced to a quadratic equation in $(\phi/M)^{2}$ with solutions:

$\displaystyle\frac{\phi^{2}}{M^{2}}$
$\displaystyle=\frac{\mu_{1}\pm\sqrt{\mu_{2}}}{\mu_{3}}$ (43)

where

$\displaystyle\mu_{1}$ $\displaystyle=\cos^{2}\theta M_{\rm re}^{4}(1+M_{\rm
re}^{2})\gamma^{5}-\cos^{3}\theta\cos^{2}\Omega(M_{\rm
re}^{2})[\cos^{2}\Omega-2(1+M_{\rm
re}^{2})]\sin\theta(1-\gamma^{2})\gamma^{4}$
$\displaystyle\quad\quad+\cos^{8}\theta\sin^{2}\Omega(1-\gamma^{2})^{2}(\gamma^{3})[\cos^{2}\Omega(1+\gamma^{2})-\gamma^{2}]-\cos^{7}\theta\cos^{2}\Omega\sin\theta(1-\gamma^{2})^{2}(\gamma^{2})[\cos^{2}\Omega(1+3\gamma^{2})-3\gamma^{2}]$
$\displaystyle\quad\quad+\cos^{4}\theta M_{\rm
re}^{2}(1-\gamma^{2})(\gamma^{3})[\cos^{4}\Omega(2\sin^{2}\theta-1)-1-2\gamma^{2}+M_{\rm
re}^{2}(1-2\gamma^{2})+\cos^{2}\Omega(1+2(1+M_{\rm re}^{2})\gamma^{2})]$
$\displaystyle\quad\quad-\cos^{6}\theta(1-\gamma^{2})(\gamma^{3})[\gamma^{4}-1+M_{\rm
re}^{2}(1-3\gamma^{2}+\gamma^{4})+\cos^{2}\Omega(1+\gamma^{2}-2\gamma^{4}+M_{\rm
re}^{2}(1+4\gamma^{2}-2\gamma^{4}))$
$\displaystyle\quad\quad+\cos^{4}\Omega((1-\gamma^{2})(2\sin^{2}\theta-\gamma^{2})+M_{\rm
re}^{2}(\gamma^{4}-2-\gamma^{2}))]-\cos^{5}\theta\cos^{2}\Omega\sin\theta(1-\gamma^{2})(\gamma^{2})[(1-\gamma^{2})(1+\sin^{2}\Omega\gamma^{2})$
$\displaystyle\quad\quad+M_{\rm
re}^{2}(6\gamma^{2}-1-2\gamma^{4}+\cos^{2}\Omega(-1-5\gamma^{2}+2\gamma^{4}))]$
(44) $\displaystyle\mu_{2}$
$\displaystyle=\cos\theta^{4}\gamma^{4}[\cos\theta\cos^{2}\Omega\sin\theta\
(1-\gamma^{2})+M_{\rm
re}^{2}\gamma-\cos^{2}\theta\sin^{4}\Omega(1-\gamma^{2})\gamma]^{2}M_{\rm
re}^{4}(\cos^{4}\Omega+4M_{\rm re}^{2})\gamma^{4}$
$\displaystyle\quad\quad-4\cos\theta\cos^{2}\Omega M_{\rm
re}^{4}\sin\theta\gamma^{3}(\cos^{2}\Omega-2+2\gamma^{2})+\cos^{8}\theta(\cos^{2}\Omega-\gamma^{2}+\gamma^{4}-\cos^{2}\Omega\gamma^{4})^{2}$
$\displaystyle\quad\quad+4\cos^{7}\theta\cos^{2}\Omega\sin\theta(1-\gamma^{2})^{2}\gamma[\cos^{2}\Omega(1+\gamma^{2})-\gamma^{2}]+4\cos^{5}\theta\cos^{2}\Omega\sin\theta(1-\gamma^{2})\gamma[(1-\gamma^{2})^{2}+M_{\rm
re}^{2}(5\gamma^{2}-2\gamma^{4}-1)$
$\displaystyle\quad\quad-\cos^{2}\Omega(2-(3-6M_{\rm
re}^{2})\gamma^{2}+(1-2M_{\rm
re}^{2})\gamma^{4})]+4\cos^{3}\theta\cos^{2}\Omega M_{\rm
re}^{2}\sin\theta\gamma[(1-\gamma^{2})(2\sin^{4}\Omega\gamma^{4}-1+(-1+4\cos^{2}\Omega)\gamma^{2})$
$\displaystyle\quad\quad+M_{\rm
re}^{2}(1-4\gamma^{2}+2\gamma^{4}+\cos^{2}\Omega(1+\gamma^{2}))]-2\cos^{2}\theta
M_{\rm re}^{2}(\gamma^{2})[(\cos^{4}\Omega)[2-3\gamma^{2}+\gamma^{4}+M_{\rm
re}^{2}(1-2\sin^{2}\theta+\gamma^{2})]$ $\displaystyle\quad\quad-2M_{\rm
re}^{2}(1-\gamma^{2})(M_{\rm re}^{2}-2\gamma^{2}-1)+\cos^{2}\Omega(M_{\rm
re}^{2}(4\gamma^{4}-1-4\gamma^{2})-(1-\gamma^{2})^{2})]+2\cos^{6}\theta(1-\gamma^{2})[(\gamma^{2})[M_{\rm
re}^{2}(1-4\gamma^{2}+2\gamma^{4})$
$\displaystyle\quad\quad-(1-\gamma^{2})^{2}]+(\cos^{2}\Omega)[1+\gamma^{2}-4\gamma^{4}+2\gamma^{6}+M_{\rm
re}^{2}(2\gamma^{2}-1+7\gamma^{4}-4\gamma^{6})]+\cos^{4}\Omega((1-\gamma^{2})[\gamma^{4}-2-(1-2\sin^{2}\theta)\gamma^{2}]$
$\displaystyle\quad\quad+M_{\rm
re}^{2}(1-4\gamma^{2}-3\gamma^{4}+2\gamma^{6}))]+(\cos^{4}\theta)[(1-\gamma^{2})^{4}+M_{\rm
re}^{4}(1-8\gamma^{2}+16\gamma^{4}-8\gamma^{6})+2M_{\rm
re}^{2}(6\gamma^{2}-1-7\gamma^{4}+2\gamma^{8})$
$\displaystyle\quad\quad+(\cos^{4}\Omega)[M_{\rm
re}^{4}(1+\gamma^{2})^{2}+(2-3\gamma^{2}+\gamma^{4})^{2}-4M_{\rm
re}^{2}(1-\gamma^{2})(1-2(2-\sin^{2}\theta)\gamma^{2}+\gamma^{6})]+2(\cos^{2}\Omega)[(1-\gamma^{2})^{3}(\gamma^{2}-2)$
$\displaystyle\quad\quad+M_{\rm
re}^{4}(3\gamma^{2}-8\gamma^{4}+4\gamma^{6}-1)+M_{\rm
re}^{2}(3-14\gamma^{2}+13\gamma^{4}+2\gamma^{6}-4\gamma^{8})]]$ (45)
$\displaystyle\mu_{3}$ $\displaystyle=[M_{\rm
re}^{2}\gamma^{2}-\cos^{2}\theta(1-\gamma^{2})(\cos^{2}\Omega(1-\gamma^{2})+\gamma^{2})][2\cos^{4}\theta\sin^{2}\Omega(1-\gamma^{2})^{2}\gamma-2\cos^{3}\theta\cos^{2}\Omega\sin\theta(1-\gamma^{2})^{2}$
$\displaystyle\quad\quad+2\cos\theta\cos^{2}\Omega M_{\rm
re}^{2}\sin\theta(1-\gamma^{2})\gamma^{2}+M_{\rm
re}^{4}\gamma^{3}+\cos^{2}\theta M_{\rm
re}^{2}\gamma(\gamma^{2}-2+\gamma^{4}+\cos^{2}\Omega(1-\gamma^{4}))]$ (46)

##  Appendix B Phase speed of the unstable modes

In the main text, we presented three figures showing the dependence of the
imaginary part of $\phi$ (i.e., the growth rate) on $M_{\rm re}$,
$\cos\theta$, $\cos\Omega$, $v_{\rm{A}}$ and $c_{\rm s}$. In this Appendix, we
provide corresponding plots for the real part (i.e., the phase speed). In
section 4, we claimed that unstable modes associated with the fast
magnetosonic wave have vanishing phase speed. Figure 5, Figure 6 and Figure 7
plot the real part of $\phi$, whose imaginary part was presented in Figure 2,
Figure 3 and Figure 4 respectively. To guide the eye, we overplot the
boundaries of the unstable region associated with the fast mode (white dashed
and dotted lines). In all three figures, non-zero phase speeds are all located
outside of the fast-mode unstable boundary, indicating that the unstable modes
associated with the fast wave indeed have vanishing phase speed. In contrast,
unstable modes associated with the slow wave (more precisely, these modes
correspond to slow waves in one fluid, and fast waves in the other fluid) have
nonzero phase speed.

Figure 5: Dependence of the real part of $\phi$ on $M_{\rm re}$ and $\Omega$
when $\cos\theta=1$. The corresponding imaginary part of the solution is in
Figure 2. Figure 6: Dependence of the real part of $\phi$ on $M_{\rm re}$ and
$\theta$ when $\cos\Omega=1$. The corresponding imaginary part of the solution
is in Figure 3. Figure 7: Dependence of the real part of $\phi$ on $M_{\rm
re}$ and $\theta$ when $\cos\Omega=0$. The corresponding imaginary part is in
Figure 4.