11institutetext: Univ Lyon, Univ Lyon1, Ens de Lyon
CNRS, Centre de Recherche Astrophysique de Lyon UMR5574
F-69230, Saint-Genis,-Laval, France

#  $\mathcal{PT}$ and anti-$\mathcal{PT}$ symmetries for astrophysical waves

Armand Leclerc<EMAIL_ADDRESS>Guillaume Laibe 11 Nicolas Perez
11

(Received September 15, 1996; accepted March 16, 1997)

###### Abstract

Context. Discrete symmetries have found numerous applications in photonics and
quantum mechanics, but remain little studied in fluid mechanics, particularly
in astrophysics.

Aims. We aim to show how $\mathcal{PT}$ and anti-$\mathcal{PT}$ symmetries
determine the behaviour of linear perturbations in a wide class of
astrophysical problems. They set the location of Exceptional Points in the
parameter space and the associated transitions to instability, and are
associated to the conservation of quadratic quantities that can be determined
explicitly.

Methods. We study several classical local problems: the gravitational
instability of isothermal spheres and thin discs, the Schwarzschild
instability, the Rayleigh-Bénard instability and acoustic waves in dust-gas
mixtures. We calculate the locations and the order of the Exceptional Points
with a method of resultant, as well as the conserved quantities in the
different regions of the parameter space using Krein theory.

Results. All problems studied here exhibit discrete symmetries, even though
Hermiticity is broken by different physical processes (self-gravity, buoyancy,
diffusion, drag). This analysis provides genuine explanations for certain
instabilities, and for the existence of regions in the parameter space where
waves do not propagate. Those correspond to breaking of $\mathcal{PT}$ and
anti-$\mathcal{PT}$ symmetries respectively. Not all instabilities are
associated to symmetry breaking (e.g. the Rayleigh-Benard instability).

###### Key Words.:

Waves – Instabilities – Methods: analytical

###### doi:

10.1103/PhysRevLett.80.5243

###### doi:

10.1063/1.532860

###### doi:

10.1103/PhysRevLett.89.270401

###### doi:

10.1103/PhysRevD.70.025001

###### doi:

10.48550/arXiv.astro-ph/9908095

###### doi:

10.1063/5.0088936

###### doi:

10.1364/OL.32.002632

###### doi:

10.1038/nphys4323

###### doi:

10.1038/nmat3495

###### doi:

10.1126/science.1258479

###### doi:

10.1038/s41566-017-0031-1

###### doi:

10.1088/1367-2630/aba38f

###### doi:

10.1103/PhysRevLett.103.093902

###### doi:

10.1103/PhysRevLett.101.080402

###### doi:

10.1111/j.1365-2966.2011.19291.x

###### doi:

10.1103/PhysRevLett.106.213901

###### doi:

10.1093/mnras/stac3653

###### doi:

10.1103/PhysRevLett.100.103904

###### doi:

10.1088/0305-4470/37/48/009

###### doi:

10.1038/nphys2927

###### doi:

10.1063/1.5088498

###### doi:

10.1103/PhysRevE.104.015215

###### doi:

10.1038/nature11298

###### doi:

10.1038/nphys1515

##  1 Introduction

A large class of astrophysical systems, such as stars or discs, are often
treated as fluids (e.g. Pringle and King 2007; Armitage 2010). The properties
and the stability of the corresponding linear modes are usually analysed with
methods that, to our knowledge, do not fully exploit the symmetries of the
perturbed system with respect to parity $\mathcal{P}$ and time $\mathcal{T}$.

Some of the most striking properties of $\mathcal{PT}$-symmetric systems
emerged in the quantum physics community since the seminal paper by Bender and
Boettcher (1998), the first non-Hermitian quantum Hamiltonian situation with a
real spectrum, which paved the way to open quantum systems. This framework has
been extended to situations exhibiting a spontaneous $\mathcal{PT}$ symmetry
breaking phase transition (Bender et al., 1999, 2002). Generalisations to
quantum field theory have been developed in Bender et al. (2004) and
considerations on the observables of such systems is found in Mostafazadeh and
Batal (2004). $\mathcal{PT}$ symmetry have been a particularly powerful tool
in optics and photonics, where it has led to the development of novel designs
with special properties. The analogy with quantum mechanics comes from the
fact that beam dynamics is governed by a Schrödinger equation, where the
optical index behaves as a complex potential, which can be designed to satisfy
$\mathcal{PT}$ symmetry (El-Ganainy et al., 2007; Makris et al., 2008).
Klaiman et al. (2008) identified a spontaneous breaking of $\mathcal{PT}$
symmetry in a waveguide. Musslimani et al. (2008) investigated the nonlinear
propagation of optical solitons in $\mathcal{PT}$-symmetric refraction
indexes. Guo et al. (2009); Rüter et al. (2010) made the first observations of
$\mathcal{PT}$ symmetry in optics with a complex index. A decade of
developments of $\mathcal{PT}$ symmetry in photonics lead to a multitude of
applications (e.g. Feng et al. 2017; El-Ganainy et al. 2018; Özdemir et al.
2019). Feng et al. (2014) use spontaneous breaking of $\mathcal{PT}$ symmetry
to design a pure single-mode laser. Peng et al. (2014) propose a design for a
low-power optical diode, where the nonlinear regime of the $\mathcal{PT}$
symmetry induces non-reciprocal propagation. Lin et al. (2011); Regensburger
et al. (2012); Feng et al. (2013) use $\mathcal{PT}$ symmetry in periodic
systems (through gratings, periodic crystals, or lattices) to design
unidirectional invisibility. $\mathcal{PT}$-symmetric degeneracies are now
used for their high resonant sensitivity in optical cavities (Hodaei et al.,
2017). Similar interest has been aroused by the use of anti-$\mathcal{PT}$
symmetry for novel optical designs (e.g. Ge and Türeci 2013; Zhang et al.
2020a, b). The topological properties of $\mathcal{PT}$-symmetric systems have
recently been investigated in the context of fluid dynamics (Fu and Qin,
2023). $\mathcal{PT}$ and anti-$\mathcal{PT}$ symmetries have only recently
been taken into account in fluid systems, to directly constrain fluid flows
and instabilities by spontaneous symmetry breaking (Qin et al., 2019; Fu and
Qin, 2020b; David et al., 2022). The $\mathcal{PT}$ symmetry framework has
been demonstrated to be particularly adapted to tackle high-dimensional
parameter space (e.g. the drift wave instabilities in tokamaks, Qin et al.
2021). The richness underlying $\mathcal{PT}$-symmetric systems has thus been
demonstrated unequivocally.

In this study, we aim to show that $\mathcal{PT}$ and anti-$\mathcal{PT}$
symmetries also control some of the properties of linear modes of
astrophysical systems, such as instabilities or absence of propagation,
stability exchange, conservation of quadratic quantities and this, for
physical processes of different origins. We first define $\mathcal{PT}$ and
anti-$\mathcal{PT}$ symmetries in Sect. 2. We then study how they condition
the outcome of the linear analysis of four canonical astrophysical systems,
namely the self-gravity instability in discs (Sect. 3, wave propagation in
gas-dust mixtures (Sect. 6), and buoyancy instability in stellar (Sect. 4) and
planetary interiors (Sect. 5). We then discuss these analyses in a broader
perspective for symmetric systems in Sect. 7.

##  2 Definitions

###  2.1 $\mathcal{PT}$ symmetry

The evolution of a monochromatic linear perturbation to the steady state of a
fluid can be written under the form of an eigenvalue problem ${H(a)X=\omega
X}$, where the operator $H$ that describes the evolution of the perturbation
depends continuously on some real parameter $a$. The complex eigenvalues are
generically denoted by the frequency $\omega$, which can be a complex number a
priori. This problem is said to be $\mathcal{PT}$-symmetric if there exists a
unitary operator $U$ such that

$UHU^{-1}=H^{*}\ ,$ (1)

where ∗ stands for complex conjugation. In this case, the eigenvalues go in
pairs: if $\omega$ is an eigenvalue of the problem, so is $\omega^{*}$.
Indeed, $HX=\omega X$ implies $H(U^{-1}X^{*})=\omega^{*}(U^{-1}X^{*})$. In
$\mathcal{PT}$-symmetric systems, unstable modes are therefore accompanied by
damped modes (imaginary parts of opposite sign). When $\omega\neq\omega^{*}$,
$\mathcal{PT}$ symmetry is said to be spontaneously broken by $H$. By defining
the symmetry operator $\mathcal{PT}\equiv U\Theta$ where $\Theta$ denotes
complex conjugation, one has $[\mathcal{PT},H]=0$. The two operators commute,
but do not share $X$ as an eigenvector for non-real eigenvalues.
Mostafazadeh (2002a, b); Zhang et al. (2020c) showed that $\mathcal{PT}$
symmetry is equivalent to another symmetry, pseudo-Hermiticity, for finite-
dimensional problems. A problem described by the operator $H$ is pseudo-
Hermitian if there exists an Hermitian operator $V$ such that

$VHV^{-1}=H^{\dagger},$ (2)

where † denotes the Hermitian conjugate with respect to the appropriate
Hermitian scalar product of the problem, which we will write $\cdot$ in the
following such that $Y_{1}\cdot Y_{2}\equiv Y_{1}^{*\top}Y_{2}$.
Let $X_{1}$ and $X_{2}$ be two eigenvectors of $H$ satisfying Eq.(2). The
eigenvalues $\omega_{1}$ and $\omega_{2}$ may take the same value: they
degenerate. When the eigenvectors also become identical, the degeneracy is
called an Exceptional Point (noted EP in the following). The linear stability
of these systems is described by Krein theory (Kreĭn, 1950; Kirillov, 2013):
in such systems, the transition between stability and instability necessarily
involves the exceptional degeneracy of stable eigenvalues. This phenomenon is
known as a Krein collision. Before the collision, the system is stable,
$\omega_{1}$ and $\omega_{2}$ are real, $\mathcal{PT}$ symmetry is said to be
unbroken and the Krein quantities defined as $X_{1,2}\cdot VX_{1,2}$ are non-
zero and of opposite signs. After the collision, the system is unstable, the
complex frequencies $\omega_{1}=\omega_{2}^{*}$ are complex conjugates of each
other, and the two Krein quantities are necessarily zero. The Krein quantity
$Y\cdot VY$ is a constant of motion of a general solution $Y(t,a)$ of the
evolution equation $i\partial_{t}Y=H(a)Y$. Indeed,

$\displaystyle i\partial_{t}(Y\cdot VY)$ $\displaystyle=$ $\displaystyle
Y\cdot(-H^{\dagger}V)Y+Y\cdot VHY$ (3) $\displaystyle=$ $\displaystyle
Y\cdot(VH-H^{\dagger}V)Y$ (4) $\displaystyle=$ $\displaystyle 0.$ (5)

Hence,

$Y\cdot VY=\mathrm{cst}.$ (6)

Moreover, for a Fourier mode $Y=\mathrm{e}^{-i\omega t}X$, this conserved
quantity is exactly zero for unstable and damped modes, and non-zero for
propagative modes. This is shown by expressing

$\displaystyle X\cdot VHX=X\cdot H^{\dagger}VX=\omega^{*}\;X\cdot VX$ (7)
$\displaystyle=\omega\;X\cdot VX.$ (8)

According to the identity above, a situation where $\omega\neq\omega^{*}$
implies that $X\cdot VX=0$.

###  2.2 Anti-$\mathcal{PT}$ symmetry

The eigenvalue problem $H\left(a\right)X=\omega X$ is said to be
anti-$\mathcal{PT}$-symmetric, or alternatively $\mathcal{CP}$-symmetric, if
there exists a unitary operator $\tilde{U}$ such that

$\tilde{U}H\tilde{U}^{-1}=-H^{*}.$ (9)

In this case, the eigenvalues come in pairs: $\omega$ and $-\omega^{*}$. They
have equal imaginary parts, and opposite real parts and describe counter-
propagating waves. In particular, $H$ being anti-$\mathcal{PT}$-symmetric is
equivalent for $iH$ to be $\mathcal{PT}$-symmetric. The same equivalence
between $\mathcal{PT}$ symmetry and pseudo-Hermiticity then applies, and
guarantees that an anti-$\mathcal{PT}$-symmetric operator $H$ is also pseudo-
chiral: there exist an Hermitian operator $\tilde{V}$ such that

$\tilde{V}H\tilde{V}^{-1}=-H^{\dagger}.$ (10)

The quantity $X\cdot\tilde{V}X$ is not strictly a Krein quantity, and is not
in general a constant of motion of a solution of $i\partial_{t}X=HX$. However,
for a Fourier mode $Y=\mathrm{e}^{-i\omega t}X$, it satisfies

$Y\cdot\tilde{V}Y=\mathrm{e}^{2\mathrm{Im}(\omega)t}\;X\cdot\tilde{V}X,$ (11)

and is constant to zero for modes with $\mathrm{Re}(\omega)\neq 0$, i.e for
modes that spontaneously break the anti-$\mathcal{PT}$ symmetry. Indeed,

$\displaystyle
X\cdot\tilde{V}HX=X\cdot(-H^{\dagger}\tilde{V})X=-\omega^{*}\;X\cdot\tilde{V}X$
(12) $\displaystyle=\omega\;X\cdot\tilde{V}X.$ (13)

One concludes that whenever $\omega\neq-\omega^{*}$, $X\cdot\tilde{V}X=0$: the
Krein quantity is zero at all times. An analogue of a Krein collision may
occur. When the anti-$\mathcal{PT}$ symmetry is unbroken, their exist a pair
of eigenvalues of $H$ denoted $\omega_{1}$ and $\omega_{2}$ which are purely
imaginary.

Consider the situation in which the parameter $a$ passes through an
Exceptional Point, so that after the degeneracy we have
$\omega_{1}=-\omega_{2}^{*}$ with non-zero real parts. Then, the
anti-$\mathcal{PT}$ symmetry is spontaneously broken. This phenomenon
corresponds to a Krein collision for the operator $iH$. It follows that the
quantity $X\cdot\tilde{V}X$ is non-zero and not conserved for non-propagative
modes, and strictly zero for propagating modes for all times.
Systems that are both $\mathcal{PT}$ and anti-$\mathcal{PT}$ symmetric are
called bi-symmetric.

The eigenvalues of bi-symmetric systems take the form of $(\omega,-\omega)$
pairs of real or purely imaginary frequencies. When $\omega$ is real, the
associated eigenvector spontaneously breaks the anti-$\mathcal{PT}$ symmetry.
Conversely, when $\omega$ is imaginary, the associated eigenvector
spontaneously breaks $\mathcal{PT}$ symmetry. Thus, one of those two
symmetries is necessarily broken for any value of the parameters $a$, and
either $X\cdot VX=0$ (broken $\mathcal{PT}$) or $X\cdot\tilde{V}X=0$ (broken
anti-$\mathcal{PT}$). In bi-symmetric systems, $X\cdot VX$ is always a
constant of motion, while $X\cdot\tilde{V}X$ is constant only when it is zero,
in the spontaneously broken anti-$\mathcal{PT}$ symmetry phase.

##  3 Gravitational instabilities

###  3.1 Jeans instability

The stability of self-gravitating objects plays a central role in the
formation of structures throughout the Universe. The simplest problem to study
would be the Jeans instability of a collapsing sphere (Jeans, 1902). We
consider local radial adiabatic perturbations of a 3D homogeneous self-
gravitating, non-rotating, inviscid sphere of constant density $\rho_{0}$ and
sound speed $c_{\rm s}$. Linear perturbations of mass and momentum
conservation as well as Poisson equation form a 4$\times$4 system that can be
written after Fourier transforms in time and space

$\omega\begin{pmatrix}\rho^{\prime}/\rho_{0}\\\\[5.0pt]
<EMAIL_ADDRESS>\<EMAIL_ADDRESS><EMAIL_ADDRESS><EMAIL_ADDRESS>(14)

This problem is invariant under both parity and time reversal. It is therefore
trivially $\mathcal{PT}$-symmetric, as can be proved by the statement
$H_{\mathrm{J}}=H_{\mathrm{J}}^{*}$. As such, $H_{\mathrm{J}}$ is pseudo-
Hermitian, with the symmetry operator
$V=\mathrm{diag}\left(\left(1-\frac{4\pi\mathcal{G}\rho_{0}}{c_{\mathrm{s}}^{2}|k|^{2}}\right),1,1,1\right).$
The Krein quantity associated to the Krein collision of the Jeans instability
is then $X\cdot
VX=\left(1-\frac{4\pi\mathcal{G}\rho_{0}}{c_{\mathrm{s}}^{2}|k|^{2}}\right)|\rho^{\prime}|^{2}+|\@vec{v}^{\prime}|^{2}$.
This quantity is non-zero for propagating waves, and zero for the unstable
modes. Complementarily, $H_{\mathrm{J}}$ is anti-$\mathcal{PT}$ symmetric and
pseudo-chiral. Indeed, it satisfies Eq.10 for
$\tilde{V}=\mathrm{diag}(-\left(1-\frac{4\pi\mathcal{G}\rho_{0}}{c_{\mathrm{s}}^{2}|k|^{2}}\right),1,1,1)$.
The quantity
$X\cdot\tilde{V}X=-\left(1-\frac{4\pi\mathcal{G}\rho_{0}}{c_{\mathrm{s}}^{2}|k|^{2}}\right)|\rho^{\prime}|^{2}+|\@vec{v}^{\prime}|^{2}$
is then zero for propagative modes and non-constant and non-zero for unstable
modes. The Krein collision occurs at an Exceptional Point for which $H$ cannot
be diagonalised. This EP corresponds to marginal stability, which is reached
at the Jeans wavenumber
$|\@vec{k}|=k_{\mathrm{J}}\equiv\sqrt{4\pi\mathcal{G}\rho_{0}}/c_{\mathrm{s}}$.
Longer wavelengths are unstable, and shorter wavelengths are propagating.

###  3.2 Toomre instability

More interesting is the case of a self-gravitating disc (Goodman and Narayan,
1988; Goodman, 2003; Bertin and Lodato, 1999). In this system, $\mathcal{P}$
and $\mathcal{T}$ are broken individually, but the overall system of linear
perturbations is still $\mathcal{PT}$-symmetric. Stabilization of large and
small scales by rotation and pressure can be sufficient to stabilize
astrophysical discs against gravitational collapse (Toomre, 1964). We discuss
here the simplest case of axisymmetric perturbations of short radial
wavelengths evolving in a razor-thin Keplerian disc. The stability of linear
perturbations is given by the celebrated Toomre criterion:
$Q\equiv\frac{c_{\mathrm{s}}\kappa}{\pi\mathcal{G}\Sigma_{0}}>1$, where
$\kappa$ denotes the epicyclic frequency of the disc and $\Sigma_{0}$ its
surface density. The razor-thin disc is integrated vertically (Armitage,
2010), to give the following set of radial linear perturbation

$\displaystyle H_{\mathrm{sgd}}\begin{pmatrix}h^{\prime}\\\ v_{r}^{\prime}\\\
v_{\phi}^{\prime}\end{pmatrix}=\left(\begin{array}[]{ccc}0&c_{\mathrm{s}}k&0\\\
-\frac{k}{c_{\mathrm{s}}}\left(\frac{2}{Q|k|}-1\right)&0&2i\\\
0&-\frac{i}{2}&0\\\ \end{array}\right)\begin{pmatrix}h^{\prime}\\\
v_{r}^{\prime}\\\
v_{\phi}^{\prime}\end{pmatrix}=\omega\begin{pmatrix}h^{\prime}\\\
v_{r}^{\prime}\\\ v_{\phi}^{\prime}\end{pmatrix},$ (18)

where we used perturbations of the enthalpy $h^{\prime}$ and the two
horizontal velocities $v_{r}^{\prime}$ and $v_{\phi}^{\prime}$ as variables,
and the orbital time $\Omega^{-1}$ and the pressure length
$c_{\mathrm{s}}\Omega^{-1}$ as units of time and length.

With $U=\mathrm{diag}(1,1,-1)$, this problem satisfies Eq. (1) and is
therefore $\mathcal{PT}$-symmetric. Rotation breaks the reflection symmetry in
the azimuthal direction and the time reversal symmetry, but not the
combination of both. It is expected that the system is pseudo-Hermitian, which
can be shown explicitly. With

$V=\left(\begin{array}[]{ccc}\frac{1}{c_{\mathrm{s}}^{2}}&0&-\frac{4ik}{c_{\mathrm{s}}Q|k|}\\\
0&1&0\\\ \frac{4ik}{c_{\mathrm{s}}Q|k|}&0&4\left(\frac{2|k|}{Q}+1\right)\\\
\end{array}\right),$ (20)

the system satisfies Eq. (2). The Krein quantity associated with this
$\mathcal{PT}$ symmetry provides the following energy partition, which applies
to the unstable mode

$\begin{split}X\cdot
VX=0=\frac{|h^{\prime}|^{2}}{c_{\mathrm{s}}^{2}}+|v_{r}^{\prime}|^{2}&+4|v_{\phi}^{\prime}|^{2}\left(1+\frac{2\pi\mathcal{G}\Sigma_{0}}{\Omega^{2}}|k|\right)\\\
&-\frac{8\pi\mathcal{G}\Sigma_{0}}{c_{\mathrm{s}}^{2}\Omega}\text{sgn}(k)\mathrm{Im}\left(hv_{\phi}^{*}\right),\end{split}$
(21)

where $k$ denotes the wavenumber in dimensional form. We illustrate this
result by computing numerically the value of this quantity through a Krein
collision, since $H_{\rm sgd}$ has degenerated eigenvalues when
$Q=\frac{2|k|}{1+k^{2}}$. These are the positions of Exceptional Points for
$H_{\mathrm{sgd}}$, which corresponds to positions of a Krein collision. The
curve of Exceptional Points separates two regions in $(Q,k)$\- space: one is
the stable region, the other is the unstable region. Figure 1 shows the EPs
curve and the Krein collision that occurs when crossing this curve.

Figure 1: Top panel : curve of Exceptional Points corresponding to marginal
stability. Maximum is reached for $Q=1$ as expected. Bottom panel : Krein
collision of the stability of the Toomre problem when varying the parameters
along the black line of the top panel. Depending on the value of $Q$, exactly
one of the two Krein quantities identified is zero, depending on which
symmetry is broken by the perturbation. The curve in green is $X\cdot VX$,
orange is $X\cdot\tilde{V}X$, solid blue is $\mathrm{Re}(\omega)$ and dotted
blue is $\mathrm{Im}(\omega)$, for a destabilizing mode.

In addition to being $\mathcal{PT}$-symmetric, $H_{\rm sgd}$ is also
anti-$\mathcal{PT}$-symmetric since

$\displaystyle\tilde{U}H_{\rm sgd}\tilde{U}^{-1}=-H_{\rm sgd}^{*},$ (22)

with $\tilde{U}=\mathrm{diag}(1,-1,-1)$. The anti-$\mathcal{PT}$ symmetry is
complementary to the $\mathcal{PT}$ symmetry, and the combination of the two
symmetries implies that the eigenvalues always come in pairs: they are always
opposite $(\omega,-\omega)$, either real or purely imaginary. They are never
general complex numbers, a generic property for bi-symmetric systems, since
one of the two symmetries is necessarily spontaneously broken.
One has

$\displaystyle\tilde{V}H_{\mathrm{sgd}}\tilde{V}^{-1}$ $\displaystyle=$
$\displaystyle-H_{\mathrm{sgd}}^{\dagger},$ (23) $\displaystyle\tilde{V}$
$\displaystyle\equiv$
$\displaystyle\left(\begin{array}[]{ccc}\frac{1}{c^{2}}&\;\;0&-\frac{4ik(Q|k|-1)}{cQ|k|}\\\
0&\;\;1&0\\\ \frac{4ik(Q|k|-1)}{cQ|k|}&\;\;0&-\frac{8k^{2}}{Q|k|}+8k^{2}-4\\\
\end{array}\right).$ (27)

Equation (27) explicitly shows that $H_{\mathrm{sgd}}$ is pseudo-chiral as
expected, and Krein theory then applies. When $iH_{\mathrm{sgd}}$ has real
eigenvalues, anti-$\mathcal{PT}$ symmetry is unbroken, and the Krein
quantities $X\cdot\tilde{V}X$ are non-zero. This holds when
$\mathcal{PT}$-symmetry is spontaneously broken, i.e in the unstable region of
parameters, where the first Krein quantity $X\cdot VX$ is zero. Conversely, in
the stable range, anti-$\mathcal{PT}$ symmetry is broken and $\mathcal{PT}$
symmetry is unbroken. In this case, $X\cdot\tilde{V}X$ is zero and $X\cdot VX$
is not equal to zero. The Krein quantity $X\cdot\tilde{V}X$ associated to the
anti-$\mathcal{PT}$-symmetry gives the following energy distribution for
propagating waves in dimensional form

$\begin{split}X\cdot\tilde{V}X=&\frac{|h^{\prime}|^{2}}{c_{\mathrm{s}}^{2}}+|v_{r}^{\prime}|^{2}+|v_{\phi}^{\prime}|^{2}\left(-4+8\left(\frac{c_{\mathrm{s}}^{2}k^{2}}{\Omega^{2}}-\frac{\pi\mathcal{G}\Sigma_{0}}{\Omega^{2}}|k|\right)\right)\\\
&-\frac{8\pi\mathcal{G}\Sigma_{0}}{c_{\mathrm{s}}^{2}\Omega}\text{sgn}(k)\left(\frac{c_{\mathrm{s}}^{2}}{\pi\mathcal{G}\Sigma_{0}}|k|-1\right)\mathrm{Im}\left(hv_{\phi}^{*}\right).\end{split}$
(28)

Hence, for each mode, either $X\cdot VX$ or $X\cdot\tilde{V}X$ is necessarily
zero: the former is zero for unstable modes and the latter is zero for
propagating waves. Both are exactly zero at the degeneracy, i.e. for
$\omega=0$, which is true for $Q=\frac{2|k|}{1+k^{2}}$. Figure 1 shows this
exchange in non-zero quantities when crossing the Krein collision. Krein
quantities simply provide energy partitions of the system in both regimes.

The matrices $V$ and $\tilde{V}$ identified above are particular cases of more
general families of matrices

$V=\left(\begin{array}[]{ccc}\frac{2c_{\rm
s}\left(-2k^{2}+\left(-1+k^{2}\right)Q\left|k\right|\right)a+k\left(-2+Q\left|k\right|\right)c}{4c_{\rm
s}^{2}kQ\left|k\right|}&b\frac{\left(kQ-2\mathrm{sgn{k}}\right)}{2c_{\rm
s}Q}&ia\\\ b\frac{\left(kQ-2\mathrm{sgn{k}}\right)}{2c_{\rm
s}Q}&\frac{1}{4}\left(2c_{\rm s}ka+c\right)&ib\\\ -ia&-ib&c\\\
\end{array}\right),$ (29)

with $a,b,c$ any real numbers such that $\det W\neq 0$. The matrix of Eq. (20)
is recovered with $a=-\frac{4k}{c_{\rm s}Q\left|k\right|}$, $b=0$ and
$c=\frac{8\left|k\right|}{Q}+2$. Similarly,

$\tilde{V}=\left(\begin{array}[]{ccc}\frac{2c_{\rm
s}\left(-2k^{2}+\left(-1+k^{2}\right)Q\left|k\right|\right)a+k\left(-2+Q\left|k\right|\right)c}{4c_{\rm
s}^{2}kQ\left|k\right|}&b\frac{i\left(kQ-2\mathrm{sgn{k}}\right)}{2c_{\rm
s}Q}&ia\\\ -b\frac{i\left(kQ-2\mathrm{sgn{k}}\right)}{2c_{\rm
s}Q}&-\frac{1}{4}\left(2c_{\rm s}ka+c\right)&b\\\ -ia&b&c\\\
\end{array}\right)$ (30)

The matrix of Eq. (27) is recovered with
$a=-\frac{4k(Q|k|-1)}{c_{\mathrm{s}}Qk}$, $b=0$ and
$c=-\frac{8}{Q}|k|+8k^{2}-4$. The same remark holds for the different problems
adressed in this study.

##  4 Buoyancy instability in stars

Consider a non-rotating star in static equilibrium, balanced by thermal
pressure and self-gravity. The steady state is spherically symmetric, the
density is stratified and decreases towards the surface. This equilibrium can
be unstable, since a perturbation can be amplified by buoyancy if the square
of the buoyancy frequency $N^{2}$ is negative, as given by the Schwarzschild
criterion (Schwarzschild, 1906).

This problem has been revisited by Leclerc et al. (2023) in inhomogeneous
media, where the system proved to be pseudo-Hermitian and pseudo-chiral.
Pseudo-Hermiticity was used to determine the Krein quantity. We will now
complete the picture by determining the Krein quantity associated to pseudo-
chirality, and show the Krein collision in this bi-symmetric problem.
Neglecting the self-gravity of the perturbations (Cowling’s approximation
Cowling 1941), we start from Eq.5 of Leclerc et al. (2023), which is the
eigenvalue equation

$\displaystyle\omega X$ $\displaystyle=$ $\displaystyle H_{\rm b}X,$ (31)
$\displaystyle H_{\rm b}$ $\displaystyle=$
$\displaystyle\begin{pmatrix}0&0&0&L_{\ell}\\\ 0&0&|N|&K_{r}-iS\\\
0&-|N|&0&0\\\ L_{\ell}&K_{r}+iS&0&0\end{pmatrix}$ (32)

where $\omega$ is the complex eigenfrequency of the perturbation, and
$X\equiv\begin{pmatrix}{v}&{w}&{\Theta}&{p}\end{pmatrix}^{\top}$ contains the
perturbation’s horizontal velocity, radial velocity, entropy and pressure
after appropriate rescaling.
$N^{2}\equiv-g\frac{\mathrm{d}{\ln\rho_{0}}}{\mathrm{d}r}-\frac{g^{2}}{c_{\mathrm{s}}^{2}}$
is the square of the buoyancy frequency which is negative here,
$S\equiv\frac{c_{\mathrm{s}}}{2g}\left(N^{2}-\frac{g^{2}}{c_{\mathrm{s}}^{2}}\right)-\frac{1}{2}\frac{\mathrm{d}{c_{\mathrm{s}}}}{\mathrm{d}r}+\frac{c_{\mathrm{s}}}{r}$
is another characteristic frequency called the buoyant-acoustic frequency,
which quantifies the coupling between $g$-modes and $p$-modes in asteroseismic
problems, $L_{\ell}=\frac{c_{\mathrm{s}}}{r}\sqrt{\ell(\ell+1)}$ is the Lamb
frequency, $K_{r}$ represents the local radial wavenumber of the wave (see
Appendix E of Leclerc et al. (2023) for details).
This problem is both pseudo-Hermitian and pseudo-chiral, and is associated
with the two matrices

$\displaystyle V$ $\displaystyle=$
$\displaystyle\mathrm{diag}(1,\;\;\;1,-1,\;\;\;1),$ (33)
$\displaystyle\tilde{V}$ $\displaystyle=$
$\displaystyle\mathrm{diag}(1,\;\;\;1,\;\;\;1,-1).$ (34)

Krein theory can therefore be applied to determine how discrete symmetries
constrain the partition of energy. $X\cdot X$ is generally not a conserved
quantity, since $H_{\rm b}$ is not Hermitian. It grows exponentially for
unstable modes and is conserved only for stable modes. However, $X\cdot VX$ is
a constant of motion for any solution of $\partial_{t}X=H_{\rm b}X$, as it has
been shown in Sect.2.1. Moreover, this constant is exactly zero only for
unstable modes. In contrast, $X\cdot\tilde{V}X$ also grows exponentially for
unstable modes, but is constantly equal to zero for stable modes. As such,

$\displaystyle X$ $\displaystyle\cdot X$
$\displaystyle=|v|^{2}+|w|^{2}+|\Theta|^{2}+|p|^{2}$
$\displaystyle\propto\mathrm{e}^{2\mathrm{Im}(\omega)t},$ (35) $\displaystyle
X$ $\displaystyle\cdot VX$
$\displaystyle=|v|^{2}+|w|^{2}-|\Theta|^{2}+|p|^{2}$
$\displaystyle=\mathrm{cst}\;\delta_{\mathrm{Im}(\omega),0},$ (36)
$\displaystyle X$ $\displaystyle\cdot\tilde{V}X$
$\displaystyle=|v|^{2}+|w|^{2}+|\Theta|^{2}-|p|^{2}$
$\displaystyle\propto\mathrm{e}^{2\mathrm{Im}(\omega)t}\delta_{\mathrm{Re}(\omega),0},$
(37)

where $\delta$ is the Kronecker delta.
For the radial modes $L_{\ell}=0$, the modes degenerate over the ring
parameterised by $K_{r}^{2}+S^{2}=-N^{2}>0$. These parameter values are
Exceptional Points, as $H_{\rm b}$ can no longer be diagonalized there. The
ring separates two ranges of parameters, the stable and the unstable part.
Figure 2 shows the Krein quantities, over the ring of EPs in which the
collision takes place. The Krein collision behaves as expected: inside the
ring, the system is unstable, and $X\cdot VX$ is zero and $X\cdot\tilde{V}X$
is non-zero. Outside the ring, the modes are stable and propagate as sound
waves, and $X\cdot VX$ is non-zero and $X\cdot\tilde{V}X$ is zero.

This structure is equivalent to the bi-symmetric problem of self-gravitating
discs discussed in Sect.3, and shown on Fig.1, even if the physical situation
is very different. The approach followed in this problem is the reverse of
that followed in section 3 for the Toomre problem: it is easier to first
establish that the system is pseudo-Hermitian and pseudo-chiral, since the
matrices $V$ and $\tilde{V}$ have $\pm 1$ on the diagonal.

Figure 2: Same as Fig. 1 for the Schwarzschild instability. Non-hermiticity is
broken in this case by buoyancy instead of self-gravity.

##  5 Rayleigh-Bénard instability

Not all instabilities correspond to spontaneous symmetry breakings. The
Rayleigh-Bénard instability (Bénard, 1900; Rayleigh, 1916) is such a case.
Interiors of rocky planets are often treated as incompressible with both
significant viscous and thermal diffusion (e.g. Bergé and Dubois 1984;
Bodenschatz et al. 2000; Brandenburg 2021). This results in a regime of
buoyancy instability that is controlled by diffusion effects and differs
fundamentally from that discussed in Sect. 4 for adiabatic perturbations in
compressible stars or planetary atmospheres. In the Rayleigh-Bénard regime,
the instability criterion is given by the value of the Rayleigh number
$\mathrm{Ra}=g\alpha\beta L^{4}/\nu\kappa$, which needs to be greater than
some critical value of order $10^{2}-10^{3}$ which depends on the boundary
conditions. This famous problem admits analytical solutions for an homogeneous
background configuration with rigid-lids isothermal boundary conditions, which
provides an illustrative example in this study for discussing the symmetries
involving diffusive effects.
The basic set of equations consists of the incompressibility condition, a
Navier-Stokes equation, and thermal heat diffusion for a 2D $(x,z)$ fluid (see
Appendix A). The $x$\- direction is invariant by translation and infinite, the
$z$\- direction is of length $L$, along which a constant temperature gradient
$\partial_{z}T_{0}=-\beta$ is directed. With dimensionless variables, the
problem for Fourier modes $\exp{(-i\omega t+ik_{z}z+ik_{x}x)}$ can be
converted into the symmetrical form

$\displaystyle\omega\begin{pmatrix}\tilde{v_{z}}\\\\[5.69046pt]
\tilde{T}\end{pmatrix}=H_{\rm RB}\begin{pmatrix}\tilde{v_{z}}\\\\[5.69046pt]
\tilde{T}\end{pmatrix},$ (38) $\displaystyle H_{\rm
RB}=-i\begin{pmatrix}k^{2}&-\frac{k_{x}}{k}\left(\frac{\mathrm{Ra}}{\mathrm{Pr}}\right)^{1/2}\\\
-\frac{k_{x}}{k}\left(\frac{\mathrm{Ra}}{\mathrm{Pr}}\right)^{1/2}&k^{2}/\mathrm{Pr}\end{pmatrix},$
(39)

where $\mathrm{Pr}=\nu/\kappa$ is the Prandtl number and
$k=\sqrt{k_{x}^{2}+k_{z}^{2}}$.
This problem is anti-$\mathcal{PT}$-symmetric and pseudo-chiral, since one
reads that

$H_{\rm RB}=-H_{\rm RB}^{*}=-H_{\rm RB}^{\dagger},$ (40)

for $\mathrm{Ra}\geq 0$. In this case, the anti-$\mathcal{PT}$ symmetry is
protected and cannot be spontaneously broken since $H_{\rm RB}$ is anti-
Hermitian (i.e. $iH_{\rm RB}$ is Hermitian). Indeed, the spectral theorem
guarantees that the eigenvalues are purely imaginary and the eigenvectors are
orthogonal to each other.
The location of EPs is given by

$\mathrm{Ra}=\left(1-\frac{(1+\mathrm{Pr})^{2}}{\mathrm{4Pr}}\right)\frac{(k_{x}^{2}+k_{z}^{2})^{3}}{k_{x}^{2}},$
(41)

which requires $\mathrm{Ra}$ to be negative. This condition can only be
fulfilled for an inverted temperature gradient, which corresponds to a stably
stratified liquid. In this case, the anti-$\mathcal{PT}$ symmetry is no longer
protected and can be broken spontaneously, which would lead to propagating,
damped internal waves.
Anti-Hermiticity ensures that the eigenvalues $\omega=i\eta$ are purely
imaginary for all $k_{x},k_{z},\mathrm{Ra},\mathrm{Pr}$. The transition to
instability is given by the change of sign $\eta<0$ to $\eta>0$. This
necessarily occurs for $\omega=0$, providing the marginal stability criterion
for the Rayleigh-Bénard instability. This condition is satisfied when

$\mathrm{Ra}=\frac{(k_{x}^{2}+k_{z}^{2})^{3}}{k_{x}^{2}},$ (42)

whose minimum for each $k_{x}$,$k_{z}$ compatible with the boundary conditions
gives the critical value of the Rayleigh number (e.g.
$\mathrm{Ra_{c}}=\frac{27}{4}\pi^{4}$ for stress-free, impenetrable and
isothermal boundaries). Figure 3 shows the regions of this problem.

Figure 3: Different regions of the Rayleigh-Bénard problem in the parameter
space. In the upper part $\mathrm{Ra}>0$, anti-$\mathcal{PT}$ symmetry is
protected, frequencies are purely imaginary, the system exhibits an
instability (black dashed) and no EPs. In the lower part $\mathrm{Ra}<0$,
anti-$\mathcal{PT}$ is no longer protected by anti-Hermiticity, and the system
exhibits a curve of EP2s (red), delimiting a region where perturbations
propagate. It is a region where anti-$\mathcal{PT}$ symmetry is spontaneously
broken.

##  6 Waves in dusty mixtures

Mixtures of pressureless dust and inviscid non-magnetised gas are generically
used to model basic hydrodynamical properties of the dusty interstellar medium
(e.g. Saffman 1962; Baines et al. 1965; Ahuja 1973; Gumerov et al. 1988).
Exchange of momentum between gas and dust is modeled as a drag force
proportional to the difference in their velocities. The dustywave problem
consist of studying the propagation of a 1D acoustic wave in an homogeneous
dusty medium. The question is whether small disturbances propagate or are
simply attenuated, since both regimes are possible (e.g. Laibe and Price 2011;
David-Cléris and Laibe 2021). The mixture is initially homogeneous and at
rest. Looking for small perturbative solutions as Fourier modes
$\exp{(i(\omega t-kx))}$, one obtains the eigenvalue equation

$\displaystyle H_{\rm dw}X=\omega X,$ (43) $\displaystyle H_{\rm
dw}\equiv\begin{pmatrix}0&\;\;\;0&\;\;\;\;k&0\\\ 0&\;\;\;0&\;\;\;\;0&k\\\
k&\;\;\;0&\;\;\;\;i\epsilon&\;-i(1-\epsilon)\\\
0&\;\;\;0&\;-i\epsilon&\;\;\;i(1-\epsilon)\end{pmatrix},$ (44)

after appropriate rescaling and choice of dimensionless parameters given in
App. B. The operator $H_{\rm dw}$ formally depends on two parameters: the
wavenumber $k$ and the dust density fraction of the mixture $\epsilon$.
$H_{\rm dw}$ satisfies Eq. (9) with $\tilde{U}=\mathrm{diag}(1,1,-1,-1)$,
meaning that $H_{\rm dw}$ is anti-$\mathcal{PT}$-symmetric. This property
stems from the reflection symmetry $x\mapsto-x$ of the fluid in physical
space, as it has been shown generically by David et al. (2022). Hence, the
eigenvalues of $H_{\rm dw}$ are either imaginary or consist of pairs of
complex numbers with opposite real parts. In parameter space, the transition
between propagating and non-propagative regions must then be an Exceptional
Points at which the solution of two counter-propagating waves degenerates.

The characterization of these EPs is carried out according to the procedure
described in Delplace et al. (2021) (see App. C for details). Exceptional
Points of order 2 (or EP2s) are located on a manifold of dimension 1 in the
parameter space ($k,\epsilon$), which is a curve in a plane. Values of
($k,\epsilon$) for which an eigenvalue has a multiplicity of two are
characterized by the fact that the characteristic polynomial $P(X)$ of $H_{\rm
dw}$ has a root of multiplicity two. Such a root is therefore also a root of
its derivative $P^{\prime}(X)$, which means that the resultant
$R_{1}\equiv\mathcal{R}(P,P^{\prime})$ between the two polynomials is 0. The
resultant between two polynomials only cancels out if they have a common root.
Since $H_{\rm dw}$ is anti-$\mathcal{PT}$-symmetric, $R_{1}$ is a real
quantity. The values of ($k,\epsilon$), for which $R_{1}=0$, therefore define
a curve. We calculate

$R_{1}=-4k^{6}(\epsilon-1)^{2}\bigg{(}1+k^{4}-\epsilon+k^{2}\big{(}2+9\epsilon(-1+\frac{3}{4}\epsilon)\big{)}\bigg{)},$
(45)

and show on Fig.4 the curve $R_{1}=0$ (top panel, red). This curves provide a
simple alternative derivation of the result of David-Cléris and Laibe (2021).
Non-propagating waves exist for $\epsilon>8/9$, which is indeed the minimum of
the curve of EPs.

Interestingly, in addition to the EP2s, the above method reveals the existence
of Exceptional Points of order 3 for which three eigenvectors merge. These
EP3s form a fold of dimension 0, which are two points in ($k,\epsilon$) space.
These two points are located at

$k=\pm\frac{1}{\sqrt{3}},\quad\epsilon=\frac{8}{9},$ (46)

which are the points at which the curve of the EP2s loses its regularity. At
these EP3s, $H_{\rm dw}$ has only two eigenvectors corresponding to the two
eigenvalues 0 and $\frac{-i}{3}$. We note that these EP3s carry topological
charges, also called winding numbers $W_{3}$,

$W_{3}\Big{|}_{\mathrm{EP3}_{1}}=+1,\quad W_{3}\Big{|}_{\mathrm{EP3}_{2}}=-1,$
(47)

where $\mathrm{EP3}_{1}$ is the point
$(k,\epsilon)=(+\frac{1}{\sqrt{3}},\frac{8}{9})$ and $\mathrm{EP3}_{2}$ is the
point $(k,\epsilon)=(-\frac{1}{\sqrt{3}},\frac{8}{9})$ (see App. D). The
question of connecting these topological charges to edge modes or particular
global modes remains open (e.g. Leclerc et al. 2023).

According to the discussion in Sect. 2.2, $H_{\rm dw}$ must be pseudo-chiral.
The operator $\tilde{V}$ can be found explicitly

$\displaystyle\tilde{V}$ $\displaystyle=$
$\displaystyle\left(\begin{array}[]{cccc}\epsilon&\epsilon-1&ik&0\\\
\epsilon-1&\frac{(\epsilon-1)^{2}}{\epsilon}&0&-\frac{ik(\epsilon-1)}{\epsilon}\\\
-ik&0&0&0\\\ 0&\frac{ik(\epsilon-1)}{\epsilon}&0&0\\\ \end{array}\right).$
(52)

As has been shown, pseudo-chirality implies that if $\omega\neq-\omega^{*}$,
the quantity $X\cdot\tilde{V}X$ must be zero. This condition is fulfilled when
sound propagates. Although the system is not conservative, it still possesses
an associated Krein quantity $X\cdot\tilde{V}X$

$\displaystyle X\cdot\tilde{V}X$ $\displaystyle=$
$\displaystyle\frac{1-\epsilon}{\epsilon}|\rho_{\mathrm{d}}|^{2}+\frac{\epsilon}{1-\epsilon}|\rho_{\mathrm{g}}|^{2}$
$\displaystyle-2\left(\mathrm{Re}(\rho_{\mathrm{g}}\rho_{\mathrm{d}})-k\frac{\rho_{0}}{c_{\mathrm{s}}}\mathrm{Im}(\rho_{\mathrm{g}}v_{\mathrm{g}}^{*}+\rho_{\mathrm{d}}v_{\mathrm{d}}^{*})\right),$

where $\rho_{0}=\rho_{\mathrm{g},0}+\rho_{\mathrm{d},0}$ is the total
background density. $X\cdot\tilde{V}X$ must be zero for propagating waves and
non-zero for non-propagating waves. Figure 4 shows the numerical confirmation
of this result. $X\cdot\tilde{V}X$ is only a constant of motion if it is zero,
i.e. only for sound waves. Another technique, where the anti-$\mathcal{PT}$
operator is diagonalized, can be used to extract constants of motions in the
unbroken phase, as it is done on the Kelvin-Helmholtz instability by Qin et
al. (2019) (see App. E for details).
Finally, it should be noted that in a monofluid description of the mixture,
the variables are $\rho=\rho_{\mathrm{g}}+\rho_{\mathrm{d}}$,
$v=\rho_{\mathrm{g}}v_{\mathrm{g}}+\rho_{\mathrm{d}}v_{\mathrm{d}}$,
$\frac{\rho_{\mathrm{d}}}{\rho_{\mathrm{g}}}$ and $\Delta
v=v_{\mathrm{d}}-v_{\mathrm{g}}$ (Laibe and Price, 2014). If the perturbed
quantities are denoted by ′, the Krein quantity that cancels out for
propagating perturbations is

$\begin{split}&X\cdot\tilde{V}_{\mathrm{mono}}X=\\\
&\frac{(1-\epsilon)^{3}}{\epsilon}\left|\left(\frac{\rho_{\mathrm{d}}}{\rho_{\mathrm{g}}}\right)^{\prime}\right|^{2}+2k\mathrm{Im}\left((1-\epsilon)^{2}\left(\frac{\rho_{\mathrm{d}}}{\rho_{\mathrm{g}}}\right)^{\prime}\Delta
v^{\prime*}+\frac{\rho^{\prime}}{\rho_{0}}v^{\prime*}\right).\end{split}$

It is associated with the symmetry operator

$\tilde{V}_{\mathrm{mono}}=\left(\begin{array}[]{cccc}0&0&-ik&0\\\
0&\frac{(\epsilon-1)^{3}}{\epsilon}&0&-ik(\epsilon-1)^{2}\\\ ik&0&0&0\\\
0&ik(\epsilon-1)^{2}&0&0\\\ \end{array}\right).$ (55)

Figure 4: Same as Fig. 1 for the Dustywave problem. For ${k=0.4}$, the Krein
quantity $X\cdot\tilde{V}X$ is represented as a function of $\epsilon$. This
Krein quantity provides a quadratic constant of motion for dusty sound waves
in the propagating phase (low $\epsilon$).

Squire and Hopkins (2018) have shown that when grains are additionally
streaming though the gas, a generic resonant drag instability develops in the
mixture. The system is neither $\mathcal{PT}$ nor anti-$\mathcal{PT}$
symmetric and Krein theory can not be applied directly. The instability comes
from a more complex mechanism where a resonance occurs between an
anti-$\mathcal{PT}$ symmetric correction due to drag to a leading
$\mathcal{PT}$ symmetric perturbation that arises from the streaming (e.g.
Zhuravlev 2019; Magnan et al. 2024). Analysis of these problems from the point
of view of discrete symmetries requires further specific studies.

##  7 Gravitational instability of a dusty disc

To our knowledge, there are still no theoretical predictions about what
happens when systems with different discrete symmetries are combined. To
illustrate this, we consider a final example, drawing on the analysis of
Longarini et al. (2022), which deals with the gravitational instability of a
dusty razor-thin disc. In the context of planet formation, the aim is to
quantify the ways in which dust can favor the local collapse of gas or even
clump itself into planetary embryos. This problem is a combination of the two
problems presented in Sect. 3 and Sect. 6 and as such, is of order 6. With the
notations used above, the matrix of linear perturbation in the local shearing
box is

$\displaystyle H_{\rm sgd}=$
$\displaystyle\left(\begin{array}[]{cccccc}0&ik\Sigma_{g}&0&0&0&0\\\
-ik\left(\frac{2\pi\mathcal{G}}{|k|}-\frac{c_{g}^{2}}{\Sigma_{g}}\right)&\frac{\epsilon}{t_{s}}&-2\Omega&-\frac{2i\pi\mathcal{G}k}{|k|}&-\frac{\epsilon}{t_{s}}&0\\\
0&-2B&\frac{\epsilon}{t_{s}}&0&0&-\frac{\epsilon}{t_{s}}\\\
0&0&0&0&ik\Sigma_{d}&0\\\
-\frac{2i\pi\mathcal{G}k}{|k|}&-\frac{1}{t_{s}}&0&-ik\left(\frac{2\pi\mathcal{G}}{|k|}-\frac{c_{d}^{2}}{\Sigma_{d}}\right)&\frac{1}{t_{s}}&-2\Omega\\\
0&0&-\frac{1}{t_{s}}&0&-2B&\frac{1}{t_{s}}\\\ \end{array}\right),$ (62)

where $\Omega$ is the orbital frequency and $B$ is the local Oort parameter.
$c_{\rm g}$ and $c_{\rm d}$ denote the sound speeds of the gas and dust
respectively. Time and lengths are rescaled to the stopping time $t_{\rm s}$
and the stopping length $c_{\rm g}t_{\rm s}$. The characteristic polynomial of
the system is of order 5, so that the roots $\omega\left(k\right)$ cannot be
determined analytically. In the regime of weak drag, Longarini et al. (2022)
determines the marginal stability numerically.

On the other hand, this problem is $\mathcal{PT}$-symmetric, with the operator
${U=\mathrm{diag}(-1,1,1,-1,1,1)}$. The marginal stability curve is therefore
a curve of EPs that corresponds to the spontaneous breaking of $\mathcal{PT}$
symmetry. From the analysis above, its exact expression is

$\mathcal{R}_{\rm sgd}(P,P^{\prime})=0,$ (63)

a polynomial of high order with respect to the parameters of the problem. The
explicit expression of Eq. (63) extends over two pages and is given in the
worksheet in the Acknowledgements. Hence, from a general perspective, the
$\mathcal{PT}$-symmetric characterization of a marginal stability criterion
allows to derive directly its analytical expression, albeit cumbersomely,
without approximations relying on asymptotically strong or weak drag regimes.
Following Qin et al. (2021), this astrophysical example shows that the
analysis of discrete symmetries therefore provides a powerful way to make
analytical predictions about the marginal stability of symmetric systems with
high dimensionalities.

Figure 5: Exceptional Points of the dusty self-gravitating disk are a
consequence of the spontaneous breaking of $\mathcal{PT}$ symmetry. The
equation $R\left(P,P^{\prime}\right)=0$ is obtained analytically with a
resultant method, reproducing the results of Longarini et al. (2022).

##  8 Conclusion

For several astrophysical objects, the properties of the small perturbations
around the equilibrium are controlled by discrete symmetries, such as the
$\mathcal{PT}$ symmetry or the anti-$\mathcal{PT}$ symmetry. For example, we
show that such symmetries are relevant for

  * •

The stability of isothermal spheres and thin self-gravitating disc: the Jeans
and the are Toomre problem are bi-symmetric due to self-gravity.

  * •

The propagation of waves and the onset of convection in stratified fluids for
compressible and adiabatic perturbations: the Schwarzschild problem is bi-
symmetric due to buoyancy.

  * •

The propagation of waves and the onset of convection in stratified fluids for
incompressible and diffusive perturbations: the Rayleigh-Bénard problem is
anti-$\mathcal{PT}$-symmetric due to diffusion, as well as anti-Hermitian.

  * •

The propagation of a wave in a dust-gas mixture: the dustywave problem is
anti-$\mathcal{PT}$ symmetric due to drag.

In these systems, the parameter space can be divided into different regions,
in which the discrete symmetries are either broken or not. These regions are
separated by Exceptional Points (EPs) for which the system of linear
perturbation can no longer be diagonalised.

From a methodology perspective, the analysis of discrete symmetries allows
when relevant to determine directly physical quantities that may be harder to
obtain alternatively:

  * •

The location of the EPs in the parameter space can be determined analytically
from the symmetries without knowing the dispersion relation. This method is
particularly effective for systems with high dimensions, where the dispersion
relations may be untractable for analytical techniques, as in the case of the
dusty self-gravitating disc.

  * •

Krein theory then makes it possible to directly determine the exchange of
stability when EPs are crossed, as well as the partition of energy enforced by
Krein invariants associated with the discrete symmetries, even for dissipative
systems.

Some of these EPs can be of high order and carry topological charges, as this
is the case for the dustywave problem. Note that not all physical systems
allow such discrete symmetries. In particular, their spontaneous breaking is
not a necessary condition for the system to be linearly unstable, as the
Rayleigh-Bénard problem shows. We have focused here on the analysis of local
stability, which applies to homogeneous systems. The study of global modes for
inhomogeneous systems requires an extension of the Wigner-Weyl framework to
non-Hermitian systems (e.g. Onuki 2020; Perez et al. 2021, 2022; Perez 2022;
Leclerc et al. 2023).

## Acknowledgements

We thank T. David–Cléris, E. Lynch, C. Longarini, P. Delplace, S. Labrosse and
Y. Ricard for stimulating discussions and comments. We thank the conference
Star@Lyon for stimulation of this study. GL, AL and NP acknowledge the support
of the European Research Council (ERC) CoG project PODCAST No. 864965. AL is
funded by Contrat Doctoral Spécifique Normalien. We have used Mathematica
(Wolfram Research Inc., 2021), worksheets can be found at
https://github.com/ArmandLeclerc/PTsymInAstroWaves.

## References

  * Ahuja [1973] Avtar S Ahuja.  Wave equation and propagation parameters for sound propagation in suspensions.  _Journal of Applied Physics_ , 44(11):4863–4868, 1973.
  * Armitage [2010] Philip J. Armitage.  _Astrophysics of Planet Formation_.  Cambridge University Press, 2010.
  * Baines et al. [1965] MJ Baines, IP Williams, AS Asebiomo, and RL Agacy.  Resistance to the motion of a small sphere moving through a gas.  _Monthly Notices of the Royal Astronomical Society_ , 130(1):63–74, 1965.
  * Bénard [1900] Henri Bénard.  Les tourbillons cellulaires dans une nappe liquide.  _Revue Gen. Sci. Pure Appl._ , 11:1261–1271, 1900.
  * Bender and Boettcher [1998] Carl M. Bender and Stefan Boettcher.  Real spectra in non-hermitian hamiltonians having $pt$ symmetry.  _Phys. Rev. Lett._ , 80:5243–5246, Jun 1998.  .  URL https://link.aps.org/doi/10.1103/PhysRevLett.80.5243.
  * Bender et al. [1999] Carl M. Bender, Stefan Boettcher, and Peter N. Meisinger.  PT-symmetric quantum mechanics.  _Journal of Mathematical Physics_ , 40(5):2201–2229, May 1999.  ISSN 0022-2488.  .
  * Bender et al. [2002] Carl M. Bender, Dorje C. Brody, and Hugh F. Jones.  Complex Extension of Quantum Mechanics.  _Physical Review Letters_ , 89(27):270401, December 2002.  .
  * Bender et al. [2004] Carl M. Bender, Dorje C. Brody, and Hugh F. Jones.  Extension of $\mathcal{}PT{}$-symmetric quantum mechanics to quantum field theory with cubic interaction.  _Physical Review D_ , 70(2):025001, July 2004\.  .
  * Bergé and Dubois [1984] P Bergé and M Dubois.  Rayleigh-bénard convection.  _Contemporary Physics_ , 25(6):535–582, 1984\.
  * Bertin and Lodato [1999] G. Bertin and G. Lodato.  A class of self-gravitating accretion disks.  _A &A_, 350:694–704, October 1999.  .
  * Bodenschatz et al. [2000] Eberhard Bodenschatz, Werner Pesch, and Guenter Ahlers.  Recent developments in rayleigh-bénard convection.  _Annual review of fluid mechanics_ , 32(1):709–778, 2000.
  * Brandenburg [2021] Axel Brandenburg.  Lecture 4: Rayleigh–bénard problem.  http://norlx51.nordita.org/~brandenb/teach/AdvAstroFluids/4_Convection/notes.pdf, 2021\.
  * Cowling [1941] Thomas G Cowling.  The non-radial oscillations of polytropic stars.  _Monthly Notices of the Royal Astronomical Society, Vol. 101, p. 367_ , 101:367, 1941.
  * David et al. [2022] Tomos W. David, Pierre Delplace, and Antoine Venaille.  How do discrete symmetries shape the stability of geophysical flows?  _Physics of Fluids_ , 34(5):056605, May 2022.  ISSN 1070-6631.  .
  * David-Cléris and Laibe [2021] Timothée David-Cléris and Guillaume Laibe.  Large dust fractions can prevent the propagation of soundwaves.  _Monthly Notices of the Royal Astronomical Society_ , 504(2):2889–2894, 2021.
  * Delplace et al. [2021] Pierre Delplace, Tsuneya Yoshida, and Yasuhiro Hatsugai.  Symmetry-protected multifold exceptional points and their topological characterization.  _Physical Review Letters_ , 127(18):186602, 2021\.
  * El-Ganainy et al. [2007] R. El-Ganainy, K. G. Makris, D. N. Christodoulides, and Ziad H. Musslimani.  Theory of coupled optical PT-symmetric structures.  _Optics Letters_ , 32(17):2632–2634, September 2007.  ISSN 1539-4794.  .
  * El-Ganainy et al. [2018] Ramy El-Ganainy, Konstantinos G. Makris, Mercedeh Khajavikhan, Ziad H. Musslimani, Stefan Rotter, and Demetrios N. Christodoulides.  Non-Hermitian physics and PT symmetry.  _Nature Physics_ , 14(1):11–19, January 2018\.  ISSN 1745-2481.  .
  * Feng et al. [2013] Liang Feng, Ye-Long Xu, William S. Fegadolli, Ming-Hui Lu, José E. B. Oliveira, Vilson R. Almeida, Yan-Feng Chen, and Axel Scherer.  Experimental demonstration of a unidirectional reflectionless parity-time metamaterial at optical frequencies.  _Nature Materials_ , 12(2):108–113, February 2013\.  ISSN 1476-4660.  .
  * Feng et al. [2014] Liang Feng, Zi Jing Wong, Ren-Min Ma, Yuan Wang, and Xiang Zhang.  Single-mode laser by parity-time symmetry breaking.  _Science_ , 346(6212):972–975, November 2014\.  .
  * Feng et al. [2017] Liang Feng, Ramy El-Ganainy, and Li Ge.  Non-Hermitian photonics based on parity–time symmetry.  _Nature Photonics_ , 11(12):752–762, December 2017.  ISSN 1749-4893.  .
  * Fu and Qin [2020a] Yichen Fu and Hong Qin.  The physics of spontaneous parity-time symmetry breaking in the kelvin–helmholtz instability.  _New Journal of Physics_ , 22(8):083040, 2020a.
  * Fu and Qin [2020b] Yichen Fu and Hong Qin.  The physics of spontaneous parity-time symmetry breaking in the Kelvin–Helmholtz instability.  _New Journal of Physics_ , 22(8):083040, August 2020b.  ISSN 1367-2630.  .
  * Fu and Qin [2023] Yichen Fu and Hong Qin.  Topological modes and spectral flows in inhomogeneous PT-symmetric continuous media, September 2023.
  * Ge and Türeci [2013] Li Ge and Hakan E Türeci.  Antisymmetric pt-photonic structures with balanced positive-and negative-index materials.  _Physical Review A_ , 88(5):053810, 2013.
  * Goodman [2003] Jeremy Goodman.  Self-gravity and quasi-stellar object discs.  _Monthly Notices of the Royal Astronomical Society_ , 339(4):937–948, 2003.
  * Goodman and Narayan [1988] Jeremy Goodman and Ramesh Narayan.  The stability of accretion tori–iii. the effect of self-gravity.  _Monthly Notices of the Royal Astronomical Society_ , 231(1):97–114, 1988.
  * Gumerov et al. [1988] NA Gumerov, AI Ivandaev, and RI Nigmatulin.  Sound waves in monodisperse gas-particle or vapour-droplet mixtures.  _Journal of Fluid Mechanics_ , 193:53–74, 1988.
  * Guo et al. [2009] A. Guo, G. J. Salamo, D. Duchesne, R. Morandotti, M. Volatier-Ravat, V. Aimez, G. A. Siviloglou, and D. N. Christodoulides.  Observation of $\mathcal{}P{}\mathcal{}T{}$-Symmetry Breaking in Complex Optical Potentials.  _Physical Review Letters_ , 103(9):093902, August 2009.  .
  * Hodaei et al. [2017] Hossein Hodaei, Absar U Hassan, Steffen Wittek, Hipolito Garcia-Gracia, Ramy El-Ganainy, Demetrios N Christodoulides, and Mercedeh Khajavikhan.  Enhanced sensitivity at higher-order exceptional points.  _Nature_ , 548(7666):187–191, 2017.
  * Jeans [1902] James Hopwood Jeans.  I. the stability of a spherical nebula.  _Philosophical Transactions of the Royal Society of London. Series A, Containing Papers of a Mathematical or Physical Character_ , 199(312-320):1–53, 1902.
  * Kirillov [2013] Oleg N Kirillov.  Nonconservative stability problems of modern physics.  In _Nonconservative Stability Problems of Modern Physics_. de Gruyter, 2013.
  * Klaiman et al. [2008] Shachar Klaiman, Uwe Günther, and Nimrod Moiseyev.  Visualization of Branch Points in $\mathcal{}P{}\mathcal{}T{}$-Symmetric Waveguides.  _Physical Review Letters_ , 101(8):080402, August 2008.  .
  * Kreĭn [1950] MG Kreĭn.  A generalization of some investigations on linear differential equations with periodic coefficients.  In _Dokl. Akad. Nauk SSSR A_ , volume 73, page 445, 1950.
  * Laibe and Price [2011] Guillaume Laibe and Daniel J. Price.  dustybox and dustywave: two test problems for numerical simulations of two-fluid astrophysical dustâ€“gas mixtures.  _Monthly Notices of the Royal Astronomical Society_ , 418(3):1491–1497, 12 2011.  ISSN 0035-8711.  .  URL https://doi.org/10.1111/j.1365-2966.2011.19291.x.
  * Laibe and Price [2014] Guillaume Laibe and Daniel J Price.  Dusty gas with one fluid.  _Monthly Notices of the Royal Astronomical Society_ , 440(3):2136–2146, 2014.
  * Leclerc et al. [2023] Armand Leclerc, Lucien Jezequel, Nicolas Perez, Asmita Bhandare, Guillaume Laibe, and Pierre Delplace.  The exceptional ring of buoyancy instability in stars, 2023.
  * Lin et al. [2011] Zin Lin, Hamidreza Ramezani, Toni Eichelkraut, Tsampikos Kottos, Hui Cao, and Demetrios N. Christodoulides.  Unidirectional Invisibility Induced by $\mathcal{}P{}\mathcal{}T{}$-Symmetric Periodic Structures.  _Physical Review Letters_ , 106(21):213901, May 2011.  .
  * Longarini et al. [2022] Cristiano Longarini, Giuseppe Lodato, Giuseppe Bertin, and Philip J Armitage.  The role of the drag force in the gravitational stability of dusty planet forming disc – i. analytical theory.  _Monthly Notices of the Royal Astronomical Society_ , 519(2):2017–2029, dec 2022.  .  URL https://doi.org/10.1093%2Fmnras%2Fstac3653.
  * Magnan et al. [2024] Nathan Magnan, Tobias Heinemann, and Henrik N Latter.  A physical picture for the acoustic resonant drag instability.  _Monthly Notices of the Royal Astronomical Society_ , page stae052, 2024.
  * Makris et al. [2008] K. G. Makris, R. El-Ganainy, D. N. Christodoulides, and Z. H. Musslimani.  Beam Dynamics in $\mathcal{}P{}\mathcal{}T{}$ Symmetric Optical Lattices.  _Physical Review Letters_ , 100(10):103904, March 2008.  .
  * Mostafazadeh [2002a] Ali Mostafazadeh.  Pseudo-hermiticity versus pt symmetry: the necessary condition for the reality of the spectrum of a non-hermitian hamiltonian.  _Journal of Mathematical Physics_ , 43(1):205–214, 2002a.
  * Mostafazadeh [2002b] Ali Mostafazadeh.  Pseudo-hermiticity versus pt-symmetry iii: Equivalence of pseudo-hermiticity and the presence of antilinear symmetries.  _Journal of Mathematical Physics_ , 43(8):3944–3951, 2002b.
  * Mostafazadeh and Batal [2004] Ali Mostafazadeh and Ahmet Batal.  Physical aspects of pseudo-Hermitian and PT-symmetric quantum mechanics.  _Journal of Physics A: Mathematical and General_ , 37(48):11645, November 2004.  ISSN 0305-4470.  .
  * Musslimani et al. [2008] Ziad H Musslimani, Konstantinos G Makris, Ramy El-Ganainy, and Demetrios N Christodoulides.  Optical solitons in p t periodic potentials.  _Physical Review Letters_ , 100(3):030402, 2008\.
  * Onuki [2020] Yohei Onuki.  Quasi-local method of wave decomposition in a slowly varying medium.  _Journal of Fluid Mechanics_ , 883:A56, 2020.
  * Özdemir et al. [2019] Şahin Kaya Özdemir, Stefan Rotter, Franco Nori, and L Yang.  Parity–time symmetry and exceptional points in photonics.  _Nature materials_ , 18(8):783–798, 2019.
  * Peng et al. [2014] Bo Peng, Şahin Kaya Özdemir, Fuchuan Lei, Faraz Monifi, Mariagiovanna Gianfreda, Gui Lu Long, Shanhui Fan, Franco Nori, Carl M. Bender, and Lan Yang.  Parity–time-symmetric whispering-gallery microcavities.  _Nature Physics_ , 10(5):394–398, May 2014.  ISSN 1745-2481.  .
  * Perez [2022] Nicolas Perez.  _Topological waves in geophysical and astrophysical fluids_.  PhD thesis, Ecole normale supérieure de lyon-ENS LYON, 2022.
  * Perez et al. [2021] Nicolas Perez, Pierre Delplace, and Antoine Venaille.  Manifestation of the berry curvature in geophysical ray tracing.  _Proceedings of the Royal Society A_ , 477(2248):20200844, 2021.
  * Perez et al. [2022] Nicolas Perez, Pierre Delplace, and Antoine Venaille.  Unidirectional modes induced by nontraditional coriolis force in stratified fluids.  _Physical Review Letters_ , 128(18):184501, 2022\.
  * Pringle and King [2007] James E Pringle and Andrew King.  _Astrophysical flows_.  Cambridge University Press, 2007.
  * Qin et al. [2019] Hong Qin, Ruili Zhang, Alexander S. Glasser, and Jianyuan Xiao.  Kelvin-Helmholtz instability is the result of parity-time symmetry breaking.  _Physics of Plasmas_ , 26(3):032102, March 2019\.  ISSN 1070-664X.  .
  * Qin et al. [2021] Hong Qin, Yichen Fu, Alexander S. Glasser, and Asher Yahalom.  Spontaneous and explicit parity-time-symmetry breaking in drift-wave instabilities.  _Physical Review E_ , 104(1):015215, July 2021\.  .
  * Rayleigh [1916] Lord Rayleigh.  Lix. on convection currents in a horizontal layer of fluid, when the higher temperature is on the under side.  _The London, Edinburgh, and Dublin Philosophical Magazine and Journal of Science_ , 32(192):529–546, 1916.
  * Regensburger et al. [2012] Alois Regensburger, Christoph Bersch, Mohammad-Ali Miri, Georgy Onishchukov, Demetrios N. Christodoulides, and Ulf Peschel.  Parity–time synthetic photonic lattices.  _Nature_ , 488(7410):167–171, August 2012.  ISSN 1476-4687.  .
  * Rüter et al. [2010] Christian E. Rüter, Konstantinos G. Makris, Ramy El-Ganainy, Demetrios N. Christodoulides, Mordechai Segev, and Detlef Kip.  Observation of parity–time symmetry in optics.  _Nature Physics_ , 6(3):192–195, March 2010.  ISSN 1745-2481.  .
  * Saffman [1962] PG Saffman.  On the stability of laminar flow of a dusty gas.  _Journal of fluid mechanics_ , 13(1):120–128, 1962.
  * Schwarzschild [1906] K. Schwarzschild.  On the equilibrium of the Sun’s atmosphere.  _Nachrichten von der Königlichen Gesellschaft der Wissenschaften zu Göttingen. Math.-phys. Klasse_ , 195:41–53, January 1906.
  * Squire and Hopkins [2018] Jonathan Squire and Philip F Hopkins.  Resonant drag instability of grains streaming in fluids.  _The Astrophysical Journal Letters_ , 856(1):L15, 2018.
  * Toomre [1964] Alar Toomre.  On the gravitational stability of a disk of stars.  _Astrophysical Journal, vol. 139, p. 1217-1238 (1964)._ , 139:1217–1238, 1964.
  * Wolfram Research Inc. [2021] Wolfram Research Inc.  Mathematica, Version 13.0.0, 2021.  URL https://www.wolfram.com/mathematica.  Champaign, IL.
  * Zhang et al. [2020a] Fangxing Zhang, Yaming Feng, Xianfeng Chen, Li Ge, and Wenjie Wan.  Synthetic anti-pt symmetry in a single microcavity.  _Physical review letters_ , 124(5):053901, 2020a.
  * Zhang et al. [2020b] Huilai Zhang, Ran Huang, Sheng-Dian Zhang, Ying Li, Cheng-Wei Qiu, Franco Nori, and Hui Jing.  Breaking anti-pt symmetry by spinning a resonator.  _Nano Letters_ , 20(10):7594–7599, 2020b.
  * Zhang et al. [2020c] Ruili Zhang, Hong Qin, and Jianyuan Xiao.  Pt-symmetry entails pseudo-hermiticity regardless of diagonalizability.  _Journal of Mathematical Physics_ , 61(1):012101, 2020c.
  * Zhuravlev [2019] VV Zhuravlev.  On the nature of the resonant drag instability of dust streaming in protoplanetary disc.  _Monthly Notices of the Royal Astronomical Society_ , 489(3):3850–3869, 2019.

##  Appendix A Rayleigh-Bénard: linear perturbations

Consider a two-dimensional fluid parametrised by the coordinates
$\left(x,z\right)$. The fluid is at rest in a constant gravity field
$\mathbf{g}=-g\mathbf{u}_{z}$. It has constant density $\rho_{0}$ and the
temperature profile has a constant gradient $T(z)=T_{\mathrm{h}}-\beta z$ with
$\beta=(T_{\mathrm{h}}-T_{\mathrm{c}})/L$. The coefficient of thermal
expansion of the fluid is $\alpha=-\rho_{0}^{-1}\left(\partial\rho/\partial
T\right)_{P}$. The dynamical viscosity and thermal diffusivity of the fluid
are denoted by $\nu$ and $\kappa$ respectively. The linear system describing
the evolution of small perturbations in the Boussinesq approximation around
this background state are

$\displaystyle\partial_{x}v_{x}+\partial_{z}v_{z}$ $\displaystyle=$
$\displaystyle 0,$ (64) $\displaystyle(\partial_{t}-\nu\Delta)v_{x}$
$\displaystyle=$ $\displaystyle-\frac{1}{\rho_{0}}\partial_{x}p$ (65)
$\displaystyle(\partial_{t}-\nu\Delta)v_{z}$ $\displaystyle=$
$\displaystyle-\frac{1}{\rho_{0}}\partial_{z}p+\alpha gT,$ (66)
$\displaystyle(\partial_{t}-\kappa\Delta)T$ $\displaystyle=$
$\displaystyle\beta v_{z}.$ (67)

Applying double curl, and enforcing the incompressibility equation
$\nabla\cdot\mathbf{v}=0$, the equation for $u_{x}$ decouples from the other
two equations, leaving

$\displaystyle(\partial_{t}-\nu\Delta)\Delta v_{z}$ $\displaystyle=$
$\displaystyle\alpha g\partial_{xx}T,$ (68)
$\displaystyle(\partial_{t}-\kappa\Delta)T$ $\displaystyle=$
$\displaystyle\beta v_{z}.$ (69)

For plane wave solutions of the form $\exp{(-i\omega t+ik_{z}z+ik_{x}x)}$,
Eqs. 68 – 69 become

$\displaystyle\left(-i\omega-\nu(k_{x}^{2}+k_{z}^{2})\right)(k_{x}^{2}+k_{z}^{2})v_{z}$
$\displaystyle=$ $\displaystyle-\alpha gk_{x}^{2}T,$ (70)
$\displaystyle\left(-i\omega-\kappa(k_{x}^{2}+k_{z}^{2})\right)T$
$\displaystyle=$ $\displaystyle\beta v_{z}.$ (71)

Rescaling the time by $L^{2}/\nu$, the lengths by $L$ and the temperature with
$\beta L$ leads to the introduction of the Rayleigh and Prandtl numbers,
$\mathrm{Ra}=g\alpha\beta L^{4}/\nu\kappa$ and $\mathrm{Pr}=\nu/\kappa$. The
following rescaling of the equations
$T\mapsto\frac{k_{x}}{k}\left(\frac{\mathrm{Ra}}{\mathrm{Pr}}\right)^{1/2}T$
with $k=\sqrt{k_{x}^{2}+k_{z}^{2}}$ provides a system of eigenvalue that has
the symetric form

$\displaystyle\omega X$ $\displaystyle=$ $\displaystyle H_{\rm RB}X,$ (72)
$\displaystyle H_{\rm RB}$ $\displaystyle=$
$\displaystyle-i\begin{pmatrix}k^{2}&-\frac{k_{x}}{k}\left(\frac{\mathrm{Ra}}{\mathrm{Pr}}\right)^{1/2}\\\
-\frac{k_{x}}{k}\left(\frac{\mathrm{Ra}}{\mathrm{Pr}}\right)^{1/2}&k^{2}/\mathrm{Pr}\end{pmatrix},$
(73)

which corresponds to Eq.39. Consequently, $iH_{\rm RB}$ is Hermitian if and
only if $\mathrm{Ra}\geq 0$.

##  Appendix B dustywave: linear perturbations

Consider a dusty mixture where gas and dust are treated as two fluid
interacting with each other via a drag force. For pressureless dust in 1D, the
equations for mass and momentum conservation are

$\displaystyle(\partial_{t}+{v}_{\mathrm{g}}\partial_{x}){\rho}_{\mathrm{g}}$
$\displaystyle=$
$\displaystyle-{\rho}_{\mathrm{g}}\partial_{x}{v}_{\mathrm{g}},$ (74)
$\displaystyle(\partial_{t}+{v}_{\mathrm{d}}\partial_{x}){\rho}_{\mathrm{d}}$
$\displaystyle=$
$\displaystyle-{\rho}_{\mathrm{d}}\partial_{x}{v}_{\mathrm{d}},$ (75)
$\displaystyle(\partial_{t}+{v}_{\mathrm{g}}\partial_{x}){v}_{\mathrm{g}}$
$\displaystyle=$
$\displaystyle\frac{K}{{\rho}_{\mathrm{g}}}({v}_{\mathrm{d}}-{v}_{\mathrm{g}})-\frac{1}{{\rho}_{\mathrm{g}}}\partial_{x}{P},$
(76)
$\displaystyle(\partial_{t}+{v}_{\mathrm{d}}\partial_{x}){v}_{\mathrm{d}}$
$\displaystyle=$
$\displaystyle\frac{K}{{\rho}_{\mathrm{d}}}({v}_{\mathrm{g}}-{v}_{\mathrm{d}}),$
(77)

where $\rho_{\rm g}$, $\rho_{\rm d}$ denote the gas and the dust densities,
$v_{\rm g}$, $v_{\rm d}$ denote the gas and the dust velocities, and $K$ is
the drag coefficient. The mixture is initially homogeneous and at rest.
Writing each quantity as a background and perturbation component
$f=f_{0}+f^{\prime}$. The system of equations can be put in a symmetric form
by using the rescaled quantities
$(\rho^{\prime}_{\mathrm{g}},\rho^{\prime}_{\mathrm{d}},v^{\prime}_{\mathrm{g}},v^{\prime}_{\mathrm{d}})\mapsto
X=\begin{pmatrix}\frac{1}{\rho_{\mathrm{g},0}+\rho_{\mathrm{d},0}}\rho^{\prime}_{\mathrm{g}},&\frac{1}{\rho_{\mathrm{g},0}+\rho_{\mathrm{d},0}}\rho^{\prime}_{\mathrm{d}},&\frac{\rho_{\mathrm{g},0}}{c_{\mathrm{s}}(\rho_{\mathrm{g},0}+\rho_{\mathrm{d},0})}v^{\prime}_{\mathrm{g}},&\frac{\rho_{\mathrm{d},0}}{c_{\mathrm{s}}(\rho_{\mathrm{g},0}+\rho_{\mathrm{d},0})}v^{\prime}_{\mathrm{d}}\end{pmatrix}^{\top}$.
For isothermal pertubations
${P^{\prime}(x,t)=c_{\mathrm{s}}^{2}\rho^{\prime}_{\mathrm{g}}(x,t)}$, the
evolution of the perturbations are given by the system

$\partial_{t}X=\begin{pmatrix}0&0&-c_{\mathrm{s}}\partial_{x}&0\\\
0&0&0&-c_{\mathrm{s}}\partial_{x}\\\
-c_{\mathrm{s}}\partial_{x}&0&-\frac{K}{\rho_{\mathrm{g},0}}&+\frac{(1-\epsilon)K}{\epsilon\rho_{\mathrm{g},0}}\\\
0&0&\frac{\epsilon
K}{(1-\epsilon)\rho_{\mathrm{d},0}}&-\frac{K}{\rho_{\mathrm{d},0}}\end{pmatrix}X,$
(78)

where
$\epsilon\equiv\frac{\rho_{\mathrm{d},0}}{\rho_{\mathrm{g},0}+\rho_{\mathrm{d},0}}$
is the dust fraction, and
$1-\epsilon=\frac{\rho_{\mathrm{g},0}}{\rho_{\mathrm{g},0}+\rho_{\mathrm{d},0}}$
is the gas fraction. The barycentric stopping time
$t_{\mathrm{s}}\equiv\frac{\rho_{\mathrm{g},0}\rho_{\mathrm{d},0}}{K(\rho_{\mathrm{g},0}+\rho_{\mathrm{d},0})}$,
corresponds to the time needed for drag to damp the differential velocity
between gas and dust. The stopping length $c_{\rm s}t_{\rm s}$ is the typical
spatial separation between the fluids. Using $t_{\mathrm{s}}$ and
$c_{\mathrm{s}}t_{\mathrm{s}}$ as units of time and length, we obtain the
following dimensionless system of equations

$\displaystyle H_{\rm dw}X$ $\displaystyle=$ $\displaystyle i\partial_{t}X,$
(79) $\displaystyle H_{\rm dw}$ $\displaystyle\equiv$
$\displaystyle\begin{pmatrix}0&\;\;\;0&\;\;\;\;k&0\\\ 0&\;\;\;0&\;\;\;\;0&k\\\
k&\;\;\;0&\;\;\;\;i\epsilon&\;-i(1-\epsilon)\\\
0&\;\;\;0&\;-i\epsilon&\;\;\;i(1-\epsilon)\end{pmatrix}.$ (80)

##  Appendix C  EPs from resultants

Consider two polynomials $Q1\left(X\right)$ and $Q2\left(X\right)$ defined by

$\displaystyle Q_{1}\left(X\right)$ $\displaystyle=$
$\displaystyle{\displaystyle\sum_{i=1}^{n}}a_{i}X^{i}=a_{n}{\displaystyle\prod_{i=1}^{n}}(X-\alpha_{i}),$
(81) $\displaystyle Q_{2}\left(X\right)$ $\displaystyle=$
$\displaystyle{\displaystyle\sum_{j=1}^{m}}b_{j}X^{j}=b_{m}{\displaystyle\prod_{j=1}^{m}}(X-\beta_{j}),$
(82)

where $\alpha_{i}$ and $\beta_{j}$ are their roots. The resultant
$\mathcal{R}(Q_{1},Q_{2})$ of these two polynomials is defined as

$\mathcal{R}(Q_{1},Q_{2})\equiv\det S(Q_{1},Q_{2}),$ (83)

where $S(Q_{1},Q_{2})$ is the Sylvester matrix of $Q_{1}$ and $Q_{2}$. It
satisfies

$\mathcal{R}(Q_{1},Q_{2})=a_{n}^{m}b_{m}^{n}{\displaystyle\prod_{i,j}}(\alpha_{i}-\beta_{j}),$
(84)

which cancels if and only if $Q_{1}$ and $Q_{2}$ share a common root. The
quantity (83) is determined by the coefficients $a_{i}$, $b_{j}$ and does not
require to find the roots $\alpha_{i}$, $\beta_{j}$. Define $P^{\prime}$, the
derivative of $P(X)$. The cancellation of the resultant $R_{1}$ of $P$ and
$P^{\prime}$ characterize the fact that $P$ has a root with double
multiplicity. In the context of this study, Exceptional Points of order at
least 2 is the space formed by parameters for which

$R_{1}\equiv\mathcal{R}(P,P^{\prime})=0.$ (85)

##  Appendix D Dustywave: topological charges

Figure 6: Exceptional Points of the dustywave in the parameter space
($k,\epsilon$). Intersections between the red ($R_{1}=0$) and blue ($R_{2}=0$)
curves indicate positions of EP3s. Two regions are zoomed in on the two bottom
panels. The EP3s are marked by black points, with associated topological
charge of $\pm 1$.

In addition to EPs of order 2, the Dustywave is also found to have EPs of
order 3. The structure of these EP3s is folds of dimension 0, which are points
in the ($k,\epsilon$) space. These points are located by a procedure which
generalises the determination of EP2s. Noting
$R_{1}=\mathcal{R}(P,P^{\prime})$ and $R_{2}=\mathcal{R}(P,P^{\prime\prime})$
the resultants of the characteristic polynomial and its derivatives, EP3s
corresponds to points in the parameters space for which $P$,$P^{\prime}$ and
$P^{\prime\prime}$ have a common root. These points are those where both
$R_{1}$ and $R_{2}$ cancel simultaneously. For the Dustywave problem, $R_{1}$
is given by Eq. (45), and $R_{2}$ by

$R_{2}=400k^{8}-432k^{4}(\epsilon-1)+96k^{6}(14-51\epsilon+36\epsilon^{2}).$
(86)

Fig.6 shows the curves $R_{1}=0$ (red) and $R_{2}=0$ (blue). These curves
cross four times, and two of these points are EP3s. They are located at

$k=\pm\frac{1}{\sqrt{3}},\quad\epsilon=\frac{8}{9}.$ (87)

For the other two points, $P,P^{\prime\prime}$ share a root different from the
one shared by $P,P^{\prime}$. As such, they are EP2s only. At the EP3s, $H$
has two eigenvalues, 0 and $\frac{-i}{3}$. It is non-diagonalizable, since it
only has two eigenvectors,

$\displaystyle\begin{pmatrix}0\\\ 1\\\ 0\\\
0\end{pmatrix},\begin{pmatrix}-i\sqrt{3}/2\\\ 2i\sqrt{3}\\\ -1/2\\\
2\end{pmatrix}.$ (88)

The EP3s have a non-zero topological charge, as defined in [Delplace et al.,
2021]. These topological invariants are defined as the winding numbers

$W_{3}=\frac{1}{2\pi}\oint_{C}\frac{1}{R_{1}^{2}+R_{2}^{2}}(R_{1}\nabla_{\boldsymbol{p}}R_{2}-R_{2}\nabla_{\boldsymbol{p}}R_{1})\mathrm{d}\boldsymbol{p},$
(89)

for $\nabla_{\boldsymbol{p}}=\begin{pmatrix}\partial_{k}\\\
\partial_{\epsilon}\end{pmatrix}$ and $\boldsymbol{p}=\begin{pmatrix}k\\\
\epsilon\end{pmatrix}\in C$ where $C$ is any closed loop around one of the
EP3s. We compute

$\displaystyle W_{3}\Big{|}_{\mathrm{EP3}_{1}}$ $\displaystyle=$
$\displaystyle+1,$ (90) $\displaystyle W_{3}\Big{|}_{\mathrm{EP3}_{2}}$
$\displaystyle=$ $\displaystyle-1,$ (91)

where $\mathrm{EP3}_{1}$ is the point at
$(k,\epsilon)=(+\frac{1}{\sqrt{3}},\frac{8}{9})$ and $\mathrm{EP3}_{2}$ is the
point at $(k,\epsilon)=(-\frac{1}{\sqrt{3}},\frac{8}{9})$.

##  Appendix E Dustywave: Anti-$\mathcal{PT}$ operator

Define $\Theta$ the operator of complex conjugation and
$\mathcal{CP}\equiv\tilde{U}\Theta$ the operator for anti-$\mathcal{PT}$
symmetry. Equation (9) can be written under the form

$H_{\rm dw}\mathcal{CP}+\mathcal{CP}H_{\rm dw}=0.$ (92)

Eq. (92) constrains the phase between the components of the eigenvectors of
$H_{\rm dw}$, also called relations of polarisation. Consider $X_{0}$ an
eigenvector of $H_{\rm dw}$ with pure imaginary eigenvalue. Formally,

$H_{\rm dw}X_{0}=i\eta_{0}X_{0},$ (93)

with $\eta_{0}$ real. Applying the operator $\mathcal{CP}$ on Eq. (93) and the
commutation rule Eq.(92) gives the relation

$H_{\rm dw}(\mathcal{CP}X_{0})=i\eta_{0}(\mathcal{CP}X_{0}).$ (94)

If the eigenvalue is not degenerated (i.e. the parameters are chosen for the
system not to be at an Exceptional Point), $\mathcal{CP}X_{0}$ is colinear to
$X_{0}$, i.e.

$\mathcal{CP}X_{0}=\alpha X_{0}$ (95)

for some complex number $\alpha$. Hence, $H_{\rm dw}$ and $\mathcal{CP}$ share
the same eigenvectors and the anti-$\mathcal{PT}$ symmetry is said to be
unbroken. As a remark, the same reasoning on the operator $iH_{\rm dw}$ would
lead to show that $iH_{\rm dw}$ is $\mathcal{PT}$-symmetric. $\mathcal{CP}$
can be straightforwardly diagonalised. Its eigenvectors are

$\begin{pmatrix}1\\\ 0\\\ 0\\\ 0\end{pmatrix},\begin{pmatrix}0\\\ 1\\\ 0\\\
0\end{pmatrix},\begin{pmatrix}0\\\ 0\\\ i\\\
0\end{pmatrix},\begin{pmatrix}0\\\ 0\\\ 0\\\ i\end{pmatrix}.$ (96)

The condition given by Eq. (95) imposes therefore the particular polarisation
relations on non-propagating waves:

$\displaystyle\left(\frac{\rho^{\prime}_{\mathrm{g}}}{\rho^{\prime}_{\mathrm{d}}}\right)^{*}$
$\displaystyle=$
$\displaystyle+\frac{\rho^{\prime}_{\mathrm{g}}}{\rho^{\prime}_{\mathrm{d}}},$
(97)
$\displaystyle\left(\frac{\rho^{\prime}_{\mathrm{g}}}{v^{\prime}_{\mathrm{g}}}\right)^{*}$
$\displaystyle=$
$\displaystyle-\frac{\rho^{\prime}_{\mathrm{g}}}{v^{\prime}_{\mathrm{g}}},$
(98)
$\displaystyle\left(\frac{\rho^{\prime}_{\mathrm{d}}}{v^{\prime}_{\mathrm{d}}}\right)^{*}$
$\displaystyle=$
$\displaystyle-\frac{\rho^{\prime}_{\mathrm{d}}}{v^{\prime}_{\mathrm{d}}},$
(99)

and other linear combinations of Eqs. 97 – 99. Hence, the analysis of the
discrete symmetries of the problem shows with a very few amount of
calculations that gas and dust densities are in phase as well as their
velocities, while densities and velocities are in phase quadrature. Similar
relations exist between fields in the unbroken $\mathcal{PT}$ symmetry phase
of the Kelvin-Helmholtz problem discussed by Fu and Qin [2020a].

This result can be reformulated with real fields. In the unbroken
anti-$\mathcal{PT}$ symmetry phase, the perturbations are damped without
oscillating, and are of the form

$\displaystyle\rho^{\prime}_{\mathrm{g}}$ $\displaystyle\propto$
$\displaystyle\mathrm{e}^{-\eta_{0}t}\cos{(kx)},$ (100)
$\displaystyle\rho^{\prime}_{\mathrm{d}}$ $\displaystyle\propto$
$\displaystyle\mathrm{e}^{-\eta_{0}t}\cos{(kx)},$ (101) $\displaystyle
v^{\prime}_{\mathrm{g}}$ $\displaystyle\propto$
$\displaystyle\mathrm{e}^{-\eta_{0}t}\sin{(kx)},$ (102) $\displaystyle
v^{\prime}_{\mathrm{d}}$ $\displaystyle\propto$
$\displaystyle\mathrm{e}^{-\eta_{0}t}\sin{(kx)},$ (103)

which imposes for quantities to be zero when they are spatially averaged.
Formally, with

$\langle
f\rangle(t)\equiv\lim_{L\to\infty}\frac{1}{L}\int_{-L/2}^{L/2}f(x,t)\mathrm{d}x,$
(104)

one obtains for example

$\displaystyle\langle
v^{\prime}_{\mathrm{g}}\rho^{\prime}_{\mathrm{g}}\rangle$ $\displaystyle=$
$\displaystyle 0,$ (105) $\displaystyle\langle
v^{\prime}_{\mathrm{d}}\rho^{\prime}_{\mathrm{d}}\rangle$ $\displaystyle=$
$\displaystyle 0,$ (106)

for all times. On the contrary, in the phase where anti-$\mathcal{PT}$
symmetry is broken, the frequency of the perturbations have a real part. Waves
are propagating, and the polarisation relations are different than
Eqs.(97)-(99), and time-average of quantities defined in Eqs. 105 – 106 are
non zero for this phase. Only propagating waves have non-zero macroscopic mass
flux.