#1 #1

11institutetext: Leibniz-Institut für Astrophysik Potsdam (AIP), An der
Sternwarte 16, D-14482 Potsdam, Germany
11email<EMAIL_ADDRESS>22institutetext: Niels Bohr Institute, University of
Copenhagen, Blegdamsvej 17, 2100 Copenhagen, Denmark

# Thermal conductivity with bells and whistlers: suppression of the
magnetothermal instability in galaxy clusters

Lorenzo Maria Perrone 11 Thomas Berlok 2211 Christoph Pfrommer 11

In the hot and dilute intracluster medium (ICM) in galaxy clusters, kinetic
plasma instabilities that are excited at the particle gyroradius may play an
important role in the transport of heat and momentum, thus affecting the
large-scale evolution of these systems. In this paper, we continue our
investigation of the effect of whistler suppression of thermal conductivity on
the magneto-thermal instability (MTI), which may be active in the periphery of
galaxy clusters and may contribute to the observed levels of turbulence. We
use a subgrid closure for the heat flux inspired from kinetic simulations and
show that MTI turbulence with whistler suppression exhibits a critical
transition as the suppression parameter is increased: for modest suppression
of the conductivity, the turbulent velocities generated by the MTI decrease
accordingly in agreement with scaling laws found in previous studies of the
MTI. However, for suppression above a critical threshold, the MTI loses its
ability to maintain equipartition-level magnetic fields through a small-scale
dynamo (SSD), and the system enters a “death-spiral”. We show that analogous
levels of suppression of thermal conductivity with a simple model of “flat”
uniform suppression would not inhibit the dynamo. We propose a model to
explain this critical transition, and speculate that conditions in the hot ICM
are favourable to the upkeep of the dynamo. We then look at spatial
correlations and energy transfers in spectral space, and find that with
whistler suppression, high-$\beta$ regions are brought out of thermal
equilibrium while the efficiency of MTI turbulent driving at the largest
turbulent scales is reduced. Finally, we show that external turbulence
interferes with the MTI and leads to reduced levels of MTI turbulence. While
individually both external turbulence and whistler suppression reduce MTI
turbulence, we find that they can exhibit a complex interplay when acting in
conjunction, with external turbulence boosting the whistler-suppressed thermal
conductivity and even reviving a “dead” MTI. Our study illustrates how
extending magnetohydrodynamics with a simple prescription for microscale
plasma physics can lead to the formation of a complicated dynamical system and
demonstrates that further work is needed in order to bridge the gap between
micro- and macro scales in galaxy clusters.

###### Key Words.:

Galaxies: clusters: intracluster medium – Instabilities – Magnetohydrodynamics
(MHD) – Plasmas – Turbulence – Methods: numerical

##  1 Introduction

In the hot and dilute plasma that fills galaxy clusters, the ICM, transport
processes are expected to be significantly modified compared to the standard
collisional picture due to the presence of magnetic fields (Sarazin, 1988).
First, despite the relatively large plasma $\beta$ (Carilli & Taylor, 2002;
Vogt & Ensslin, 2005; Kuchar & Enßlin, 2011), charged particles gyrate around
magnetic field lines on much shorter timescales than Coulomb scattering, and
as a result transport of heat and momentum take place primarily along the
local direction of the magnetic field (Braginskii, 1965). This results in a
new class of buoyancy instabilities that eschew the stable entropy
stratification of galaxy clusters, and which can perturb the ICM respectively
in the cluster periphery (the MTI; Balbus, 2000) or in the core of cool-core
clusters (the heat-flux driven buoyancy instability, or HBI; Quataert, 2008).
The presence of the MTI may have implications for the re-acceleration of
volume-filling, relativistic electrons far in the periphery of galaxy clusters
(Cuciti et al., 2022; Beduzzi et al., 2023). Second, high-$\beta$, weakly-
collisional plasmas like the ICM are a fertile ground for microinstabilities
that grow on kinetic scales (Schekochihin & Cowley, 2006; Schekochihin et al.,
2008) and that may substantially impede particle motion. Among these, we
mention the mirror, firehose and whistler instabilities (Chandrasekhar et al.,
1958; Parker, 1958; Gary, 1993; Gary & Li, 2000).

A large body of numerical work (Kunz et al., 2014; Riquelme et al., 2016;
Melville et al., 2016) and observational evidence from the solar wind (also a
dilute and magnetized plasma like the ICM; Bale et al., 2009; Chen et al.,
2016) suggests that these microinstabilities act in such a way as to counter
the pressure anisotropy and to enhance the scattering of charged particles,
leading to a suppression of the effective viscosity and thermal conductivity
(Bale et al., 2013; Roberg-Clark et al., 2016, 2018; Komarov et al., 2018).
The implications for the thermodynamic evolution of galaxy clusters are
potentially far-reaching, and include suppressed dissipation of turbulence
(Zhuravleva et al., 2019), plasma heating (Ley et al., 2022) and explosive
magnetic field amplification (Mogavero & Schekochihin, 2014; Zhou et al.,
2023; Rappaz & Schober, 2023). Of particular importance is the impact on
thermal conductivity, which may play a role in regulating the radiative losses
in cool-core clusters (Narayan & Medvedev, 2001; Zakamska & Narayan, 2003;
Kunz et al., 2012; Jacob & Pfrommer, 2017a, b), and whose suppression may
explain certain observed features such as cold fronts (Markevitch & Vikhlinin,
2007; ZuHone et al., 2013; ZuHone & Roediger, 2016). Furthermore, the
suppression of thermal conductivity can potentially lead to an impairment of
the MTI, since the instability relies on fast conduction along magnetic field
lines and drives turbulent motions that scale with the value of thermal
diffusivity (Perrone & Latter, 2022a, b, hereafter ; ).

###  1.1 The role of microscale instabilities

Predicting the impact of microscale instabilities and the resulting
suppression of heat flux on the MTI is not straightforward, and modeling self-
consistently the saturation of microinstabilities within large-scale
simulations is currently unfeasible due to the vast scale separation. For this
reason, several studies have sought to capture the effect of individual
microinstabilities through effective models inspired by theoretical
considerations as well as by kinetic plasma simulations. For instance, Drake
et al. (2021) suggest through theoretical arguments that whistler suppression
of the heat flux could severely constrain the scales where the MTI can
efficiently grow in the periphery of galaxy clusters. However, because of the
complexity of the physics and the disparate scales involved, these
considerations warrant a more detailed investigation and leave the door open
to potential surprises. In the case of the mirror instability, for example,
magnetohydrodynamics (MHD) simulations have indicated that the overall impact
on the MTI may be comparatively small (Berlok et al., 2021), because regions
that are mirror-unstable (where conduction is suppressed) coincide with those
of fastest linear MTI growth. This demonstrates that the MTI can minimise, to
some degree, the onset and impact of deleterious micro-instabilities, and
calls into questions simple models where the thermal conductivity is reduced
by some fixed uniform factor. Finally, the existence of a kinetic counterpart
of the MTI that grows on electron kinetic scales (and possibly competes with
other high-$\beta$ microinstabilities) further complicates the picture (Xu &
Kunz, 2016).

###  1.2 The role of external turbulence

The MTI is not the only source of turbulence in galaxy clusters on large
scales: mergers and (substructure) accretion can stir turbulent motions
(Churazov et al., 2003; Schuecker et al., 2004), potentially competing with
the MTI and suppressing the instability (Sharma et al., 2009; Ruszkowski & Oh,
2010; McCourt et al., 2011), although hints of additive behaviour between the
MTI and other sources of turbulence have also been reported (Parrish et al.,
2012). The interplay between MTI and other sources of turbulence has been
previously investigated with local and global numerical setups (McCourt et
al., 2011; Parrish et al., 2012), which showed that strong external turbulence
can dominate over MTI-driven turbulence if the the eddy turnover time of the
external turbulence is much shorter than the MTI e-folding time. In practice,
however, given the broad range of scales where the MTI can grow close to its
maximum rate this is unlikely to happen at all scales up to the pressure
scale-height (McCourt et al., 2011). Some open questions remain, particularly
concerning the impact of external forcing on the overall buoyant instability
of the ICM with anisotropic thermal conduction (Sharma et al., 2009).

A number of large-scale cosmological simulations looking at anisotropic heat
conduction in galaxy clusters have also been performed (Ruszkowski & Oh, 2010;
Ruszkowski et al., 2011; ZuHone et al., 2013; Yang & Reynolds, 2016; Kannan et
al., 2017; Barnes et al., 2019), but have generally not found strong evidence
of MTI action against the turbulence produced by cosmological accretion
(Ruszkowski et al., 2011; Kannan et al., 2017; Barnes et al., 2019). However,
several questions remain on whether these works adequately resolve the broad
range of scales ($10-100\,$\mathrm{k}\mathrm{p}\mathrm{c}$$ in the periphery
of galaxy clusters) that are key to energy injection by the MTI, without which
the instability would be numerically inhibited. Similarly, high levels of
numerical heat diffusion in the low-resolution regions in the outskirts of
galaxy clusters will also artificially degrade the growth of the MTI (PL22b).
Carefully constructed large-scale cosmological simulations that adequately
resolve the MTI scales are still required in order to draw firm conclusions on
the relevance of the MTI in a cosmological context.

###  1.3 Aim and plan of the paper

In Perrone et al. (2023, hereafter PBP23), we set out to explore the impact of
the whistler instability on the MTI through small-scale Boussinesq simulations
where we employed a subgrid closure for the heat flux. We found that,
depending on the exact parameters of the local environment, the effect on the
MTI can range from a mild to a complete suppression of the instability. In
particular, we identified the shutdown of the SSD as the culprit behind the
death of the MTI in our simulations. We then focused on the interplay between
external turbulence and whistler suppression and showed that, with none or
little whistler suppression, externally-driven turbulence interferes with but
does not wipe out the MTI, and that the flow remains overall buoyantly
unstable. Contrary to expectations, external turbulence was shown to actually
play a beneficial role by reviving the MTI in strongly whistler suppressed
simulations.

In this work we complement the analysis of PBP23 along three main scientific
drivers: (a) we explain the mechanism behind the critical transition of MTI
turbulence with whistler suppression, (b) we explore in what ways whistler
suppression is different to a flat uniform suppression of thermal
conductivity, and (c) we examine the interaction of whistler-suppressed MTI
turbulence with externally-driven turbulence. Our aim is to add another
building block to the ongoing conversation of how macro- and micro-turbulence
interact in the hot dilute ICM, and to inform future theoretical and numerical
work.

This paper is structured as follows: in Sect. 2 we briefly review how the
presence of magnetic fields and the low collisionality in the ICM affect the
standard collisional picture of particle transport, and we introduce our
subgrid model for the suppression of thermal diffusivity due to whistler
waves; in Sect. 3 we detail our numerical methods, the parameters of the
simulations and the diagnostics employed; in Sect. 4 we present the first
results of our simulations, and focus on the weakening of the MTI with
whistler suppression and the appearance of a critical transition above which
MTI-driven turbulence dies entirely; in Sect. 5 we take a closer look at the
backreaction of whistler suppression on the flow and on how the MTI energy
fluxes in spectral space are affected; in Sect. 6 we re-examine the role of
external turbulence and show that it can interfere with the workings of the
MTI but at the same time be beneficial for simulations where the MTI
turbulence is dying. Finally, we discuss in what respects whistler suppression
of thermal diffusivity can be approximated by a flat uniform suppression below
the Spitzer value (Sect. 7) and we draw our conclusions in Sect. 8. The main
text is complemented by three appendices as follows: in Appendix A we rewrite
the subgrid model of whistler suppression in terms of Boussinesq variables; in
Appendix B we show through linear theory how whistler suppression can
interrupt the growth of an MTI eigenmode; finally, in Appendix C we detail and
validate through benchmarks the numerical implementation of whistler
suppression in our code.

##  2 Theoretical preliminaries

###  2.1 Collisional transport in the ICM

In the standard picture of transport processes in plasmas, transport of heat
is mediated by particle-particle Coulomb collisions (Spitzer, 1962), with
electrons (the main particle species responsible for the conduction of heat)
scattering off each other every $\lambda_{\text{mfp,e}}$, while traveling at
the thermal speed $\varv_{\text{th,e}}$. This gives the so-called Spitzer
diffusivity (Spitzer, 1962)

$\displaystyle\chi_{0}\simeq 4.98\times
10^{31}\left(\frac{T_{\mathrm{e}}}{5\leavevmode\nobreak\
$\mathrm{k}\mathrm{e}\mathrm{V}$}\right)^{5/2}\left(\frac{n_{\mathrm{e}}}{10^{-3}\leavevmode\nobreak\
$\mathrm{c}\mathrm{m}^{-3}$}\right)^{-1}$\mathrm{c}\mathrm{m}^{2}\,\mathrm{s}^{-1}$,$
(1)

(here and in the rest of the paper we express temperature in energy units),
which we have normalized for typical values in the periphery of galaxy
clusters. However, despite having subthermal magnetic fields
($\beta=p_{\mathrm{th}}/p_{\mathrm{mag}}\gtrsim 100$; Carilli & Taylor, 2002),
the ICM is strongly magnetized in the sense that the collisional mean-free
path of charged particles is much larger than their gyroradius. In fact, for
typical conditions in the ICM, the electron mean-free path is

$\displaystyle\lambda_{\text{mfp,e}}\simeq
5.83\left(\frac{n_{\mathrm{e}}}{10^{-3}\leavevmode\nobreak\
$\mathrm{c}\mathrm{m}^{-3}$}\right)^{-1}\left(\frac{T_{\mathrm{e}}}{5\leavevmode\nobreak\
$\mathrm{k}\mathrm{e}\mathrm{V}$}\right)^{2}$\mathrm{k}\mathrm{p}\mathrm{c}$,$
(2)

which is many orders of magnitude larger than their gyroradius

$\displaystyle\rho_{\text{e}}\simeq 2.73\times
10^{-14}\left(\frac{B}{2\leavevmode\nobreak\
$\mu\mathrm{G}$}\right)^{-1}\left(\frac{T_{\mathrm{e}}}{5\leavevmode\nobreak\
$\mathrm{k}\mathrm{e}\mathrm{V}$}\right)^{1/2}$\mathrm{k}\mathrm{p}\mathrm{c}$.$
(3)

This means that particles collide against each other very infrequently as they
gyrate along magnetic field lines, and transport across the local magnetic
field is suppressed. As a result, the heat flux and viscous tensor become
anisotropic with respect to the local direction of the magnetic field
(Braginskii, 1965). In particular, the anisotropic (or Braginskii) heat flux
along the magnetic field takes the form:

$\displaystyle q_{\parallel,\mathrm{s}}\simeq
n_{\mathrm{e}}\chi_{0}(\bm{b}\bm{\cdot}\bm{\nabla})T_{\mathrm{e}},$ (4)

which shows that only temperature gradients along the magnetic field can be
thermalized. This property of dilute, magnetized plasmas like the ICM may help
explain the presence of features such as cold fronts, which are characterized
by sharp temperature gradients on spatial scales comparable or shorter than
the electron mean-free path (see, e.g., Markevitch & Vikhlinin, 2007).
Simulations find that magnetic fields are typically draped around the cold
fronts and the suppression of heat conduction across the front, that arises
with an anisotropic conduction and this field configuration, has been invoked
to explain cold front observations (Lyutikov, 2006; Asai et al., 2007; Dursi &
Pfrommer, 2008; Pfrommer & Dursi, 2010), although magnetic draping alone may
not be sufficient (ZuHone et al., 2013; ZuHone & Roediger, 2016).

###  2.2 Whistler suppression of thermal conductivity

The Braginskii form of the heat flux and viscous tensor is derived under the
assumption that the plasma is collisional on macroscopic (fluid) scales. In
many astrophysical plasmas, including the ICM, this requirement is however
only marginally satisfied. As a consequence, a host of kinetic instabilities
can grow on small scales (comparable to the particle gyroradius) at a much
faster rate than any macroscopic process (see, e.g., Schekochihin et al.,
2008, and references therein). The resulting electromagnetic fluctuations can
significantly affect particle trajectories by enhancing their scattering rate,
effectively leading to a suppression of the plasma viscosity and thermal
conductivity (Bale et al., 2013; Roberg-Clark et al., 2016, 2018; Komarov et
al., 2018). In particular, the whistler instability (which can be excited in
high-$\beta$ plasmas by an electron heat flux; Levinson & Eichler 1992, or by
an electron temperature anisotropy; Gary & Wang 1996) has long been identified
as a leading mechanism for the suppression of electron transport in weakly
collisional plasmas (Levinson & Eichler, 1992; Pistinner & Eichler, 1998), and
could also appear concomitantly to ion-scale microinstabilities, such as the
mirror instability (Ley et al., 2023). Theoretical arguments and kinetic
particle-in-cell simulations showed that at saturation the whistler
instability establishes a marginal heat flux which is independent of the
temperature gradient and suppressed by a factor of $\sim 1/\beta_{\mathrm{e}}$
(where $\beta_{\mathrm{e}}$ is the electron plasma beta) compared to the value
it would attain solely due to the streaming of electrons in the absence of
electromagnetic fluctuations (Riquelme et al., 2016; Roberg-Clark et al.,
2018; Komarov et al., 2018):

$\displaystyle q_{\parallel,\mathrm{w}}\simeq
1.5n_{\mathrm{e}}m_{\mathrm{e}}\varv_{\mathrm{th},\mathrm{e}}^{3}/\beta_{\mathrm{e}}.$
(5)

Since in the ICM $\beta\gtrsim 100$, the suppression of the heat flux could
then be significant.

As it is currently unfeasible to self-consistently model the saturation of the
whistler instability in macroscopic fluid simulations of the MTI, we resolve
on adopting the following subgrid formula for the heat flux, which is inspired
by the work of Komarov et al. (2018) and smoothly interpolates between the
collisional heat flux (Eq. 4) and the marginal heat flux allowed in the
whistler-dominated regime (Eq. 5):

$\displaystyle
q_{\parallel}=\left(\frac{1}{q_{\parallel,\mathrm{s}}}+\frac{1}{q_{\parallel,\mathrm{w}}}\right)^{-1}=\frac{q_{\parallel,\mathrm{s}}}{1+\frac{1}{3}\beta_{\mathrm{e}}\lambda_{\text{mfp,e}}/L_{T,\parallel}},$
(6)

where $L_{T,\parallel}\equiv\left\lvert\bm{b}\bm{\cdot}\bm{\nabla}\ln
T_{\text{e}}\right\rvert^{-1}$ is the temperature gradient scale parallel to
the magnetic field. Since in the periphery of galaxy clusters the plasma is
relatively collisional ($\lambda_{\text{mfp,e}}/L_{T,\parallel}\lesssim 0.1$,
see Table 1), in our model we do not include saturation of the heat flux due
to the free-streaming of electrons (Cowie & McKee, 1977), but depending on the
regime of interest other interpolation formulae can be devised (see, e.g.,
Beckmann et al., 2022).

The closure for the heat flux in Equation (6) has the advantage of tying
together the dynamics on small and large scales, since the large-scale fluid
motions determine the degree of whistler suppression through the plasma
$\beta$, and are themselves affected by the resulting parallel heat flux. In
doing so, it implicitly assumes that, due to the vast separation of temporal
and spatial scales, the saturation of the whistler instability is effectively
instantaneous with respect to the large-scale MTI motions. While this seems
plausible on the macroscopic fluid scales evolved in our simulations, we note
that the MTI is not necessarily confined to those scales but exists all the
way down to electron kinetic scales, where the whistler instability lives (Xu
& Kunz, 2016). At this point, the simple model in Eq. (6) is harder to
justify, since the two instabilities could interfere with one another.
Finally, our model neglects the combined impact of other microscale
instabilities, such as the mirror and firehose, which could also change the
picture (Komarov et al., 2018, although individually they are likely
subdominant compared to the whistler instability).

###  2.3 Theory of MTI turbulence in galaxy clusters’ periphery

In this and the next section we briefly review the theory of MTI turbulence
and the Boussinesq approximation. Some readers may choose to skip ahead to
Sect. 2.5.

Since the seminal work of Balbus (2000) it was recognized that, in addition to
modifying particle transport on small scales, anisotropic thermal conduction
can also affect the ICM dynamics on large scales, triggering a new class of
buoyancy instabilities. In fact, with heat conduction along magnetic field
lines, the hydrodynamic stability to convection is determined by the sign of
the temperature gradient (Balbus, 2000), rather than the entropy gradient
(Schwarzschild, 1906). Under these conditions, two new instabilities arise:
the HBI (Quataert, 2008), which destabilizes plasmas with a positive
temperature gradient with respect to gravity and may be relevant in the
central region of cool-core clusters, and the MTI, which can destabilize the
periphery of galaxy clusters (Parrish & Stone, 2007; McCourt et al., 2011;
Kunz, 2011), and is the main focus of this paper. The MTI operates by
extracting energy from the background temperature gradient via conduction of
heat along magnetic field lines, which are frozen-in and dragged along with
the fluid. As a result, a displaced fluid element remains thermally connected
with the plasma at its original location, and to get into pressure balance
with its new environment expands (if rising) or contracts (if sinking),
becoming buoyantly unstable.

In PL22a and PL22b it was shown that the MTI produces a state of buoyancy-
driven turbulence sustained by energy drawn from the background temperature
gradient and composed of density and velocity fluctuations over a wide range
of scales, with turbulent properties at saturation that follow clear power
laws. In the stably stratified regime typical of the ICM, the rms velocities
obey

$\displaystyle
u_{\text{rms}}^{2}\propto\dfrac{\chi\omega_{\mathrm{T}}^{3}}{N^{2}},$ (7)

where $\omega_{\mathrm{T}}$ is the MTI frequency, i.e., the maximum growth
rate of the instability, which also corresponds to the typical timescale of
MTI turbulence, and $N$ is the Brunt-Väisälä frequency, which determines the
buoyant response of the system in a stably-stratified environment. Locally (at
a given spherical radius $R_{0}$), they are given by:

$\displaystyle\omega_{\mathrm{T}}^{2}=-g_{0}\left.\frac{\partial\ln
T}{\partial R}\right|_{0},\quad
N^{2}=\frac{g_{0}}{\gamma}\left.\frac{\partial\ln p\rho^{-\gamma}}{\partial
R}\right|_{0},\qquad$ (8)

where $\rho,p,T$ are the fluid density, pressure and temperature,
respectively, $g_{0}=-(1/\rho)\mathrm{d}p/\mathrm{d}R|_{0}$ is the local
gravitational acceleration (here and in Eq. 8 the derivatives are computed at
$R_{0}$), and $\gamma$ is the adiabatic index. In the periphery of galaxy
clusters their values are comparable, and of the order of
$(600\,$\mathrm{M}\mathrm{y}\mathrm{r}$)^{-1}\simeq 5.27\times
10^{-17}\,$\mathrm{s}^{-1}$$ (PL22a). The largest scale excited by MTI
turbulence (the peak of the kinetic energy spectrum, or the integral scale) is
set by a balance between energy injection by the MTI and the stabilizing
entropy stratification – hence the moniker of “buoyancy scale” – and on
theoretical grounds was predicted to scale as111Note, however, that in the
three-dimensional simulations of Perrone & Latter (2022b) the scaling of
$l_{\mathrm{i}}$ with $\chi$ was found to be somewhat shallower than $1/2$, a
discrepancy which was attributed to the limited scale separation between
injection and viscous dissipation.

$\displaystyle
l_{\mathrm{i}}\propto\dfrac{(\chi\omega_{\mathrm{T}})^{1/2}}{N}.$ (9)

A second important scale in MTI turbulence is the conduction length

$\displaystyle l_{\chi}\equiv(\chi/\omega_{\mathrm{T}})^{1/2},$ (10)

which represents how far heat can diffuse on timescales of
$\omega_{\mathrm{T}}^{-1}$, the turnover time of MTI turbulence. The integral
and conduction scales differ by a factor of $\omega_{\mathrm{T}}/N$, which
however is of order unity in the periphery of galaxy clusters, and the two are
thus roughly comparable. Turbulence on scales below the integral scale is
characterized by vertically elongated eddies excited by the buoyancy force.
These eddies are efficient at transporting heat through advection and they
drive an SSD to approximate equipartition with the kinetic energy. On scales
$\gtrsim l_{\mathrm{i}}$, the stabilizing effect of the entropy stratification
is stronger, and the turbulent motions become relatively isotropic. On even
longer scales, vertical motions are effectively suppressed by the stable
entropy stratification. Finally, combining the scalings of the rms velocities
with the integral scale allowed PL22b to define a Reynolds number of MTI
turbulence $\text{Re}_{\mathrm{i}}$ as

$\displaystyle\text{Re}_{\mathrm{i}}\propto\dfrac{\chi\omega_{\mathrm{T}}^{2}}{\nu
N^{2}}.$ (11)

For typical ICM conditions in the periphery of galaxy clusters and assuming
unsuppressed to moderately suppressed thermal diffusivity with respect to the
Spitzer value (Eq. 1), the MTI is capable of driving turbulent motions of
several hundred kilometers per second, on scales of $\sim
100\,$\mathrm{k}\mathrm{p}\mathrm{c}$$, thus potentially contributing to the
observed levels of turbulence (PL22b).

###  2.4 Governing equations

In this work we adopt the Boussinesq approximation to simulate the subsonic
and small-scale dynamics of fluctuations within a quasi-spherical cluster
characterised by a background density $\rho_{0}$, pressure $p_{0}$, and
temperature $T_{0}$. By small scales, we here mean on scales much smaller than
the typical pressure (or density) scale-height in the periphery of galaxy
clusters, typically several hundred $\mathrm{k}\mathrm{p}\mathrm{c}$. The main
advantage of the Boussinesq approximation is that the fast-propagating sound
waves are screened out of the fluid equations, leading to a considerable
computational speed-up. As a result, the fluid is in pressure equilibrium at
all times, and temperature and density fluctuations are related by

$\displaystyle\frac{\delta T}{T_{0}}=-\frac{\delta\rho}{\rho_{0}}.$ (12)

Our model describes a small Cartesian block of ICM plasma located at a given
spherical radius $R_{0}$, where the vertical direction is along the radius and
antiparallel to the gravitational acceleration $\bm{g}=-g_{0}\bm{e}_{z}$. The
background radial stratification in both entropy and temperature is included
via the (constant-in-$z$) squared Brunt-Väisälä, $N^{2}$, and squared MTI,
$\omega_{\mathrm{T}}^{2}$, frequencies. Their values are computed at the
reference radius, $R_{0}$, using Eq. (8). For further discussion on the
specifics and limitations of the Boussinesq approximation, see PL22a.

The Boussinesq variables evolved in our simulations are the velocity
fluctuation $\bm{u}$, the magnetic field $\bm{B}$ (rescaled by
$\sqrt{4\pi\rho_{0}}$ so that it has the dimensions of a velocity), and the
buoyancy variable $\theta$

$\displaystyle\theta\equiv g_{0}\frac{\delta\rho}{\rho_{0}},$ (13)

(in units of $\mathrm{c}\mathrm{m}\,\mathrm{s}^{-2}$). The appropriate
equations are derived in PL22a, and are shown hereafter including a spatially-
varying anisotropic thermal diffusivity and a driving random external
acceleration $\bm{f}$ in the velocity equation. They read:

$\displaystyle\bm{\nabla}\bm{\cdot}\bm{u}=\bm{\nabla}\bm{\cdot}\bm{B}$
$\displaystyle=0,$ (14)
$\displaystyle\left(\partial_{t}+\bm{u}\bm{\cdot}\bm{\nabla}\right)\bm{u}$
$\displaystyle=-\frac{\bm{\nabla}p_{\text{tot}}}{\rho_{0}}-\theta\bm{e}_{z}+\left(\bm{B}\bm{\cdot}\bm{\nabla}\right)\bm{B}+\nu\bm{\nabla}^{2}\bm{u}+\bm{f},$
(15)
$\displaystyle\left(\partial_{t}+\bm{u}\bm{\cdot}\bm{\nabla}\right)\bm{B}$
$\displaystyle=\left(\bm{B}\bm{\cdot}\bm{\nabla}\right)\bm{u}+\eta\bm{\nabla}^{2}\bm{B},$
(16)
$\displaystyle\left(\partial_{t}+\bm{u}\bm{\cdot}\bm{\nabla}\right)\theta$
$\displaystyle=N^{2}u_{z}+\bm{\nabla}\bm{\cdot}\left[\chi\bm{b}\left(\bm{b}\bm{\cdot}\bm{\nabla}\right)\theta\right]+\omega_{\mathrm{T}}^{2}\bm{\nabla}\bm{\cdot}\left(\chi\bm{b}b_{z}\right)$
$\displaystyle+\chi_{\text{hyp}}\bm{\nabla}^{6}\theta,$ (17)

where $p_{\text{tot}}$ is the sum of thermal and magnetic pressure; $\nu$ and
$\eta$ are the fluid viscosity and magnetic resistivity (both in units of
$\mathrm{c}\mathrm{m}^{2}\,\mathrm{s}^{-1}$), while in Eq. (17) the thermal
diffusivity $\chi$ (also in $\mathrm{c}\mathrm{m}^{2}\,\mathrm{s}^{-1}$)
absorbs a factor of $(\gamma-1)/\gamma$. Note that in the buoyancy equation we
have also included a sixth order diffusion operator where $\chi_{\text{hyp}}$
is the corresponding hyperdiffusivity (in units of
$\mathrm{c}\mathrm{m}^{6}\,\mathrm{s}^{-1}$). We refer to Appendix C for a
detailed description of our implementation of whistler suppression and a
discussion of why hyperdiffusivity is required to regularize the smallest
scales when whistler suppression causes strong spatial variations in the
thermal diffusivity.

###  2.5 Boussinesq model of whistler suppression

In this section we introduce the final form of the whistler-suppressed thermal
diffusivity that we employ in our numerical model. In terms of Boussinesq
variables, Eq. (6) can be rewritten as follows (see Appendix A for further
discussion):

$\displaystyle\chi\equiv\dfrac{\chi_{0}}{1+\frac{1}{3}\alpha\tilde{\beta}\sigma},$
(18)

where we introduced the following three nondimensional quantities:

$\displaystyle\alpha\equiv\frac{1}{\gamma}\left(\frac{\lambda_{\text{mfp,e}}}{H}\right)\left(\frac{H^{2}}{\chi_{0}/\omega_{\mathrm{T}}}\right),$
(19)
$\displaystyle\tilde{\beta}\equiv\left(\frac{2\varv_{\mathrm{cond}}^{2}}{\varv_{\mathrm{A}}^{2}}\right),$
(20) $\displaystyle\sigma\equiv\left\lvert
b_{z}+\frac{\bm{b}\bm{\cdot}\bm{\nabla}\theta}{\omega_{\mathrm{T}}^{2}}\right\rvert.$
(21)

In Eq. (19) the $\alpha$ parameter encodes the physical information of the
volume of ICM plasma that we simulate. More precisely, it represents the
combined effect of plasma collisionality (given by the ratio
$\lambda_{\text{mfp,e}}/H$, the Knudsen number, where $H$ is the pressure
scale-height) and of the efficiency of thermal conduction on large (cluster-
size) scales over an MTI e-folding time, as expressed by the squared ratio of
$H$ to the unsuppressed conduction length (see Eq. 10). This parameter is
small when the plasma has large collisionality or is ineffective at conducting
heat. Conversely, it is large when the plasma becomes less and less
collisional, or when heat transport is efficient on large scales. In our
subgrid model for whistler suppression, $\alpha$ is the only free parameter
and is fixed at the beginning of each run. The second quantity in Eq. (20) is
a modified plasma $\tilde{\beta}$ that plays the same role as the usual plasma
beta but with a “conduction speed” $\varv_{\mathrm{cond}}$, defined as

$\displaystyle\varv_{\mathrm{cond}}\equiv(\chi_{0}\omega_{\mathrm{T}})^{1/2},$
(22)

which replaces the ion thermal velocity. As we show below, the two betas are
roughly comparable. Finally, in Eq. (21) the isothermality parameter $\sigma$
is a dimensionless quantity that represents the ratio of the background to the
total (background $+$ perturbation) temperature gradient parallel to the
magnetic field. It can be understood as follows: assuming for simplicity
purely radial magnetic fields, $\sigma$ is equal to unity at quiescence (no
temperature fluctuations), while in a turbulent medium it is greater (smaller)
than unity if the gradient of the local temperature fluctuation is aligned
(anti-aligned) with the background gradient. When the two gradients exactly
cancel each other ($\sigma=0$), the plasma is locally isothermal. In MTI
turbulence at saturation $\sigma$ is smaller than unity, since the MTI
effectively enforces isothermality along magnetic field lines (particularly so
on scales below the conduction length $l_{\chi}$, see PL22b, ). Much like
$\tilde{\beta}$ is the Boussinesq counterpart of the usual plasma $\beta$, the
product $\alpha\sigma$ plays the role of the ratio
$\lambda_{\text{mfp,e}}/L_{T,\parallel}$ in Eq. (6).

###  2.6 Reference parameters in the ICM

We now estimate reference values for the suppression parameter $\alpha$, and
proceed to quantify the impact of whistler suppression of thermal diffusivity
on the MTI and its saturated turbulent state. These estimates will then inform
the choice of parameters of our simulations (Sect. 3.2). Using the collisional
mean-free path (Eq. 2) and the Spitzer diffusivity (Eq. 1) we estimate the
suppression parameter $\alpha$ to be of the order of

$\displaystyle\alpha\simeq
0.01\left(\frac{H}{300\,$\mathrm{k}\mathrm{p}\mathrm{c}$}\right)\left(\frac{\omega_{\mathrm{T}}}{(700\,$\mathrm{M}\mathrm{y}\mathrm{r}$)^{-1}}\right)\left(\frac{T_{\mathrm{e}}}{5\,$\mathrm{k}\mathrm{e}\mathrm{V}$}\right)^{-1/2}.$
(23)

To evaluate $H$ we use the best-fit profiles for the X-COP sample reported by
Ghirardini et al. (2019), and take a fiducial radius (where $H$ is to be
computed) of the order of several hundred $\mathrm{k}\mathrm{p}\mathrm{c}$. By
doing so, we find that the realistic range of $\alpha$ allowed by observations
is 0.01–0.05, where the lower estimate is obtained from Eq. (23) with the
values given, while the higher estimate is calculated assuming
$H=600\,$\mathrm{k}\mathrm{p}\mathrm{c}$$,
$\omega_{\mathrm{T}}=(600\,$\mathrm{M}\mathrm{y}\mathrm{r}$)^{-1}$ and
$T_{\mathrm{e}}=1\,$\mathrm{k}\mathrm{e}\mathrm{V}$$. It is important to note
that $\alpha$ can have a large spread in values between different regions of
the cluster periphery, and across different clusters. Nevertheless, we believe
that the quoted numbers represent a realistic order-of-magnitude estimate of
$\alpha$ in the ICM.

The modified $\tilde{\beta}$ can be similarly estimated using Eq. (22) and the
Spitzer diffusivity (Eq. 1),

$\displaystyle\tilde{\beta}\simeq
24\left(\frac{\omega_{\mathrm{T}}}{(700\,$\mathrm{M}\mathrm{y}\mathrm{r}$)^{-1}}\right)\left(\frac{T_{\mathrm{e}}}{5\,$\mathrm{k}\mathrm{e}\mathrm{V}$}\right)^{5/2}\left(\frac{B}{2\,\mu$\mathrm{G}\,$}\right)^{-2},$
(24)

which could vary by up to an order of magnitude ($\tilde{\beta}\simeq 2-24$,
where the lower bound was computed assuming
$T_{\mathrm{e}}=1\,$\mathrm{k}\mathrm{e}\mathrm{V}$,B=1\,$\mu\mathrm{G}$$) due
to the strong dependence on the temperature and on the magnetic field
strength. Assuming typical conditions in the periphery of galaxy clusters
($n_{\mathrm{e}}=10^{-3}$\mathrm{c}\mathrm{m}^{-3}$,T_{\mathrm{e}}=5\,$\mathrm{k}\mathrm{e}\mathrm{V}$,B=2\,$\mu\mathrm{G}$$),
$\tilde{\beta}$ is approximately $1/2$ of the electron plasma $\beta$ (see
Table 1). Since for MTI turbulence $\sigma$ is less than unity, we can expect
that the overall suppression of thermal diffusivity by whistler waves is
relatively modest, with an average diffusivity of the order of $71-99\%$ of
the full Spitzer value (assuming conservatively $\sigma=1$), which is in
agreement with the estimates presented in Komarov et al. (2018). While this
level of suppression may not appear particularly significant, we note that
these estimates do not take into account the nonlinear backreaction of
whistler suppression onto the MTI turbulence: as we will show in the next
sections, even mild levels of suppression can have a strong impact on MTI-
driven turbulence. A summary of all the definitions and the numerical values
of the quantities introduced in this Section can be found in Table 1.

Table 1: Physical parameters of the ICM in the periphery of galaxy clusters. For simplicity we assume quasi-neutrality ($n_{\mathrm{e}}\simeq n_{\mathrm{i}}$) and equal temperature between ions and electrons ($T_{\mathrm{e}}\simeq T_{\mathrm{i}}\simeq T$). The numbers quoted have been computed taking $n_{\mathrm{e}}=10^{-3}$\mathrm{c}\mathrm{m}^{-3}$,T_{\mathrm{e}}=5\,$\mathrm{k}\mathrm{e}\mathrm{V}$,B=2\,$\mu\mathrm{G}$$, except for the suppression parameter $\alpha$ and the modified plasma beta $\tilde{\beta}$, for which we quote the realistic range of values that they could take in the periphery of galaxy clusters as discussed in the main text. The values of $\omega_{\mathrm{T}}$ and $N$ are taken from PL22a, using data from Ghirardini et al. (2019). The ‡ symbol means that the quantity has been evaluated assuming diffusivity equal to the Spitzer value. To compute $H$ we assume a power-law profile with slope $-2$ to $-3$ and evaluate it at $R=1\,$\mathrm{M}\mathrm{p}\mathrm{c}$$.  Quantity |  Definition |  Value
---|---|---
$\chi_{0}$ |  Eq. (1) |  $1.65\times 10^{2}\;$\mathrm{k}\mathrm{p}\mathrm{c}^{2}\,\mathrm{M}\mathrm{y}\mathrm{r}^{-1}$$
$\omega_{\mathrm{T}}$ |  Eq. (8) |  $(600\,$\mathrm{M}\mathrm{y}\mathrm{r}$)^{-1}\simeq 5.3\times 10^{-17}$\mathrm{s}^{-1}$$
$N$ |  Eq. (8) |  $\simeq\omega_{\mathrm{T}}$
$\rho_{\text{e}}$ |  Eq. (3) |  $2.73\times 10^{-14}\,$\mathrm{k}\mathrm{p}\mathrm{c}$$
$\lambda_{\text{mfp,e}}$ |  Eq. (2) |  $5.83\,$\mathrm{k}\mathrm{p}\mathrm{c}$$
$l_{\chi}$ |  Eq. (10)‡ |  $320\,$\mathrm{k}\mathrm{p}\mathrm{c}$$
$H$ |  $\lvert d\ln p/dR\rvert^{-1}$ |  $300-600\,$\mathrm{k}\mathrm{p}\mathrm{c}$$
$c_{\mathrm{s}}$ |  $(\gamma T/\mu m_{\mathrm{p}})^{1/2}$ |  $890\,\mu^{-1/2}\,$\mathrm{k}\mathrm{m}\,\mathrm{s}^{-1}$$
$\varv_{\mathrm{cond}}$ |  Eq. (22)‡ |  $510\,$\mathrm{k}\mathrm{m}\,\mathrm{s}^{-1}$$
$\varv_{\mathrm{A}}$ |  $B/(4\pi\rho)^{1/2}$ |  $140\,$\mathrm{k}\mathrm{m}\,\mathrm{s}^{-1}$$
$\beta_{\mathrm{e}}$ |  $8\pi n_{\mathrm{e}}T_{\mathrm{e}}/B^{2}$ |  $50$
$\mathrm{Kn}$ |  $\lambda_{\text{mfp,e}}/H$ |  $0.01-0.02$
$\alpha$ |  Eq. (19) |  $0.01-0.05$
$\tilde{\beta}$ |  Eq. (20) |  $2-24$
$\sigma$ |  Eq. (21) |  $\lesssim 1$ (in MTI turbulence)

##  3 Methods

In this section we describe the numerical methods that we have used, alongside
the diagnostics, the initial conditions and the parameters of our simulations.

###  3.1 SNOOPY code

The simulations presented in this work have been performed with a modified
version of SNOOPY, a pseudo-spectral three-dimensional (3D) MHD code (Lesur,
2015), where we included anisotropic heat conduction and implemented our
closure for the whistler-suppressed heat diffusivity (Eq. 18). In the SNOOPY
code, a Fourier transform is applied to the MHD equations, and the Fourier
coefficients of each wavenumber
$\bm{k}=k_{x}\bm{e}_{x}+k_{y}\bm{e}_{y}+k_{z}\bm{e}_{z}$ are then advanced in
time with a three-step Runge-Kutta algorithm. In addition to employing a
standard $2/3$ de-aliasing rule to all the physical variables, we decide to
filter the whistler-suppressed diffusivity through a two-point harmonic
averaging and to apply a hyperdiffusion operator to the buoyancy equation to
avoid Gibbs oscillations at the grid-scale, see Eq. (17). For further details
on the numerical implementation, see Appendix C.

####  3.1.1 Diagnostics

We use a variety of diagnostics to explore the impact of whistler suppression
on the MTI, both in real and in spectral space. We track the time evolution of
the kinetic ($K$), magnetic ($M$) and potential ($U$) specific energies, as
follows:

$\displaystyle K=\frac{1}{2}\langle u^{2}\rangle,\hskip
10.00002ptM=\frac{1}{2}\langle B^{2}\rangle,\hskip
10.00002ptU=\frac{1}{2N^{2}}\langle\theta^{2}\rangle,$ (25)

where the angle brackets denote the volume average (we remind the reader that
here and in the next sections the magnetic field is rescaled by
$\sqrt{4\pi\rho_{0}}$ so that it has units of a velocity). We similarly track
the mean suppression of thermal diffusivity, $\langle\chi/\chi_{0}\rangle$. In
addition to these volume-averaged quantities, we also follow the time-
evolution of the energy contained at each Fourier wavenumber
$k=\lvert\bm{k}\rvert$ (the so-called “shell-integrated” power spectra), as
well as the rates of energy injection and removal at each given $k$ by all the
various physical processes (e.g. fluid advection, magnetic tension,
dissipative processes, and so forth) corresponding to the terms in Eqs.
(15)-(17). For instance, the evolution of the potential energy contained at
each wavenumber is described by

$\displaystyle\frac{d}{dt}\left(\frac{1}{2N^{2}}\lvert\hat{\theta}(\bm{k},t)\rvert^{2}\right)=\frac{1}{N^{2}}\sum_{i}\Re\left\\{\hat{\theta}^{*}(\bm{k},t)\mathcal{F}\left[w_{i}\right]\right\\}$
(26)

where $\mathcal{F}$ is the Fourier transform operator (also denoted by a hat
$\hat{\leavevmode\nobreak\ }$), the asterisk denotes the complex conjugate,
and $w_{i}$ are the terms corresponding to the fluid advection
($-\bm{u}\bm{\cdot}\bm{\nabla}\theta$), the buoyancy force ($N^{2}u_{z}$), the
anisotropic dissipation
($\bm{\nabla}\bm{\cdot}\left[\chi\bm{b}\left(\bm{b}\bm{\cdot}\bm{\nabla}\right)\theta\right]$),
the energy injection rate from the background temperature gradient due to the
MTI
($\omega_{\mathrm{T}}^{2}\bm{\nabla}\bm{\cdot}\left(\chi\bm{b}b_{z}\right)$),
and the hyperdissipation ($\chi_{\text{hyp}}\bm{\nabla}^{6}\theta$). The
shell-averaged potential spectral energy density $E_{U}(k)$ is then computed
as

$E_{U}(k,t)dk=\sum_{k\leq\lvert\bm{k}^{\prime}\rvert<k+dk}\frac{1}{2N^{2}}\lvert\hat{\theta}(\bm{k}^{\prime},t)\rvert^{2},$
(27)

and similarly for the kinetic and magnetic spectral energy densities.

###  3.2 Initial conditions and parameters

####  3.2.1 2D simulations

We perform four two-dimensional (2D) runs with increasingly stronger whistler
suppression ($\alpha=4\times 10^{-5},4\times 10^{-4},2\times 10^{-3},4\times
10^{-3}$) and compare them to a reference unsuppressed ($\alpha=0$) MTI run.
The reason why we first focus on 2D runs is that despite the different
saturation mechanism (which in 2D involves an inverse energy cascade to large
scales) some key features of MTI turbulence are similar between 2D and 3D (;
). One notable exception is the absence of an SSD in 2D (Zeldovich, 1957),
which will be shown to play a crucial role in the whistler suppression of MTI
turbulence. For this reason we believe that 2D simulations constitute a useful
benchmark which helps us investigate and isolate the main features of whistler
suppression and the impact of the SSD.

Note that in the 2D runs, the values of the suppression parameter $\alpha$ are
lower than what was estimated in Sect. 2.5. This is to obviate the fact that
in 2D without an SSD, the magnetic energy at saturation (which affects the
level of whistler suppression according to Eq. 18) is not pinned at
equipartition with the kinetic energy, but depends on the initial strength of
the imposed (weak) uniform magnetic field. Boosting the initial uniform
$\bm{B}$ field to match the equipartition strength is also problematic, since
in 2D strong fields lead to formation of structures at the box size (PL22b).
To get around this issue, for our 2D simulations only we decide to compensate
the unphysically larger value of $\tilde{\beta}$ with a correspondingly lower
$\alpha$, such that their product (and therefore the overall suppression of
thermal conduction, see Eq. 18) roughly coincides with the levels estimated in
Sect. 2.5. We believe that this strategy allows us to obtain a more realistic
picture of whistler suppression on the MTI, and note that this issue will be
fully overcome with our 3D simulations of Sect. 3.2.2.

All runs have a resolution of $1024^{2}$ and aspect ratio unity, thermal
diffusivity $\chi_{0}/(L^{2}\omega_{\mathrm{T}})=10^{-3}$ with a thermal
Prandtl number $\text{Pr}=\nu/\chi=0.01$, and a magnetic Prandtl number
$\text{Pm}=\nu/\eta=1$. The hyperdiffusivity is $\chi_{\mathrm{hyp}}/(\Delta
h^{4}\chi_{0})=5.2\times 10^{-3}$, where $\Delta h=L/1024$ the grid size. The
runs are initialized with a mean uniform horizontal magnetic field with
strength $B_{0}/(L\omega_{\mathrm{T}})=10^{-4}$, and random velocity
fluctuations of amplitude $|\delta u_{i}|/(L\omega_{\mathrm{T}})=10^{-4}$, for
$i=x,z$. The magnetic field lines are initially isothermal ($\sigma=0$) and
the thermal diffusivity does therefore not feel any whistler suppression at
time $t=0$.

Figure 1: Time evolution of our MTI runs with suppression of thermal
diffusivity in 2D (left column), and in 3D with whistler and flat suppression
(middle and right column, respectively). From top to bottom: volume-averaged
kinetic, magnetic energy and average suppression factor as a function of time.
In 3D with whistler suppression, MTI turbulence decays if the suppression
parameter $\alpha$ is large enough, while in 2D (where there can be no SSD),
or in 3D with flat uniform suppression (with or without an SSD) MTI turbulence
never dies completely. Left panel: the cyan vertical line at $t=150$ indicates
the time of the snapshots shown in Fig. 4. The inset in the top panel zooms in
on the simulation with $\alpha=4\times 10^{-3}$ (gold solid line) at late
times, and shows periodic oscillations of frequency $\omega=N/\sqrt{5}$
associated with a dominant large-scale $g$ mode $(k_{x},k_{z})=(1,2)\times
2\pi/L$. Middle and right panels: the cyan vertical line at $t=186$ indicates
the time when we turn on suppression of thermal diffusivity in our reference
run.

####  3.2.2 3D simulations

In addition to our 2D runs, we perform a number of 3D simulations with
whistler suppression. This allows us to examine a more realistic regime where
the magnetic field strength is kept self-consistently at equipartition by the
SSD, and feeds back into the suppression of thermal conductivity.

We perform six 3D runs with whistler suppression and increasing $\alpha$
between $[4\times 10^{-3},0.4]$, and compare the results to a reference 3D run
without suppression ($\alpha=0$). The values of $\alpha$ are chosen so as to
sample the range that can be realistically expected in the periphery of galaxy
clusters (see Sect. 2.5). All our 3D runs have a resolution of $288^{3}$,
aspect ratio of unity, and $\chi_{0}/(L^{2}\omega_{\mathrm{T}})=0.0025$. To
reflect the hierarchy between thermal, viscous and resistive diffusivities in
the ICM, we choose explicit viscosity and resistivity such that
$\text{Pr}\simeq 0.042$ and $\text{Pm}=4$. The hyperdiffusivity is
$\chi_{\mathrm{hyp}}/(\Delta h^{4}\chi_{0})=4.8\times 10^{-3}$.

In order to save computational time, we simulate the linear phase of the
instability and the subsequent kinematic and nonlinear dynamo phase only once
in our reference run. Then, once it reaches saturation, we choose one snapshot
and use it as a starting point for our series of runs with whistler
suppression. The reference run is initialized with random velocity
fluctuations of amplitude $|\delta u|/(L\omega_{\mathrm{T}})=10^{-4}$ in all
directions. Contrary to our 2D runs, we do not include a mean uniform
horizontal magnetic field, but start with a helicoidal field varying with $z$
pointing in the $x-y$ plane: $\bm{B}=B_{0}\left[\sin(2\pi
z/L)\bm{e}_{x}+\cos(2\pi z/L)\bm{e}_{y}\right]$, with
$B_{0}/(L\omega_{\mathrm{T}})=10^{-4}$.

We stress that the lack of a net magnetic flux in our 3D simulations means
that – if not replenished by an SSD – magnetic fields are now allowed to
resistively decay to negligible values. This introduces a new equilibrium
state available to the system where magnetic fields are zero, thermal
conduction completely suppressed and the MTI rendered ineffective (a “dead”
state). Assessing whether and under what conditions the system reaches this
dead state will be one of the main topics of the following sections. We decide
not to run 3D simulations with a (weak) net magnetic flux because, if they are
able to sustain an SSD, we expect that their final state at saturation will
show little difference from an equivalent no-net flux run (as shown in PL22b);
if they cannot sustain an SSD, on the other hand, then the residual thermal
diffusivity after the turbulent magnetic fields decay will depend on the
initial arbitrary weak flux, and will not be representative of the physical
conditions in the periphery of galaxy clusters.

###  3.3 Simulations with external forcing

To study the impact of additional sources of turbulence to the MTI with
whistler suppression, we run a number of MTI simulations where we add an
external forcing to the momentum equation. In this paper we show the results
of two MTI runs with external forcing, one without whistler suppression
($\alpha=0$), and another with strong whistler suppression ($\alpha=0.12$).
For comparison, we also run a forcing-only run (without the MTI) with
isotropic thermal conductivity equal to $\chi_{0}$ and no whistler
suppression, while keeping all other parameters unchanged.

These simulations are the same as in PBP23, and are initialized from a
snapshot at late times of the corresponding 3D run without external forcing
introduced in Sect. 3.2.2. We report here the numerical details for
completeness. We implement a delta-correlated in time, solenoidal, white-noise
forcing in Fourier space by drawing six random numbers at each timestep for
every wavevector $\bm{k}$ in a specified range (we choose
$|\bm{k}|=\\{4,5,6\\}\times 2\pi/L$): three random numbers determine an
auxiliary random vector used to enforce the divergence-free condition
($\bm{k}\bm{\cdot}\hat{\bm{f}}_{\bm{k}}=0$), and three multiply the complex
phases of the forcing vector $\hat{\bm{f}}_{\bm{k}}$. The amplitude of
$\hat{\bm{f}}_{\bm{k}}$, instead, is constant and the same for all the forced
wavevectors. The resulting forcing has zero-mean helicity but allows for small
fluctuations around zero.

###  3.4 Runs with flat uniform suppression

In addition to whistler-suppressed thermal diffusivity, we also run a number
of simulations with a flat uniform suppression below the Spitzer value:

$\displaystyle\chi=f\chi_{0},$ (28)

where $f<1$ is a constant. Comparing the two models will allow us to clearly
identify the features of whistler suppression and what physics might be missed
by employing a simple flat uniform suppression.

##  4 Weakening of the MTI and critical suppression

In this section, we show that both in 2D and in 3D simulations the MTI
turbulence is weakened as the suppression parameter $\alpha$ increases (Sect.
4.1). However, while in 2D we observe a gradual decrease in the turbulent
levels as we sweep through the range of $\alpha$, in 3D there is a sharp
transition at a critical value of the suppression parameter, beyond which the
MTI turbulence dies (Sect. 4.2). We relate this behaviour to the impairment of
the SSD, and we propose a simple model to predict the critical $\alpha$ (Sect.
4.3).

###  4.1 Time evolution of the MTI with whistler suppression

In Fig. 1 we show the time-series of the volume-averaged kinetic, potential
and magnetic energy and the average suppression factor of our whistler-
suppressed runs in 2D (left column), 3D with whistler suppression (middle
column), and 3D with uniform suppression (right column).

In 2D (where we evolve the MTI with whistler suppression throughout the linear
phase of the instability), we observe that simulations with larger $\alpha$
show a shorter phase of exponential growth, followed by a secondary growth
phase until the system finally saturates. We also note that the average
suppression factor decreases sharply during the linear phase of the
instability before eventually stabilizing at somewhat larger values. This is
due to the fact that as the instability grows the magnetic field line is
brought further away from isothermality ($\sigma$ increases). However, as the
MTI approaches the nonlinear phase and the vertical magnetic perturbation
grows stronger, the total magnetic energy increases (see second row of Fig. 1,
left column), $\tilde{\beta}$ decreases and the whistler suppression becomes
less effective. For further details concerning the nonlinear saturation of the
MTI with whistler suppression, see Appendix B.

After the instability has saturated, we observe that in both 2D and 3D with
whistler suppression and with the lowest $\alpha$ there is little difference
compared to the reference MTI run, and the final energies are similar. With
higher suppression parameter the behaviour between 2D and 3D drastically
differs: while in 2D the final saturated energies progressively decrease as
$\alpha$ increases, in 3D this gradual trend holds only for the three runs
with the lowest suppression parameter ($\alpha=4\times 10^{-3},0.012,0.024$),
beyond which the system undergoes a sharp transition and the turbulence dies
out. This complete switch-off of turbulence is particularly evident in the run
with $\alpha=0.04$ (dark orange), which experiences a slow decline over
hundreds of dynamical times before suddenly decaying. This behaviour is
visible across all the volume-averaged quantities tracked in Fig. 1 and
further increasing $\alpha$ only accelerates the decay. This is in contrast to
simulations with a model of flat uniform suppression of thermal conductivity
(third column of Fig. 1), where irrespective of whether the turbulence is
sufficient to drive an SSD ($\chi/\chi_{0}=0.40,0.68$) or not
($\chi/\chi_{0}=0.33$), the MTI turbulence never dies completely but is only
reduced.

We interpret the radically different behaviour with small changes to the
$\alpha$ parameter as a sign of a critical transition in whistler-suppressed
MTI turbulence. This is a feature of many nonlinear systems and we speculate
that two dynamical fixed points exist for low values of $\alpha$,
corresponding to the dead state and the MTI-turbulent state. As $\alpha$ is
increased above a critical value the MTI turbulent state ceases to exist, and
the system moves onto the dead state. A more in-depth analysis of this
critical transition will be carried out in Sect. 4.3. This picture may change
substantially in simulations with a (weak) net magnetic flux, which would set
a lower limit for the thermal conductivity, thus removing the dead state. We
speculate that, depending on the strength of the initial uniform magnetic
field, the simulations may either become effectively dead (if the residual
thermal diffusivity is too low for the MTI to grow), or possibly exhibit a
cyclic behaviour of MTI growth and decay.

Before proceeding with the analysis of the turbulence at saturation, we note
that in the 2D simulations with largest $\alpha$ the system starts exhibiting
an oscillatory pattern at late times, most clearly visible in the kinetic
energy of Fig. 1, left column. We explain this behaviour with the emergence of
a large-scale $g$ mode on top of the smaller-scale MTI turbulence, which we
are able to track over a few cycles in the inset in Fig. 1 (top-left panel).
In this particular case the $g$ mode has wavevector $(k_{x},k_{z})=(1,2)\times
2\pi/L$, but in different runs we detected other orientations. The appearance
of large-scale $g$ modes in MTI simulations is in itself not surprising, and
indeed in 2D the instability saturates by exciting large-scale $g$ modes that
block the inverse cascade of kinetic energy and redistribute it in the form of
density fluctuations (PL22a). The difference here is that with strong whistler
suppression $g$ modes end up dominating energetically over the smaller-scale
turbulence. Interestingly, we observe the emergence of large-scales $g$ modes
also in the decaying 3D runs, possibly excited by the weak residual MTI energy
injection, meaning that the dead solution is not completely dead. However, we
believe this phenomenon to be only transient given the ongoing decrease of the
thermal diffusivity which further quenches MTI injection.

Figure 2: Time-averaged kinetic energy versus mean suppression factor for our
2D simulations (left panel) and 3D simulations (right panel). The points are
color-coded as in Fig. 1. In 2D the time-averages have been performed over a
window of $t=[130,170]$ to reduce contamination by the large-scale $g$ modes.
In 3D we show time-averages for the three surviving runs only, while the three
decaying runs are shown as solid lines (the full extent of the decay is shown
in the inset). In both panels the black solid line represents a relation of
linear proportionality. The linear proportionality between between
conductivity and turbulence levels found in PL22a and PL22b is seen to remain
valid in an average sense in simulations including whistler-suppression of
thermal conduction.

###  4.2 Correlation of kinetic energy with thermal diffusivity

We now take a closer look at how the suppression of MTI turbulence relates to
the volume-averaged $\chi/\chi_{0}$. In Fig. 2 we plot the time-averaged
kinetic energy versus the suppression factor for the 2D runs (left panel) and
the 3D runs (right panel). Note that in 3D we compute the time-averages
(denoted by $\langle...\rangle_{t}$) only for the three surviving runs
($\alpha=4\times 10^{-3},0.012,0.024$) and we instead show the correlation of
kinetic energy with thermal diffusivity as a function of time for the three
runs that die ($\alpha=0.04,0.12,0.4$). The different datapoints (each
corresponding to a value of $\alpha$) closely follow a linear trend with the
averaged $\chi/\chi_{0}$ both in 2D and in 3D. Interestingly, in 3D this trend
applies not only to the runs that adjust to the new turbulent saturated state,
but also to the runs that die as they decay. In other words, despite the
strong spatial inhomogeneity of the effective thermal diffusivity
$\chi/\chi_{0}$ (which will be the subject of Sect. 5), in our simulations the
turbulent kinetic energy appears to scale proportionally with the volume-
averaged suppression factor, extending the basic scaling relationship of MTI
turbulence between $K$ and $\chi$ in the whistler-suppressed regime. This
finding lends weight to the idea that, at least as a first approximation, the
impact of whistler suppression of thermal conductivity on the MTI can be
described by a simple model with a flat uniform suppression. However, as we
show in the next Section, a flat uniform suppression cannot explain the
appearance of a critical transition in 3D MTI turbulence that we witnessed in
Fig. 1.

Figure 3: Top panel: time-averaged turbulent kinetic energies, Reynolds number
at the integral scale, and average suppression factor for our 3D runs as a
function of $\alpha$. The values are normalized to those of the unsuppressed
run. The leftmost point, corresponding to $\alpha=0$, has been placed at
$\alpha=2\times 10^{-3}$ for visualization purposes. The three right-most data
points represent decaying runs, thus we do not show error bars. The light-gray
shaded region indicates the approximate range of the $\alpha$ parameter for
typical ICM parameters (see Table 1), while the dark-gray shaded region
represents the predicted critical $\alpha$ where transition from MTI
turbulence to dead state occurs. Bottom panel: model for critical $\alpha$
with whistler suppression (see text). The blue band indicates the
equipartition value of $\tilde{\beta}$ with the MTI turbulence (Eq. 36). The
red band divides the $\alpha-\tilde{\beta}$ space in two regions: where an SSD
is possible (below the curve), and where the SSD is not sustainable (hatched
region), see Eq. (39). The two bands cross approximately at $\alpha\simeq
0.038-0.046$ above which MTI turbulence dies. The gold symbols show the
saturated values of $\tilde{\beta}$ for the non-decaying runs.

###  4.3 Impairment/shutdown of the SSD in 3D and the “death-spiral”

In this section we propose a simple toy model to explain the critical
transition in 3D between weakened MTI turbulence with whistler suppression and
the dead state. Our model rests on the assumption that at saturation MTI-
driven turbulence pins kinetic and magnetic turbulent fluctuations
approximately at equipartition

(30)

provided that the criterion for an SSD (magnetic Reynolds number
$\mathrm{Rm}\equiv\mathrm{Pm}\times\mathrm{Re}_{\mathrm{i}}$ larger than a
critical value, PL22b, ) is fulfilled:

(32)

More precisely, our hypothesis is that the critical transition emerges when
these two conditions cease to be satisfied: if magnetic fields cannot be kept
at equipartition levels by an SSD, their decay leads to further suppression of
thermal conductivity through its dependence on the plasma $\beta$. Stronger
whistler suppression results in weaker MTI turbulence, and an even faster
decay of magnetic energy, producing a self-reinforcing cycle which leads to
the dead state and which we call the “death-spiral”.

The previous argument can be made more quantitative using the two conditions
(Eqs. 30-32) combined with the MTI scaling laws found in ;  and summarized in
Sect. 2.3, allowing us to estimate the critical value of $\alpha$. In the
following, we make the simplifying assumption that the effective volume-
averaged suppression factor $\langle\chi/\chi_{0}\rangle$ can be approximately
expressed in terms of appropriate volume averages of $\tilde{\beta}$ and
$\sigma$ as

$\displaystyle\left\langle\dfrac{\chi}{\chi_{0}}\right\rangle\approx\dfrac{1}{1+\dfrac{\phi}{3}\alpha\langle\tilde{\beta}\rangle\langle\sigma\rangle},$
(33)

where $\phi$ is a parameter that takes into account the nonlinear dependence
of $\langle\chi/\chi_{0}\rangle$ on $\tilde{\beta}$ and $\sigma$, as well as
their statistical correlation (see Sect. 5.1). We test Eq. (33) comparing the
predicted value of $\langle\chi/\chi_{0}\rangle$ with the actual simulation
data from the three surviving 3D runs, and find that the two sets of values
match reasonably well for $\phi\simeq 3.0-3.7$ (the average relative error is
$\sim 3\%$). With this caveat in mind, we drop the angle-bracket notation for
the rest of the section and assume that we are always referring to volume-
averaged quantities.

Using the MTI scaling law in Eq. (7), the equipartition condition (Eq. 30) can
be rewritten in terms of the specific magnetic energy as

$\displaystyle\frac{M}{\varv_{\text{cond}}^{2}}=\text{const.}\times\left\langle\dfrac{\chi}{\chi_{0}}\right\rangle\dfrac{\omega_{\mathrm{T}}^{2}}{N^{2}},$
(34)

where the constant of proportionality can be read directly from the value of
$M/\varv_{\text{cond}}^{2}$ at saturation for our reference 3D unsuppressed
run in Fig. 1, for which we have $M/\varv_{\text{cond}}^{2}\simeq 0.043$, and
thus $\text{const.}\simeq 0.043N^{2}/\omega_{T}^{2}=0.0043$. Replacing
$\langle\chi/\chi_{0}\rangle$ in the above equation with Eq. (33), and solving
for $\tilde{\beta}$ the equipartition condition can be rewritten as

(36)

where we used the equivalence $1/\tilde{\beta}\equiv
M/\varv_{\mathrm{cond}}^{2}$. The physical meaning of Eq. (36) is that, given
a certain value of $\alpha$ (and of effective suppression of thermal
diffusivity), the value of $\tilde{\beta}$ changes accordingly to maintain
equipartition with the turbulent kinetic energy.

A second independent constraint on $\tilde{\beta}$ is provided by the
criterion for the SSD which sets a lower limit on the magnetic Reynolds number
below which the dynamo cannot be sustained. This constraint originates from
the fact that in whistler-suppressed MTI turbulence the Reynolds number at the
integral scale $\text{Re}_{\mathrm{i}}$ also scales approximately linearly
with the average suppression factor (Fig. 3, top panel):

$\displaystyle\text{Re}_{\mathrm{i}}\simeq\text{Re}_{\mathrm{i}0}\times\left\langle\dfrac{\chi}{\chi_{0}}\right\rangle,$
(37)

where $\text{Re}_{\mathrm{i}0}$ denotes the Reynolds numbers of the reference
unsuppressed run. Solving for $\tilde{\beta}$ and using Eq. (37), the
criterion for the SSD in Eq. (32) can be recast as follows:

(39)

where $\mathrm{Rm}_{0}\equiv\mathrm{Pm}\mathrm{Re}_{i0}$ is the unsuppressed
magnetic Reynolds number of the plasma. The physical meaning of Eq. (39) is
that in MTI turbulence with whistler suppression, only magnetic fields of
strength above a certain $\alpha$-dependent threshold (equivalently, low
enough plasma $\tilde{\beta}$) can be sustained self-consistently by the MTI-
driven SSD.

If we plot the two conditions in Eqs. (36) and (39) (lower panel of Fig. 3,
assuming a value of $\phi$ between $3.0-3.7$), we see that for very low values
of $\alpha$, the equipartition $\tilde{\beta}$ (solid blue line) is well below
the SSD threshold (solid red curve), and thus equipartition-level magnetic
fields can be sustained by the SSD in presence of moderate suppression, as we
observe in our 3D runs (golden hexagons in Fig. 3). As $\alpha$ is increased,
however, the equipartition $\tilde{\beta}$ progressively nears the SSD
threshold, and eventually crosses it: at this point, equipartition-level
magnetic fields become too weak to sustain themselves, and the system enters
the death-spiral. From our model, we estimate this critical parameter
$\alpha_{\mathrm{c}}$ to be of the order of $\alpha_{\mathrm{c}}\simeq
0.038-0.046$ (depending on the value of the parameter $\phi$)222In drawing the
lower panel of Fig. 3 we have assumed a volume-averaged value of $\sigma=0.5$
which is constant with $\alpha$. This is in reasonable agreement with the
results of our three surviving 3D runs (where the average $\sigma$ is
$\left[0.52,0.48,0.46\right]$ for increasing $\alpha$, respectively); using
instead $\sigma=0.4$ the corresponding critical value is
$\alpha_{\mathrm{c}}\simeq 0.048-0.059$. , which is consistent with the
findings from our numerical simulations (where
$\alpha_{\mathrm{c}}\in\left(\left.0.024,0.04\right]\right.$), and lies at the
upper end of the values of $\alpha$ that can realistically be expected in the
periphery of galaxy clusters (see Sect. 2.5). Note that the critical value of
$\alpha$ for the impairment of the SSD is a function of $\mathrm{Rm}_{0}$: the
higher the magnetic Reynolds number, the higher the allowed values of
$\tilde{\beta}$ that can be self-sustained by a magnetic dynamo (see Eq. 39).
In the ICM we expect $\mathrm{Rm}_{0}$ to be very large, which suggests that
even with significant whistler suppression, the MTI could in theory still be
able to drive an SSD saturating at low-amplitudes. In this case, however,
other sources of turbulence would likely take over the MTI, a scenario that we
explore more in detail in Sect. 6.

###  4.4 Differences with flat uniform suppression

We now summarize our findings on the MTI with whistler suppression and compare
to a simple model of uniform suppression of thermal conductivity. As we saw in
Sect. 4.2, both models show direct relationship between the properties of the
turbulence and the level of suppression of thermal diffusivity, whether that
is a fixed ratio uniform in space, or a volume-averaged quantity. In this
respect, the two models do not differ significantly, and as a result – given
an average estimate of whistler suppression below the Spitzer value in the ICM
– it is in principle possible to determine the resulting level of MTI
turbulence. However, care must be taken when extrapolating to the cluster’s
periphery (or when running simulations with the simpler model of flat uniform
suppression) because, as we have seen in Sect. 4.3, for values of the
suppression parameter $\alpha$ above a certain critical threshold the SSD
shuts down and the MTI decays. Note that this death spiral is fundamentally a
nonlinear feedback intrinsic to our model of whistler suppression and would
not take place with flat uniform suppression. We see this difference most
clearly in two aspects: first of all, although with large enough flat
suppression the SSD also shuts down, this has no impact on the level of MTI
turbulence, which continues unperturbed (PL22b); second, from Fig. 1 we saw
that a simulation with flat suppression to $40\%\chi_{0}$ can still sustain an
SSD and keep magnetic fields at equipartition, while a whistler-suppressed run
with initial volume-averaged suppression of $\simeq 50\%\chi_{0}$ (as for the
run with $\alpha=0.12$) does not.

##  5 Modifications to the flow with whistler suppression

In our discussion so far we have mostly concerned ourselves with volume-
averaged quantities, such as the correlation between the turbulent kinetic
energy and the volume-averaged suppression of thermal diffusivity. However,
our model of whistler suppression is a nonlinear function of quantities (such
as the modified plasma $\tilde{\beta}$) that can strongly vary in space and
across scales, and therefore we expect this variation to be reflected in the
suppression of thermal conductivity as well. This will be shown in Sect. 5.1.
We then take a closer look at the spatial correlations of MTI turbulence and
we look at how whistler suppression affects the energy distribution and the
MTI energy transfers in spectral space. We find that high-$\beta$ regions are
progressively brought out of thermal equilibrium as $\alpha$ increases (Sect.
5.2). In Sect. 5.3 instead we observe that – while whistler suppression
reduces the efficiency of MTI turbulent driving – this effect appears to be
mostly limited to the largest turbulent scales.

Figure 4: Multiplot of the 2D runs at saturation ordered by increasing
suppression of thermal diffusivity (from left to right). From top to bottom:
magnitude of the velocity field, density fluctuation, Boussinesq plasma beta,
Boussinesq parallel temperature length, point-wise thermal diffusion
suppression fraction. Note that the colour bar is in common for each row.
Higher suppression parameter, $\alpha$, leads to weaker MTI turbulence, weaker
magnetic fields (higher $\tilde{\beta}$), less isothermal fieldlines (higher
$\sigma$) and lower average thermal conductivity $\chi$.

###  5.1 Morphology of the flow and spatial correlations

To visualize how the MTI turbulence reacts to increasingly stronger whistler
suppression, in Fig. 4 we show snapshots of the 2D runs after they have
reached saturation (time $t=150$). We focus here on 2D simulations because
they allow us to explore a wider range of suppression parameter $\alpha$ than
in 3D, where the differences between the surviving runs are less pronounced.
From a visual analysis of the turbulent velocity and temperature field (first
and second row) we can immediately see that higher $\alpha$ lead to a
progressive weakening of the MTI turbulence. This is in spite of the fact that
the suppression of thermal diffusion (last row) is not uniform across the
domain but is highly inhomogeneous and tends to follow the morphology of the
magnetic field (third row). Note, however, that even in the two runs with
strongest whistler suppression there are still regions with thermal
diffusivity higher than $50\%$ of its reference value, amounting to 11% and 4%
of the whole domain, respectively. While the isothermality parameter $\sigma$
undergoes a rather dramatic change as $\alpha$ is increased, its role in
determining the amount of whistler suppression is less clear. Finally, from
the velocity amplitude and the $\tilde{\beta}$ field (first and third rows),
one can see that the largest scale of the turbulence (the integral scale) also
decreases with $\alpha$. We observe a similar behaviour in our 3D surviving
runs.

To better understand the interplay between $\tilde{\beta}$ and $\sigma$, and
to assess to what extent the two quantities are spatially correlated or anti-
correlated, we employ a variety of diagnostics, including level plots, as well
as one- and two-dimensional probability distributions.

In Fig. 5 we show a zoomed-in view of the suppressed thermal diffusivity of
one 2D run ($\alpha=2\times 10^{-3}$, from Fig. 4) where we separately show in
the inset the regions with “low” ($\chi/\chi_{0}<0.05$, left column), and
“high” diffusivity ($\chi/\chi_{0}>0.2$, right column) and plot the values of
$\tilde{\beta}$ and $\sigma$. As we can see, the morphology of the low- and
high-diffusivity regions is rather different: while low-$\chi$ regions appear
as “holes”, areas with high-$\chi$ mostly look filamentary. A similar picture
can be observed in 3D as well. We relate this behaviour to the intermittency
of magnetic fields in MHD turbulence, with magnetic field lines (shown in the
second row of the inset) concentrated in thin bundles near regions of strong
magnetic fields. As a result, with whistler suppression thermal conductivity
remains high in these thin bundles while it drops significantly in between.
Moreover, from Fig. 5 we see an interesting interplay emerge between the
modified plasma $\tilde{\beta}$ and the isothermality parameter $\sigma$, as
these quantities appear to be positively correlated: where $\tilde{\beta}$ is
high, also $\sigma$ is consistently above $1$, which indicates that the plasma
is locally out of isothermality; conversely, in regions with strong magnetic
field (low $\tilde{\beta}$) the isothermality parameter $\sigma$ is typically
less than unity, and the plasma close to isothermality. This positive
correlation, however, is the result of two combined effects, which highlights
the nonlinearity of the model: in fact, with whistler suppression the thermal
diffusivity first drops in regions with high $\tilde{\beta}$, leading to an
increase in $\sigma$ as the plasma cannot thermalize as effectively as before.
This increase in $\sigma$, however, further reduces the thermal diffusivity
(Eq. 18), and will therefore lead the plasma farther away from isothermality.
These effects end up reinforcing each other, and produce the tight
correlations shown in Fig. 5.

Figure 5: Zoom-in view of the snapshot of the 2D run with $\alpha=2\times
10^{-3}$ at time $t=150$ showing the spatial correlation between
$\tilde{\beta}$ and $\sigma$. In the inset we separately show regions with low
(left column) and high thermal diffusivity (right column) by cutting out
regions where $\chi/\chi_{0}$ is above $0.05$ and below $0.2$, respectively.
In the second row of the inset, left column, the magnetic field lines are
shown in blue. With whistler suppression, regions of high thermal conductivity
are shaped as thin bundles near filaments of strong magnetic fields. Figure 6:
Probability distribution as a function of $\chi/\chi_{0}$ and $\tilde{\beta}$
in our 3D runs with whistler suppression at saturation. The colour bar
represents the probability density of finding a value in the infinitesimal
interval of size $\mathrm{d}\log(\chi/\chi_{0})\,\mathrm{d}\log\tilde{\beta}$.
The two dash-dotted lines represent the two limits of Spitzer (in grey) and
whistler-suppressed thermal diffusivity (in blue). By increasing the
suppression parameter, $\alpha$, larger parts of the domain enter the
$1/\tilde{\beta}$ whistler-suppressed regime.  Figure 7: Distribution of the
data in the $\tilde{\beta}$–$\sigma$ plane for our 3D runs at saturation. The
colour bar represents the probability density of finding a value in the
infinitesimal interval of size
$\mathrm{d}\log\tilde{\beta}\;\mathrm{d}\log\sigma$. The cyan and gold curves
are isocontours of the probability distribution function that enclose $33\%$
and $67\%$ of the box volume, respectively. The solid gray contour lines in
the plots with whistler suppression show the amount of suppression of the
thermal diffusivity (for any given line, left of it are regions with
diffusivity larger than the number quoted, in units of $\chi_{0}$). Whistler-
suppression leads to the formation of a high-$\beta$ lobe consisting of points
that are not isothermal along field lines. This feature is not captured in the
reference simulations with flat or no suppression.

###  5.2 Probability distributions

To complement our visual analysis of spatial correlations, we compute a
variety probability distribution functions for our 3D runs with whistler
suppression.

It is instructive to first look at the distribution of the thermal
diffusivities as a function of the modified plasma $\tilde{\beta}$ in the 3D
surviving runs (Fig. 6), where we can clearly see the transition between the
Spitzer (at low $\tilde{\beta}$) and the whistler-suppressed regime (at high
$\tilde{\beta}$) that is encoded in Eq. (18). The turnover value of
$\tilde{\beta}$, where the transition between the two regimes happens (the
“knee”), can be easily estimated by weighing the two terms in the denominator
of Eq. (18). This yields $\tilde{\beta}_{\mathrm{turnover}}\approx
3/(\alpha\sigma)$, which shows that it decreases as the suppression parameter
$\alpha$ increases. As a result, for higher $\alpha$ a larger proportion of
the domain enters the whistler-suppressed regime. The vertical spread in
$\chi/\chi_{0}$ at fixed $\tilde{\beta}$ is caused by the different values of
the isothermality parameter $\sigma$.

We then show in Fig. 7 the two-dimensional probability density of finding a
cell with given $(\tilde{\beta},\sigma)$ in our 3D surviving runs, while in
Fig. 8 we separately show the marginal distributions along $\chi,\beta$, and
$\sigma$ (also in 3D). Without whistler suppression (Fig. 7, left-most panel),
values are distributed in the $\tilde{\beta}$–$\sigma$ plane to resemble a
triangle, with the vast majority of the volume ($>95\%$) lying below the
$\sigma=1$ line, meaning that the MTI effectively enforces isothermality along
magnetic field lines. The values with the highest occurrence are distributed
in the oblong-shaped region defined by the isocontours around
$\tilde{\beta}\sim 10^{2}$ and $\sigma\in[10^{-1},10^{0}]$. Increasing
$\alpha$ we notice a progressive shift of the high-$\tilde{\beta}$ part of the
distribution towards larger values of $\sigma$ (a “lobe”), with a few data
points reaching as far as $\sigma=10$. However, even when $\alpha=0.024$ the
fraction of the volume far away from isothermality ($\sigma>1$) is no more
than $9\%$. Interestingly, if we compute the two-dimensional probability
density of the local turbulent velocities and suppression of thermal
diffusivity, we do not notice any positive correlations, contrary to what we
observe at the volume-average level.

The fact that only high-$\beta$ regions deviate from isothermality is because
it is precisely these regions that are impacted most by whistler suppression
when the $\alpha$ parameter is increased (see the gray $\chi/\chi_{0}$ contour
lines in Fig. 7). When thermal diffusivity is strongly suppressed locally,
thermalization along magnetic field lines stops being efficient and $\sigma$
increases. This behaviour is a direct consequence of whistler suppression, and
does not take place in simulations with a flat uniform suppression of thermal
diffusivity. To show this, we plot in Fig. 7 (right-most panel) the
probability distribution of a run with flat suppression where the (spatially-
uniform) suppression coefficient is chosen to be equal to the volume-averaged
suppression factor of the run with $\alpha=0.024$ (second-last panel). With
flat suppression, the distribution of values in the $\tilde{\beta}$–$\sigma$
plane is surprisingly similar to the reference run, except for an overall
shift towards larger $\tilde{\beta}$ in order to maintain equipartition with a
weaker turbulence, and we observe no lobe out of isothermality.

This picture is further supported by the marginal distributions shown in Fig.
8, which show the progressive emergence of fat tails for $\sigma>1$ in the
one-dimensional distribution (bottom left panel), while for $\sigma<1$ the
distribution is practically unchanged. At the same time, the
$\tilde{\beta}$-distribution undergoes an overall shift towards larger values
(weaker magnetic fields, top panel). With flat uniform suppression, instead,
the marginal distribution over $\sigma$ is identical to that of the reference
unsuppressed run, while the marginal distribution over $\tilde{\beta}$ is
shifted towards larger values. We thus distinguish between a “local” and
“global” response of the turbulence to whistler suppression, where the local
response involves the inhibition of isothermality in those regions where
$\tilde{\beta}$ is large, while the global response is the overall decrease in
turbulence strength and magnetic energy in the system.

Figure 8: Marginal distributions of the 3D runs with respect to the modified
plasma $\tilde{\beta}$ (top), the isothermality parameter $\sigma$ (bottom-
left), and the effective thermal diffusivity $\chi/\chi_{0}$ (bottom-right).
Stronger whistler suppression decreases the overall magnetic energy levels,
and brings part of the domain out of isothermality.

###  5.3 Dynamics in Spectral Space

(a)

(b)

Figure 9: Left: compensated spectral energy densities (i.e., multiplied by the
wave number $k$) of the kinetic (top) and magnetic fluctuations (bottom) in
our 3D runs at saturation (averaged over hundreds of dynamical times). The
$y$-axis is in linear scale to highlight the differences between the runs. For
comparison, we show in the insets the kinetic and magnetic (uncompensated)
power spectra in log-log scale. Right: MTI energy injection (top) and balance
between MTI injection and anisotropic thermal dissipation (bottom) in spectral
space in 3D simulations. In both panels the spectral transfers are multiplied
by the wave number $k$ to highlight the important features. Despite the
reduction in MTI energy injection due to increasingly strong whistler
suppression, the final balance between injection and dissipation remains
mostly unchanged over a wide range of scales below the conduction length
$l_{\chi}$.

We now look at the impact of whistler suppression on thermal conductivity
across scales, using spectral energy densities (power spectra) and transfer
terms in spectral space. In Fig. 9(a) we show the kinetic and magnetic power
spectra as a function of $k$, the wavenumber, for the reference MTI run and
the three surviving 3D runs with whistler suppression. To aid visual
inspection, we compensate the spectra multiplying by $k$ and use a linear
scale in the $y$-axis. For comparison, the same uncompensated spectral energy
densities in log-log scale can be found in the insets of Fig. 9(a). Looking at
the kinetic power spectrum (top panel), we find that increasing $\alpha$
mostly affects the energy contained at the largest scales, while the small
scales are hardly affected. This behaviour is reminiscent of the MTI without
whistler suppression, where decreasing the (spatially-uniform) thermal
diffusivity leads to a decrease of energy at the largest scales, with the
integral scale moving to larger $k$ (PL22b). Note, however, that with whistler
suppression this shift is hardly noticeable given that the runs shown only
have a limited amount of suppression of thermal conductivity
($\langle\chi/\chi_{0}\rangle\gtrsim 68\%$). Looking at the magnetic power
spectrum, on the other hand, we notice that the decrease in energy is more
evenly distributed across wavenumbers.

In Fig. 9(a) we also show the spectra of the run with flat uniform suppression
and $\chi/\chi_{0}=68\%$, and see that they align rather well with the
equivalent whistler-suppressed run with $\alpha=0.024$ (particularly so in the
kinetic power spectrum), suggesting a similar dynamics in spectral space. We
can get additional insight by looking at the transfers of energy in spectral
space. We focus on two key processes of the MTI, namely the injection of
energy from the background temperature gradient by the anisotropic
conductivity and the net injection plus anisotropic dissipation, which is a
more meaningful measure of how quickly energy is eventually added or removed
to the system. The behaviour of these two terms (shown in Fig. 9(b), top and
bottom panels, respectively) is particularly revealing: while they
individually show a clear downward trend as $\alpha$ is increased (meaning
that the energy injection and dissipation become less efficient both at large
and small scales), their sum is instead rather constant, except at the largest
scales comparable to the conduction length $l_{\chi}$. If we now look at the
run with flat uniform suppression and compare it with the run with
$\alpha=0.024$ (that has the same volume-averaged thermal diffusivity), we
notice that the energy injection efficiency looks clearly different: at
intermediate scales ($kL\lesssim 100$) MTI forcing is more efficient with
whistler suppression than with flat uniform suppression, while at small scales
the situation is reversed. At the largest scales the two are very similar.
Interestingly, these differences again mostly vanish once we look at the net
effect of injection plus dissipation.

###  5.4 Differences with flat uniform suppression

In the last sections we carried out a detailed study of the statistical
properties of MTI with whistler suppression, and looked at how the spectral
dynamics is affected. The comparison with a model of flat uniform suppression
helped us identify some key features of whistler suppression, which we
classified as “local” (the bringing of high-$\tilde{\beta}$ regions out of
isothermality) or “global” (the overall decrease of magnetic field strength).
Both of these outcomes influence the level of whistler suppression itself, and
therein lies the nonlinearity of the process. The differences in the spectral
space are more subtle: despite the clear change in the MTI forcing efficiency,
the net contribution of injection and dissipation due to anisotropic diffusion
shows striking similarities with simulations with flat suppression, which
helps explain why the kinetic and magnetic power spectra are also very much
alike.

Figure 10: Time evolution of volume-averaged kinetic, magnetic energy and
average suppression of thermal diffusivity in 3D with and without external
forcing. The zero point on the $x$-axis denotes the time when external forcing
was activated in the runs. External forcing drives its own SSD and in so doing
revives previously dead MTI simulations with strong whistler suppression.

##  6 3D MTI turbulence with external forcing and whistler suppression

Figure 11: Volume-averaged and spectral energy densities for the 3D runs with
and without external forcing at saturation. From left to right: bar plot of
the volume-averaged turbulent energies, compensated (i.e., multiplied by $k$),
and uncompensated spectral energy densities. The color legend is in common
between the bar plots and the power spectra. The pink hashed bars in the left
column mark the sum of the energy levels of the MTI-only and forcing-only
runs. The peak in the kinetic energy near $l_{\chi}$ indicates the wave
numbers excited by the external forcing. The turbulent driving of the MTI and
the external forcing do not add up exactly for the kinetic energy but nearly
so for the magnetic energy. Moreover, the MTI excites velocity fluctuations
over a wide range of scales below $l_{\chi}$.

We now turn our attention on the impact of external sources of turbulence on
MTI turbulence with whistler suppression. As we show, external forcing can
actually have a beneficial effect, reviving simulations that had entered the
death spiral because the suppression parameter $\alpha$ is above the critical
threshold for the SSD (Sect. 6.1). We then look more in detail at the
interplay between MTI and externally-driven turbulence, trying to distinguish
signatures of the MTI (with and without whistler suppression) in simulations
with external forcing, which may have implications for observations. We find
that the presence of the MTI is clearly recognizable in spectral space, where
we observe an excess of energy in the turbulent fluctuations at small scales,
compared to a forcing-only simulation with similar volume-averaged kinetic and
magnetic energies (Sect. 6.2). Interestingly, we also find that the inclusion
of externally-driven turbulence does not significantly affect the statistical
properties of whistler suppression, except for an overall increase of the
magnetic energy due to the additional SSD effect (Sect. 6.3).

###  6.1 Revival of the “dead” MTI with external forcing

In Fig. 10 we show the time evolution of the volume averages of kinetic
energy, magnetic energy and thermal conductivity for four representative 3D
runs (also discussed in PBP23): the MTI reference run (without whistler
suppression), a forcing-only run with isotropic thermal conductivity without
the MTI, a run with both MTI and external forcing (with no whistler
suppression), and finally a run with MTI, external forcing and whistler
suppression ($\alpha=0.12$) which was firmly in the “death-spiral” regime (see
Fig. 1). We choose the level of external forcing so as to drive turbulence at
saturation of similar amplitude and at similar scale as the MTI reference run
without whistler suppression (grey starred line). Our choice is supported by
previous estimates that the strength of MTI-driven turbulence might be roughly
comparable to the observed levels in the periphery of galaxy clusters (PL22b).
For more details on the forcing mechanism refer to Sect. 3.3.

By construction, the forcing-only run (solid red line) has saturated kinetic
energy levels comparable to those of the MTI (no suppression) reference run
before forcing ($0.19$ vs. $0.16$, in units of $\varv_{\mathrm{cond}}^{2}$),
and similar magnetic energy ($0.045$ vs. $0.040$). Due to the absence of the
MTI, on the other hand, in the forcing-only run the potential energy is
substantially lower ($1.3\times 10^{-3}$ vs. $0.51$, not shown). Combining
(unsuppressed) MTI and external forcing leads unsurprisingly to larger levels
of kinetic and magnetic energies than either the MTI-only or forcing-only
runs, although the external and the MTI turbulent kinetic energies do not
exactly sum up (contrary to what was observed by Parrish et al., 2012).
Interestingly, the magnetic energy at saturation appears to be closer to the
sum of the MTI-only and forcing-only runs.

The most remarkable consequence of adding external forcing, however, can be
seen in the run with strong whistler suppression ($\alpha=0.12$, solid orange
line), where the MTI turbulence was previously decaying: after switching on
external forcing, both turbulent kinetic and magnetic energies quickly
increase and overshoot the levels of the forcing-only run (first and second
rows of Fig. 10), particularly so for the magnetic energy. Similarly, the
volume-averaged thermal diffusivity returns to very high levels, around $67\%$
of the unsuppressed Spitzer value (third row). We interpret this evidence as a
clear sign that the MTI has been successfully revived by the addition of
external forcing, and now contributes to the observed turbulent levels (this
conclusion is also supported by the positive sign of the buoyancy power,
indicating instability, see PBP23, ). The reason for this phenomenon is that
the externally-driven turbulence can also drive its own SSD and, contrary to
what we examined in Sect. 4.3, this process is not affected by how much
suppression of thermal diffusivity there is. Since turbulent magnetic fields
are now propped up by the external forcing, whistler suppression of thermal
diffusivity cannot go below a certain value (determined by the equipartition
value of $\tilde{\beta}$ with the external turbulence). If the residual
thermal diffusivity is high enough, the MTI has a chance to make a comeback
and do its work on top of the external turbulence, which is indeed what we
observe in Fig. 10. Fixing the external forcing, the final turbulent levels of
the MTI with both whistler suppression and external forcing will depend on the
exact value of the suppression parameter, $\alpha$. However, we can
realistically expect that they will span the range between the turbulent
levels of the unsuppressed MTI run with forcing (in the limit of low
$\alpha$), and those of the forcing-only run without the MTI (if $\alpha$ is
so high that the residual thermal diffusivity does not allow the MTI to grow
anymore).

This behaviour is particularly interesting because it constitutes a first
example of a positive interplay between external forcing and the MTI with
whistler suppression, which is clearly non-trivial and need not be entirely
inimical. We can speculate that two limiting regimes may exist: if the
external turbulence is too weak, then the thermal diffusivity cannot be
boosted to a sufficiently large value to allow the MTI to become again
unstable, thus no revival could take place; if the external turbulence is too
strong, however, the MTI can be boosted only up to a certain extent
($\chi/\chi_{0}$ cannot grow larger than 1), and the external turbulence may
then again dominate the turbulence so that the MTI contribution is negligible
(although this should be thought of as a scale-dependent statement, McCourt et
al., 2011). Therefore it appears that a “goldilock” scenario may arise when
the external forcing and the MTI can both operate at comparable strength,
which is precisely the regime we studied and that we consider realistic for
the periphery of galaxy clusters. In order to make conclusive statements,
however, a more detailed exploration of the different forcing amplitudes and
driving scales is required.

###  6.2 Signature of the MTI in forced-turbulence simulations

We now look more in detail at the energy distribution in spectral space at
saturation of our runs with forcing. In fact, despite the similarities between
the volume-averaged energies of the forcing-only and MTI-only runs, their
power spectra differ noticeably. In Fig. 11 we plot the kinetic, magnetic and
potential power spectra (from top to bottom) both in log-log scale (right
column) and in linear scale after multiplying by the wave number $k$, to
highlight the different features (middle column). As we can see from Fig. 11
(top row), the kinetic power spectra of the unsuppressed MTI-only run presents
stronger turbulence at all wavelengths below and above the forcing scale
compared to the forcing-only run. This result underlines the fact that the
energy injection provided by the MTI works very differently from external
forcing: the MTI, in fact, drives the turbulence over a wider range of scales
below the conduction length $l_{\chi}$, whereas in our forcing model all the
energy is dumped in a small set of wavelengths (which results in a very
visible spike near $l_{\chi}$ in the kinetic power spectrum). Because of the
wide-range injection by the MTI, identifying a proper inertial range – defined
by the absence of both injection and dissipation – becomes problematic, and
one should expect deviations from the classic Kolmogorov scaling.
Interestingly, we notice that the magnetic power spectra are much less
affected by the change in the driving, with the MTI-only and forcing-only
lines closely tracing each other over a wide range of scales (middle panel).

Combining the (unsuppressed) MTI and external forcing we observe stronger
velocity fluctuations at all scales (this is consistent with earlier results
by McCourt et al., 2011) whereas the magnetic power spectrum appears to be
boosted by a multiplicative factor, while its shape remains relatively
unchanged. Finally, an intermediate behaviour is observed in the simulation
with whistler suppression and forcing, where the kinetic and magnetic energy
levels are below the run with (unsuppressed) MTI and external forcing, but
above those in either the MTI-only or forcing-only runs. Interestingly, this
trend is reversed in the potential energy power spectra (bottom panel), where
the largest thermal fluctuations are achieved in the MTI-only run without
whistler suppression, while both external forcing and whistler suppression
lead to a decrease in their amplitude, particularly at large scales.

Figure 12: Probability densities and marginal distributions with and without
external forcing and with and without whistler-suppression. The overall
magnetic field strength mainly depends on whether external forcing is
included, while the degree of isothermality is almost exclusively affected by
whistler-suppression. From bottom-left clockwise: two-dimensional probability
density as a function of $\tilde{\beta}$ and $\sigma$ for the MTI run with
whistler suppression $\alpha=0.12$ and external forcing, marginal distribution
over $\tilde{\beta}$, marginal distribution over $\chi/\chi_{0}$, and marginal
distribution over the isothermality parameter $\sigma$. In the two-dimensional
plot the solid purple line denotes the isocontour of $10^{-4}$ for the MTI run
without forcing and $\alpha=0.024$, which shows approximately the same values
of the volume-averaged $\chi/\chi_{0}$.

###  6.3 Statistics of whistler suppression unaffected by external forcing

In addition to volume-averaged and spectral diagnostics, we are also
interested in how the statistics of MTI turbulence with whistler suppression
are affected by the external forcing. In Fig. 12 we plot the two-dimensional
probability density as a function of $\tilde{\beta}$ and $\sigma$ for the MTI
run with strong whistler suppression ($\alpha=0.12$) and external forcing,
together with the marginal distributions over $\tilde{\beta}$, $\sigma$, and
$\chi/\chi_{0}$. We notice that the addition of external forcing does not
significantly affect the probability distribution, which is qualitatively very
similar to the runs without forcing analyzed in Sect. 5.2: the values are
distributed within this roughly triangular shape, with a prominent lobe at
high-$\tilde{\beta}$ representing portions of the domain out of isothermality,
whereas the majority of the volume is clustered in the region with
$\sigma\in[10^{-1},10^{0}]$ and $\tilde{\beta}\sim 10^{1}-10^{2}$.

In the bottom-left plot of Fig. 12 we also show the isocontour of the two-
dimensional probability (level $10^{-4}$, solid purple line) of the run with
MTI without forcing and $\alpha=0.024$: this run has similar volume-averaged
suppression $\langle\chi/\chi_{0}\rangle$ as the MTI run with strong whistler
suppression ($\alpha=0.12$) and external forcing ($68\%$ vs. $64\%$,
respectively). Thus, it offers us the possibility to compare two very
different scenarios: one with weak whistler suppression and no external
forcing, while the other with strong whistler suppression and forcing, but in
the end with similar volume-averaged thermal diffusivities. Interestingly, the
main difference between the two runs is an overall shift of the distribution
towards lower $\tilde{\beta}$ (higher magnetic fields) when external forcing
is included.

For a more quantitative comparison, we can look at the marginal distributions
of Fig. 12, where we compare two runs with MTI and external forcing
($\alpha=0$ and $\alpha=0.12$) with two MTI-only runs ($\alpha=0$ and
$\alpha=0.024$). From the marginal distribution over $\tilde{\beta}$ (top-
left), we notice again how stronger levels of turbulence (due to external
turbulence or lower $\alpha$) result in an overall shift of the probability
distribution towards lower values of $\tilde{\beta}$, while the shape remains
relatively unchanged. This behaviour is similar to what we observed in Fig. 8,
and is due to a stronger driving of the SSD. We distinguish a split in the
marginal distribution over the isothermality parameter $\sigma$ (bottom-right)
between the runs with and without whistler suppression. In the runs without
whistler suppression, the distributions fall off rapidly for values above
$\sigma=1$, and in the runs with whistler suppression, the distributions
develop fat tails for $\sigma>1$. This behaviour suggests that the level of
isothermality enforced by the MTI is not influenced by external forcing, and
that it only depends on the effective thermal diffusivity. Further support to
this picture is given in the top-right panel, where we plot the marginal
distribution of the two runs with whistler suppression with respect to
$\chi/\chi_{0}$: there we see confirmation that not only the volume-averaged
$\langle\chi/\chi_{0}\rangle$ between the two runs is similar, but also the
frequency with which the individual values occur.

In conclusion, there seems to be a level of degeneracy between the suppression
parameter $\alpha$ and the magnetic field strength (set by total turbulence
levels in MTI simulations with external forcing) where external forcing can
effectively compensate for larger values of $\alpha$, and produce similar
levels of thermal diffusivity. This similarity is manifested in certain
aspects (as the distribution over the isothermality parameter $\sigma$), but
not others: with external forcing, the total magnetic energy is higher,
whereas the potential energy and the buoyancy power (the rate at which energy
is added by the buoyancy force) are both lower ($0.25$ vs. $0.43$ and $0.002$
vs. $0.007$, in units of $\varv_{\mathrm{cond}}^{2}$ and $\epsilon_{0}$,
respectively).

##  7 Discussion: can we model the effect of whistler suppression on the MTI
as a “flat” decrease in thermal conductivity?

We pick up the threads of these last sections and summarize our findings on
whether it is appropriate or not to model whistler suppression numerically as
a “flat” decrease in thermal conductivity, and what physics one might miss by
doing so. We stress that this summary is necessarily focused on the
differences seen in Boussinesq MTI-driven turbulence and that we cannot
exclude significant differences with more sophisticated physical models that
take into account, e.g., radiative physics, galaxy formation, etc.

In Sect. 4.2 we observed that for modest suppression below the Spitzer value
of thermal diffusivity, whistler suppression behaves very similarly to a flat
uniform suppression in terms of the turbulent levels and other volume-averaged
quantities. This is consistent with the theoretical model proposed in PL22a
where the turbulent root-mean-square velocities are determined by the MTI
energy injection rate and the largest scale excited by the turbulence, both of
which depend on $\chi$. A similar picture emerges from looking at the power
spectra in Sect. 5.3, which look remarkably similar to those computed for
simulations with a flat uniform suppression. Despite some noticeable
differences in the MTI energy injection and anisotropic dissipation with
whistler suppression, these differences are averaged out such that the net
energy injection rate in wave number space is very close to the uniform
suppression model.

Figure 13: Representation of the “death-spiral” mechanism with whistler
suppression. When the suppression parameter $\alpha$ is increased beyond a
critical value, the turbulence weakens and the magnetic field strength
decreases to maintain equipartition. As a result whistler suppression gets
stronger and the turbulence weakens further.

One new notable feature of whistler suppression is that regions with high
plasma $\beta$ are “locally” brought out of isothermality (Sect. 5.2), which
is an effect that would not be observed with a flat uniform suppression. The
most striking piece of new physics introduced by whistler suppression,
however, is that when the level of suppression becomes too large the SSD is
shut down and the system enters into what we called the “death-spiral”. An
illustration of the death-spiral can be found in Fig. 13 and can be summarized
as follows: starting from a turbulent equilibrium state, when the suppression
parameter $\alpha$ increases beyond a critical value (either by tuning it
manually like in our local MTI simulations, or because the physical conditions
in the environment change) the volume-averaged thermal diffusivity decreases
and so does the strength of the MTI-driven turbulence in order to keep the
linear relationship between $\langle\chi\rangle$ and the turbulent kinetic
energy. Once the strength of the turbulence becomes weaker, the SSD follows
suit and the total magnetic energy adjusts to the new (lower) equipartition-
value. However (and here lies the difference between whistler suppression and
a flat uniform suppression), the lower magnetic energy directly influences the
amount of suppression of thermal diffusivity through the $1/\beta$ dependence,
which in turn leads to weaker turbulence and so forth. Even though the
suppression parameter $\alpha$ is likely on average below the critical value
for the shutdown of the SSD in the cluster periphery (Sect. 4.3), different
regions of the cluster periphery may have very different conditions (e.g.,
they could experience a merger that temporarily stirs the ICM so that the
pressure scale height locally increases, or have locally weaker magnetic
fields or stronger temperature gradients), and therefore a shutdown of the MTI
cannot be excluded a priori everywhere. By using a flat suppression model this
behaviour would be missed entirely, which underscores the importance of
adopting well-motivated subgrid models.

##  8 Summary and conclusions

In this paper, which complements and expands our previous work (PBP23), we
carried out an extensive investigation of the effects of whistler suppression
of thermal conductivity on the MTI through small-scale simulations. We adopted
a popular subgrid model for the thermal conductivity which smoothly
interpolates between the Spitzer (at low plasma $\beta$) and the whistler-
suppressed regime (at high $\beta$) and which is inspired by theoretical
arguments (Drake et al., 2021) and kinetic simulations (Roberg-Clark et al.,
2016, 2018; Komarov et al., 2018).

Our model is characterized by a dimensionless parameter ($\alpha$) that
encodes the information on the collisionality of a local patch of ICM plasma
and the efficiency of thermal conduction on cluster-size scales. By varying
$\alpha$ in the range that we can realistically expect in the periphery of
galaxy clusters, we were able to study the consequences of weaker (low
$\alpha$) or stronger suppression (large $\alpha$) of thermal conductivity on
the MTI, which relies on effective conduction of heat along magnetic field
lines, and that drives turbulence proportional to the thermal diffusivity.

We found that for relatively weak suppression below the Spitzer value, MTI
turbulence with whistler suppression appears phenomenologically rather similar
to a simpler model of “flat” (i.e. spatially uniform) suppression, in that the
turbulent energy levels decrease proportionally to the volume-averaged
suppressed thermal diffusivity (Perrone & Latter, 2022b). In a similar fashion
the turbulent power spectra also resemble their counterparts with flat
suppression. This is in spite of the fact that whistler suppression is highly
inhomogeneous in space, and that locally regions with high plasma $\beta$
(weak magnetic fields) can be brought away from isothermality, which is
instead effectively enforced by the MTI elsewhere. However, the similarities
end when the suppression parameter $\alpha$ nears a critical value (which
depends on the magnetic Reynolds number of the plasma), beyond which the MTI
enters a new regime, the “death-spiral”, and the turbulence abruptly dies.
This critical transition comes about via the suppression of the MTI-driven
SSD, and represents one of the main results of this work. While we estimate
that in the hot ICM the suppression parameter $\alpha$ is below the critical
value, the physical conditions in the periphery of galaxy clusters can vary
quite significantly in space and time, and thus simpler models of flat
suppression should be used with care.

This unique connection between MTI turbulence and whistler suppression can be
severed in the case where magnetic fields are propped up by some other
physical process than the MTI. This leads us to the second main result of this
study, which is that externally-driven turbulence (representing, e.g.,
turbulence driven by cluster mergers) can in be fact beneficial for the MTI
with whistler suppression. The external turbulence sets a minimum value of
thermal diffusivity due to its own SSD, on top of which the MTI can then add
its contribution. In this way, simulations where the MTI was previously dying
can be be revived, and the resulting turbulent levels are greater than in
analogous forcing-only runs. This shows that under certain circumstances the
interaction between two sources of turbulence can be constructive and not
necessarily detrimental.

###  8.1 Limitations of this study

This work comes with several limitations which we can broadly categorize as
(i) intrinsic to the specific subgrid model chosen for whistler suppression,
and (ii) pertaining the application of our results to galaxy clusters. In the
first camp, the main unknown is certainly to what extent the closure for
thermal conductivity that we adopted represents an accurate description of the
effect of whistler waves on electron transport in the bulk of the ICM. There
are a number of reasons why the model in Eq. (6) can only constitute a partial
approximation of the “true” thermal conductivity in galaxy clusters. First of
all, the whistler-dominated regime for the heat flux (the $1/\beta$
suppression) has been convincingly shown only in kinetic simulations with very
strong temperature gradients (on scales not much longer than electron
gyroradii, see, e.g., Roberg-Clark et al., 2016, 2018; Komarov et al., 2018),
which are necessary to drive the whistler instability. Such small-scale
temperature gradients are unlikely to be long-lived in galaxy clusters, with
typical mean-free paths $14$ orders of magnitude larger than the electron
gyroradii. Other ways of exciting whistler waves in the ICM exist (e.g., they
could appear as side-products of mirror modes, driven by the anisotropic
pressure rather than a temperature gradient, Ley et al., 2023), but if so the
appropriateness of the closure introduced in Eq. (5) for the heat flux becomes
less clear. Other sources of uncertainty lie in the fact that the
interpolation formula in Eq. (6) implicitly assumes a hierarchy of temporal
and spatial scales between the slow, large-scale fluid motions (such as those
induced by the MTI) and the fast, micro-scale whistler instability, with the
latter instantaneously adjusting to the changing conditions of the plasma,
without interference by the large-scale turbulence. Even though this is
plausible on the fluid scales that we simulated in this study, in theory the
MTI (or its kinetic counterpart) can be excited at the same electron scales
where the whistler instability develops (Xu & Kunz, 2016). As a result, the
interplay between the two instabilities is likely to be more complex than
anticipated. Finally, there is currently no consensus on what the combined
impact of whistler, mirror and firehose instabilities on the thermal
conductivity could look like, and therefore these additional effects were also
neglected in our study.

While we believe that our choice to perform small-scale simulations using the
Boussinesq approximation has proved very fruitful, as it allowed us to gain a
solid understanding of the nonlinear interplay between the MTI and whistler
suppression of thermal diffusivity, which would have been hard to obtain
otherwise, it inevitably limits the applicability of our findings to model the
large-scale dynamics in the periphery of galaxy clusters. Specifically,
cosmological-style simulations have so far failed to show robust evidence of
the MTI, even without any sort of suppression of thermal conductivity
(Ruszkowski et al., 2011; Kannan et al., 2017; Barnes et al., 2019). While
this may be due to numerical limitations (for a detailed discussion see PL22b,
), the straight-up adoption of a subgrid closure as in Eq. (6) in cosmological
simulations clearly would not help address this problem. Granted, the
suppression of thermal conduction is relevant to other physical processes in
the ICM than just the MTI, such as the formation and survival of cold fronts
(Markevitch & Vikhlinin, 2007), and the energy balance between cooling and
heat transport (Zakamska & Narayan, 2003). However, the consequences of the
two main results of this work, the shutdown of the SSD and the beneficial
effect of external turbulence on the MTI (both very relevant to the ongoing
debate on the source of relativistic electrons in the periphery of galaxy
clusters; Cuciti et al., 2022; Beduzzi et al., 2023), would be currently very
difficult to assess within large-scale cosmological simulations.

###  8.2 Future outlook

Despite earlier assessments of its inhibition due to competing external
turbulence and suppression of thermal conductivity by microinstabilities, the
MTI has proven to be remarkably resilient, and in this work and in its
companion paper (PBP23) we have shown that one problem can alleviate the
other. This clearly shows that the MTI is a complex dynamical process that
requires detailed modelling beyond linear theory in order to draw firm
conclusions. It should also serve as a cautionary tale to carefully consider
the consequences of adopting certain subgrid models, and what new physics one
might introduce in the problem by doing so. Future work is clearly needed on
two fronts: on the small scales, with the aim to strengthen our understanding
of the impact of the whistler and other microinstabilities on the transport of
electrons, and on the large scales, focusing on thorough studies of the MTI in
a cosmological context.

## Acknowledgements

The authors acknowledge support by the European Research Council under ERC-AdG
grant PICOGAL-101019746. TB gratefully acknowledges funding from the European
Union’s Horizon Europe research and innovation programme under the Marie
Skłodowska-Curie grant agreement No 101106080. LMP thanks Henrik Latter for
helpful and generous advice on a draft version of the manuscript. Parts of the
results in this work make use of the colormaps in the CMasher package (van der
Velden, 2020).

## Data Availability

The data underlying this article will be shared on reasonable request to the
corresponding author.

## References

  * Asai et al. (2007) Asai, N., Fukuda, N., & Matsumoto, R. 2007, The Astrophysical Journal, 663, 816
  * Balbus (2000) Balbus, S. A. 2000, The Astrophysical Journal, 534, 420
  * Bale et al. (2009) Bale, S. D., Kasper, J. C., Howes, G. G., et al. 2009, Physical Review Letters, 103, 211101
  * Bale et al. (2013) Bale, S. D., Pulupa, M., Salem, C., Chen, C. H. K., & Quataert, E. 2013, The Astrophysical Journal, 769, L22
  * Barnes et al. (2019) Barnes, D. J., Kannan, R., Vogelsberger, M., et al. 2019, Monthly Notices of the Royal Astronomical Society, 488, 3003
  * Beckmann et al. (2022) Beckmann, R. S., Dubois, Y., Pellissier, A., Polles, F. L., & Olivares, V. 2022, Astronomy & Astrophysics, 666, A71
  * Beduzzi et al. (2023) Beduzzi, L., Vazza, F., Brunetti, G., et al. 2023, On the Origin of Mega Radiohalos
  * Berlok et al. (2021) Berlok, T., Quataert, E., Pessah, M. E., & Pfrommer, C. 2021, Monthly Notices of the Royal Astronomical Society, 504, 3435
  * Birdsall & Langdon (1991) Birdsall, C. K. & Langdon, A. B. 1991, Plasma Physics via Computer Simulation
  * Braginskii (1965) Braginskii, S. I. 1965, Reviews of Plasma Physics, 1, 205
  * Carilli & Taylor (2002) Carilli, C. L. & Taylor, G. B. 2002, Annual Review of Astronomy and Astrophysics, 40, 319
  * Chandrasekhar et al. (1958) Chandrasekhar, S., Kaufman, A. N., & Watson, K. M. 1958, Proceedings of the Royal Society of London. Series A, Mathematical and Physical Sciences, 245, 435
  * Chen et al. (2016) Chen, C. H. K., Matteini, L., Schekochihin, A. A., et al. 2016, The Astrophysical Journal, 825, L26
  * Churazov et al. (2003) Churazov, E., Forman, W., Jones, C., & Böhringer, H. 2003, 590, 225
  * Cowie & McKee (1977) Cowie, L. L. & McKee, C. F. 1977, The Astrophysical Journal, 211, 135
  * Cuciti et al. (2022) Cuciti, V., de Gasperin, F., Brüggen, M., et al. 2022, Nature, 609, 911
  * Drake et al. (2021) Drake, J. F., Pfrommer, C., Reynolds, C. S., et al. 2021, The Astrophysical Journal, 923, 245
  * Dursi & Pfrommer (2008) Dursi, L. J. & Pfrommer, C. 2008, The Astrophysical Journal, 677, 993
  * Gary (1993) Gary, S. P. 1993, Theory of Space Plasma Microinstabilities
  * Gary & Li (2000) Gary, S. P. & Li, H. 2000, The Astrophysical Journal, 529, 1131
  * Gary & Wang (1996) Gary, S. P. & Wang, J. 1996, Journal of Geophysical Research: Space Physics, 101, 10749
  * Ghirardini et al. (2019) Ghirardini, V., Eckert, D., Ettori, S., et al. 2019, Astronomy & Astrophysics, 621, A41
  * Jacob & Pfrommer (2017a) Jacob, S. & Pfrommer, C. 2017a, Monthly Notices of the Royal Astronomical Society, 467, 1449
  * Jacob & Pfrommer (2017b) Jacob, S. & Pfrommer, C. 2017b, Monthly Notices of the Royal Astronomical Society, 467, 1478
  * Kannan et al. (2017) Kannan, R., Vogelsberger, M., Pfrommer, C., et al. 2017, The Astrophysical Journal, 837, L18
  * Komarov et al. (2018) Komarov, S., Schekochihin, A., Churazov, E., & Spitkovsky, A. 2018, Journal of Plasma Physics, 84, 905840305
  * Kuchar & Enßlin (2011) Kuchar, P. & Enßlin, T. A. 2011, Astronomy and Astrophysics, 529, A13
  * Kunz (2011) Kunz, M. W. 2011, Monthly Notices of the Royal Astronomical Society, 417, 602
  * Kunz et al. (2012) Kunz, M. W., Bogdanovic, T., Reynolds, C. S., & Stone, J. M. 2012, The Astrophysical Journal, 754, 122
  * Kunz et al. (2014) Kunz, M. W., Schekochihin, A. A., & Stone, J. M. 2014, Physical Review Letters, 112, 205003
  * Lesur (2015) Lesur, G. 2015, Snoopy: General Purpose Spectral Solver
  * Levinson & Eichler (1992) Levinson, A. & Eichler, D. 1992, The Astrophysical Journal, 387, 212
  * Ley et al. (2023) Ley, F., Zweibel, E. G., Miller, D., & Riquelme, M. 2023, Secondary Whistler and Ion-cyclotron Instabilities Driven by Mirror Modes in Galaxy Clusters
  * Ley et al. (2022) Ley, F., Zweibel, E. G., Riquelme, M., et al. 2022, A Heating Mechanism via Magnetic Pumping in the Intracluster Medium
  * Lyutikov (2006) Lyutikov, M. 2006, Monthly Notices of the Royal Astronomical Society, 373, 73
  * Markevitch & Vikhlinin (2007) Markevitch, M. & Vikhlinin, A. 2007, 443, 1
  * McCourt et al. (2011) McCourt, M., Parrish, I. J., Sharma, P., & Quataert, E. 2011, Monthly Notices of the Royal Astronomical Society, 413, 1295
  * Melville et al. (2016) Melville, S., Schekochihin, A. A., & Kunz, M. W. 2016, Monthly Notices of the Royal Astronomical Society, 459, 2701
  * Mogavero & Schekochihin (2014) Mogavero, F. & Schekochihin, A. A. 2014, Monthly Notices of the Royal Astronomical Society, 440, 3226
  * Narayan & Medvedev (2001) Narayan, R. & Medvedev, M. V. 2001, The Astrophysical Journal, 562, L129
  * Parker (1958) Parker, E. N. 1958, Physical Review, 109, 1874
  * Parrish et al. (2012) Parrish, I. J., McCourt, M., Quataert, E., & Sharma, P. 2012, Monthly Notices of the Royal Astronomical Society: Letters, 419, L29
  * Parrish & Stone (2007) Parrish, I. J. & Stone, J. M. 2007, The Astrophysical Journal, 664, 135
  * Perrone et al. (2023) Perrone, L. M., Berlok, T., & Pfrommer, C. 2023, Does the Magneto-Thermal Instability Survive Whistler-Suppression of Thermal Conductivity in Galaxy Clusters?
  * Perrone & Latter (2022a) Perrone, L. M. & Latter, H. 2022a, Monthly Notices of the Royal Astronomical Society, 513, 4605
  * Perrone & Latter (2022b) Perrone, L. M. & Latter, H. 2022b, Monthly Notices of the Royal Astronomical Society, 513, 4625
  * Pfrommer & Dursi (2010) Pfrommer, C. & Dursi, J. 2010, Nature Physics, 6, 520
  * Pistinner & Eichler (1998) Pistinner, S. L. & Eichler, D. 1998, Monthly Notices of the Royal Astronomical Society, 301, 49
  * Quataert (2008) Quataert, E. 2008, The Astrophysical Journal, 673, 758
  * Rappaz & Schober (2023) Rappaz, Y. & Schober, J. 2023, The Effect of Pressure-Anisotropy-Driven Kinetic Instabilities on Magnetic Field Amplification in Galaxy Clusters
  * Riquelme et al. (2016) Riquelme, M., Quataert, E., & Verscharen, D. 2016, The Astrophysical Journal, 824, 123
  * Roberg-Clark et al. (2016) Roberg-Clark, G. T., Drake, J. F., Reynolds, C. S., & Swisdak, M. 2016, The Astrophysical Journal, 830, L9
  * Roberg-Clark et al. (2018) Roberg-Clark, G. T., Drake, J. F., Reynolds, C. S., & Swisdak, M. 2018, Physical Review Letters, 120, 035101
  * Ruszkowski et al. (2011) Ruszkowski, M., Lee, D., Brüggen, M., Parrish, I., & Oh, S. P. 2011, The Astrophysical Journal, 740, 81
  * Ruszkowski & Oh (2010) Ruszkowski, M. & Oh, S. P. 2010, The Astrophysical Journal, 713, 1332
  * Sarazin (1988) Sarazin, C. L. 1988, X-Ray Emission from Clusters of Galaxies
  * Schekochihin & Cowley (2006) Schekochihin, A. A. & Cowley, S. C. 2006, Physics of Plasmas, 13, 056501
  * Schekochihin et al. (2008) Schekochihin, A. A., Cowley, S. C., Kulsrud, R. M., Rosin, M. S., & Heinemann, T. 2008, Physical Review Letters, 100, 081301
  * Schuecker et al. (2004) Schuecker, P., Finoguenov, A., Miniati, F., Böhringer, H., & Briel, U. G. 2004, Astronomy & Astrophysics, 426, 387
  * Schwarzschild (1906) Schwarzschild, K. 1906, Nachrichten von der Gesellschaft der Wissenschaften zu Göttingen, Mathematisch-Physikalische Klasse, 1906, 41
  * Sharma et al. (2009) Sharma, P., Chandran, B. D. G., Quataert, E., & Parrish, I. J. 2009, AIP Conference Proceedings, 1201, 363
  * Sharma & Hammett (2007) Sharma, P. & Hammett, G. W. 2007, Journal of Computational Physics, 227, 123
  * Spitzer (1962) Spitzer, L. 1962, Physics of Fully Ionized Gases
  * van der Velden (2020) van der Velden, E. 2020, The Journal of Open Source Software, 5, 2004
  * Vogt & Ensslin (2005) Vogt, C. & Ensslin, T. A. 2005, Astronomy & Astrophysics, 434, 67
  * Xu & Kunz (2016) Xu, R. & Kunz, M. W. 2016, Journal of Plasma Physics, 82, 905820507
  * Yang & Reynolds (2016) Yang, H.-Y. K. & Reynolds, C. S. 2016, The Astrophysical Journal, 818, 181
  * Zakamska & Narayan (2003) Zakamska, N. L. & Narayan, R. 2003, The Astrophysical Journal, 582, 162
  * Zeldovich (1957) Zeldovich, Y. B. 1957, Soviet Physics–JETP [translation of Zhurnal Eksperimentalnoi i Teoreticheskoi Fiziki], 4, 460
  * Zhou et al. (2023) Zhou, M., Zhdankin, V., Kunz, M. W., Loureiro, N. F., & Uzdensky, D. A. 2023, Magnetogenesis in a Collisionless Plasma: From Weibel Instability to Turbulent Dynamo
  * Zhuravleva et al. (2019) Zhuravleva, I., Churazov, E., Schekochihin, A. A., et al. 2019, Nature Astronomy, 3, 832
  * ZuHone et al. (2013) ZuHone, J. A., Markevitch, M., Ruszkowski, M., & Lee, D. 2013, The Astrophysical Journal, 762, 69
  * ZuHone & Roediger (2016) ZuHone, J. A. & Roediger, E. 2016, Journal of Plasma Physics, 82

##  Appendix A Whistler-suppressed diffusivity in SNOOPY

In this Appendix we show how to rewrite the analytical closure for whistler-
suppressed diffusivity in terms of the Boussinesq variables. Consistent with
the MHD approximation, we consider a two-species, quasi-neutral plasma
($n_{\mathrm{i}}\simeq n_{\mathrm{e}}$) with equal ion and electron
temperatures ($T_{\mathrm{i}}\simeq T_{\mathrm{e}}\simeq T$).

The temperature gradient parallel to the magnetic field $L_{T,\parallel}$ in
Eq. (6) can be split into two terms: the contributions coming from the
background profile ($\mathrm{d}T/\mathrm{d}z|_{0}$) plus the contribution of
the local temperature perturbation ($\bm{b}\bm{\cdot}\bm{\nabla}\delta T$) as
follows:

$\displaystyle L_{T,\parallel}^{-1}$
$\displaystyle\equiv\left\lvert\frac{\bm{b}\bm{\cdot}\bm{\nabla}T}{T}\right\rvert\simeq\left\lvert
b_{z}\left(\frac{\mathrm{d}\ln
T}{\mathrm{d}z}\right)_{0}+(\bm{b}\bm{\cdot}\bm{\nabla})\frac{\delta
T}{T_{0}}\right\rvert$
$\displaystyle=\frac{\omega_{\mathrm{T}}^{2}}{g_{0}}\left\lvert
b_{z}+\omega_{\mathrm{T}}^{-2}\bm{b}\bm{\cdot}\bm{\nabla}\theta\right\rvert,$
(40)

where we have used the definition of the local MTI frequency in Eq. (8)
together with the definition of the buoyancy variable $\theta$ (Eq. 13), and
where $g_{0}$ is the local gravitational acceleration. Note that in writing
Eq. (40) we have neglected terms that are higher-order than
$\mathcal{O}(\delta T/T_{0})$, consistently with the Boussinesq approximation.
Further, we assume hydrostatic equilibrium and write the pressure scale-height
as $H=c_{\mathrm{s}}^{2}/g_{0}$, where $c_{\mathrm{s}}^{2}=\gamma
T/m_{\mathrm{p}}$ is the adiabatic sound speed for an hydrogenic plasma
($m_{\mathrm{p}}$ is the proton mass). With these definitions, Eq. (6) can be
rewritten as follows:

$\displaystyle\chi=\dfrac{\chi_{0}}{1+\dfrac{1}{3}\dfrac{1}{\gamma}\left(\dfrac{\lambda_{\text{mfp,e}}}{H}\right)\left(\dfrac{H^{2}}{\chi_{0}/\omega_{\mathrm{T}}}\right)\left(\dfrac{2\chi_{0}\omega_{\mathrm{T}}}{\varv_{\mathrm{A}}^{2}}\right)\left\lvert
b_{z}+\dfrac{\bm{b}\bm{\cdot}\bm{\nabla}\theta}{\omega_{\mathrm{T}}^{2}}\right\rvert}$
(41)

where we wrote the (electron) plasma beta as
$\beta=2c_{\mathrm{s}}^{2}/(\gamma\varv_{\mathrm{A}}^{2})$. Replacing the
definitions for the modified plasma $\tilde{\beta}$ (Eq. 20), the suppression
parameter $\alpha$ (Eq. 19), and the isothermality parameter (Eq. 21) we then
obtain Equation (18).

##  Appendix B Interruption of MTI growth with whistler suppression

In this Appendix we show how the presence of whistler suppression of thermal
diffusivity can interrupt the exponential growth of the MTI modes. As we
discuss below, the effect of whistler suppression can be more or less severe,
depending on the wave number of the MTI eigenmode, with modes growing on
scales much smaller than the conduction length $l_{\chi}$ being the least
affected. We consider the onset of the MTI with whistler suppression from an
initial state of equilibrium with $\bm{u}_{0}=\mathbf{0}$,
$\bm{B}_{0}=B_{0}\bm{e}_{x}$, $\theta_{0}=0$ (with constant gravity pointing
in the $-\bm{e}_{z}$ direction), and looking at small perturbations
$\delta\bm{u}^{*},\delta\bm{B}^{*},\delta\theta^{*}\sim\mathcal{O}(\epsilon)$,
where $\epsilon$ is a small parameter. The starred quantities have been
nondimensionalized using $l_{\chi,0}=(\chi_{0}/\omega_{\mathrm{T}})^{1/2}$ and
$\omega_{\mathrm{T}}$ as the units of length and time, respectively, as
follows:

$\displaystyle\delta\bm{u}^{*}\equiv\frac{\delta\bm{u}}{l_{\chi,0}\omega_{\mathrm{T}}},\hskip
20.00003pt\delta\bm{B}^{*}\equiv\frac{\delta\bm{B}}{B_{0}},\hskip
20.00003pt\delta\theta^{*}\equiv\frac{\delta\theta}{l_{\chi,0}\omega_{\mathrm{T}}^{2}}.$
(42)

For ease of notation, throughout this section we drop the “0” subscript in
$l_{\chi,0}$, with the understanding that it always refers to the conduction
length computed using the unsuppressed diffusivity $\chi_{0}$. At first order
in $\epsilon$, Eq. (18) then reads:

$\displaystyle\chi=\chi_{0}\left[1-\frac{1}{3}\alpha\left(\frac{2\chi_{0}\omega_{\mathrm{T}}}{B_{0}^{2}}\right)\left\lvert\delta
B_{z}^{*}+\partial_{x}^{*}\theta^{*}\right\rvert+\mathcal{O}(\epsilon^{2})\right],$
(43)

where $\partial_{x}^{*}$ is the partial derivative with respect to the
nondimensionalized spatial variable variable $x^{*}\equiv x/l_{\chi}$. Because
the magnetic field lines are initially isothermal, at zeroth order in
$\epsilon$ the heat diffusivity is just $\chi_{0}$ and, consequently, the
linearized buoyancy equation (Eq. 17) is the same as that of the MTI without
whistler suppression. As the MTI eigenmode grows exponentially, however, the
magnetic field lines are progressively brought farther away from isothermality
(quantified by the term $\left\lvert\delta
B_{z}^{*}+\partial_{x}^{*}\theta^{*}\right\rvert$ in Eq. 43), thus reducing