# Tailed skyrmions – an obscure branch of magnetic solitons

Vladyslav M. Kuchkin<EMAIL_ADDRESS>Science Institute, University of Iceland,
107 Reykjavík, Iceland  Peter Grünberg Institute and Institute for Advanced
Simulation, Forschungszentrum Jülich and JARA, 52425 Jülich, Germany
Department of Physics, RWTH Aachen University, 52056 Aachen, Germany  Nikolai
S. Kiselev  Peter Grünberg Institute and Institute for Advanced Simulation,
Forschungszentrum Jülich and JARA, 52425 Jülich, Germany  Pavel F. Bessarab
Science Institute, University of Iceland, 107 Reykjavík, Iceland  Department
of Physics and Electrical Engineering, Linnaeus University, SE-39231 Kalmar,
Sweden

###### Abstract

We report tailed skyrmions – a new class of stable soliton solutions of the 2D
chiral magnet model. Tailed skyrmions have elongated shapes and emerge in a
narrow range of fields near the transition between the spin spirals and the
saturated state. We analyze the stability range of these solutions in terms of
external magnetic field and magnetocrystalline anisotropy. Minimum energy
paths and the homotopies (continuous transitions) between tailed skyrmions of
the same topological charge have been calculated using the geodesic nudged
elastic bands method. The discovery of tailed skyrmions extends the diversity
of already-known solutions illustrated by complex morphology solitons, such as
tailed skyrmion bags with and without chiral kinks.

##  I Introduction

Chiral magnets represent a unique class of materials where the competition
between local interaction, in particular, Heisenberg exchange and
Dzyaloshinskii-Moriya interaction [1, 2] (DMI), give rise to a vast variety of
topological magnetic solitons such as $k\pi$-skyrmions [3, 4, 5, 6], skyrmion
bags [7, 8], and skyrmions with chiral kinks [10, 9]. Note that the physical
systems that exhibit such diversity of solitons possessing arbitrary
topological charges and morphology are very rare in nature.

The primary parameter that distinguishes magnetic solitons is the topological
charge, which defines the homotopy class of a particular solution:

$Q=\dfrac{1}{4\pi}\int\\!\mathbf{n}\cdot\left(\partial_{\mathrm{x}}\mathbf{n}\times\partial_{\mathrm{y}}\mathbf{n}\right)\,\mathrm{d}x\mathrm{d}y,$
(1)

where $\mathbf{n}$ is the magnetization unit vector field.

The solutions with identical $Q$ can be continuously (without the appearance
of magnetic singularities) transformed into each other [7]. We refer to such
transformations as homotopies, which can be defined as follows. For two
magnetic textures $\mathbf{n}_{1}(\mathbf{r})$ and
$\mathbf{n}_{2}(\mathbf{r})$ with identical $Q$, there is a unity parametric
vector field $\mathbf{n}(\mathbf{r};\tau)$, $\tau\in[0,1]$, which is
continuously differentiable with respect of $\mathbf{r}$ and $\tau$ and obeys
the boundary conditions $\mathbf{n}(\mathbf{r};0)=\mathbf{n}_{1}(\mathbf{r})$,
$\mathbf{n}(\mathbf{r};1)=\mathbf{n}_{1}(\mathbf{r})$. Homotopy between any
two configurations $\mathbf{n}_{1}(\mathbf{r})$ and
$\mathbf{n}_{2}(\mathbf{r})$ is not defined uniquely as there are infinitely
many ways to introduce the parameter $\tau$. A constraint can be introduced by
considering only those homotopies that represent minimum energy paths (MEPs),
i.e. the paths lying lowermost on the energy surface of the system. Homotopies
satisfying the MEP condition are also physically relevant as they define
thermal stability of the magnetic configurations within harmonic rate
theories. For solitons with different $Q$, homotopies are, by definition,
impossible.

The representative and the most well-studied soliton in chiral magnets –
$\pi$-skyrmion with $Q=-1$ – is shown in Figure 1A-B. $\pi$-skyrmions were
theoretically predicted more than thirty years ago by Bogdanov and Yablonskii
[3, 4, 5]. Nowadays, there is a long list of magnetic materials where
$\pi$-skyrmions were observed experimentally [11].

Another type of axially symmetric solitons called $k\pi$-skyrmions ($k>1$) was
presented by Bogdanov and Hubert [6]. Despite the morphological diversity of
$k\pi$-skyrmions, they all belong to only two homotopy classes, with $Q=0$ and
$Q=-1$ for even and odd $k$, respectively.

The non-axially symmetric solutions with arbitrary topological charge, also
known as skyrmion bags [7, 8] and skyrmions with chiral kinks [10, 9] were
found only recently. The skyrmion bags with positive topological charge have
recently been observed in FeGe plates using transmission electron microscopy
[12]. Co-existing skyrmion and its antiparticle (antiskyrmion) i.e. soliton
with two chiral kinks were observed experimentally in FeGe thin platelet [13].

In this article, we introduce a new class of chiral magnetic skyrmions we
refer to as tailed skyrmions. We demonstrate rich diversity of tailed
skyrmions and present a comprehensive study of their static properties,
stability and finite-temperature dynamics. Tailed skyrmions can be thought of
as various elongated skyrmions stabilized in the narrow range of external
magnetic fields close to the elliptical instability field of an ordinary
$\pi$-skyrmion. Some examples of tailed skyrmions are shown in Fig. 1. To our
knowledge, this type of solitons has not been reported in the literature.

##  II Model

We consider the two-dimensional (2D) micromagnetic model of a chiral magnet
containing three energy terms:

$\mathcal{E}=\int\\{w_{\mathrm{ex}}(\textbf{n})+w_{\mathrm{D}}(\textbf{n})+w_{\mathrm{U}}(\textbf{n})\\}\,l\mathrm{d}x\mathrm{d}y,$
(2)

where $\textbf{n}=\textbf{M}/M_{\mathrm{s}}$ is the magnetization unit vector
field which is uniform across the film thickness $l$, $M_{\mathrm{s}}$ is the
saturation magnetization,
$w_{\mathrm{ex}}(\textbf{n})=\mathcal{A}\left|\nabla\mathbf{n}\right|^{2}$ is
the Heisenberg exchange interaction and
$w_{\mathrm{U}}(\textbf{n})=M_{\mathrm{s}}B_{\mathrm{ext}}\left(1-n_{\mathrm{z}}\right)+\mathcal{K}\left(1-n_{\mathrm{z}}^{2}\right)$
is the potential term containing the Zeeman interaction and the easy-
axis/easy-plane anisotropy. The external magnetic field is perpendicular to
the plane of the film,
$\textbf{B}_{\mathrm{ext}}\parallel\textbf{e}_{\mathrm{z}}$.

Figure 1:  Tailed skyrmions. (A) The magnetization vector field of the
$\pi$-skyrmion. (B) Top view of a $\pi$-skyrmion, where colors encode the
magnetization direction. (C)-(E) Color-coded spin texture of tailed skyrmions
with $Q=-1$. (F) Stability diagram of tailed skyrmions in terms of the
magnetic field, $h$, and magnetocrystalline anisotropy, $u$. (G) Minimum
energy path between skyrmions depicted in (B)-(E). All simulations were
performed on the $4L_{\mathrm{D}}\times 4L_{\mathrm{D}}$ domain.

The DMI term $w_{\mathrm{D}}(\mathbf{n})=\mathcal{D}w(\mathbf{n})$ is defined
by combinations of Lifshitz invariants [3],
$\Lambda_{ij}^{(k)}\\!=\\!n_{i}\partial_{k}n_{j}\\!-\\!n_{j}\partial_{k}n_{i}$.
The results presented below are valid for systems with Bloch-type modulations
[14], Néel-type modulations [15, 16, 17] as well as for crystals with D2d or
S4 point group symmetry [3]. Without loss of generality, we assume Bloch-type
DMI,
$w(\mathbf{n})\\!=\\!\Lambda_{\mathrm{zy}}^{(\mathrm{x})}\\!+\\!\Lambda_{\mathrm{xz}}^{(\mathrm{y})}\\!$,
in our calculations. By introducing the characteristic size of chiral
modulations $L_{\mathrm{D}}=4\pi\mathcal{A}/\mathcal{D}$ and the
characteristic magnetic field
$B_{\mathrm{D}}=\mathcal{D}^{2}/(2M_{\mathrm{s}}\mathcal{A})$, we reduce the
number of independent parameters to two, namely the dimensionless external
magnetic field $h=B_{\mathrm{ext}}/B_{\mathrm{D}}$ and anisotropy parameter
$u=\mathcal{K}/(M_{\mathrm{s}}B_{\mathrm{D}})$.

##  III Methods

This work involves several numerical methods. Stable magnetic configurations
are obtained via direct energy minimization using the conjugate gradient
method with stereographic projections of the $\mathbf{n}$-field. For details,
see Supplementary materials in Ref. [18]. MEPs) between the stable states are
calculated using the geodesic nudged elastic band (GNEB) method [19, 20].
Dynamical properties of tailed skyrmions are investigated via the numerical
integration of the Landau-Lifshitz-Gilbert (LLG) equation [21] using the well-
established semi-implicit method from Ref. [22].

The calculations were performed and double-checked with various software. In
particular, we used Excalibur [23] and Mumax3 [24] for energy minimization and
calculation of stability ranges. MEPs were calculated using the Spirit code
[25]. We used high-accuracy finite difference schemes in the energy
minimization to remove numerical artifacts in the calculation of derivatives
in Eq. (2). A typical mesh density in our simulations is 64 nodes per
$L_{\mathrm{D}}$. The details can be found in Refs. [7, 9].

For the MEP calculations and stochastic LLG simulations, we used an effective
spin-lattice model of nearest neighbors. All calculations were performed with
periodic boundary conditions in the film plane. For visualization of the
magnetization, we use the standard color code as explained in Figure 1A and B.

##  IV Initial guesses

To obtain statically stable solutions for tailed skyrmions in Mumax3, we
constructed initial guesses for the $\mathbf{n}$-field by placing domains with
magnetization antiparallel to the field, $\mathbf{n}=(0,0,-1)$, into the state
with magnetization pointing along the field, $\mathbf{n}=(0,0,1)$. In the case
of axially symmetric skyrmions, e.g., $\pi$-skyrmion in Figure 1A and B, we
use a circular shape domain with magnetization down. For tailed skyrmions
shown in Figure 1C-E, the domains with down magnetization should approximately
mimic the shape of the skyrmion, which can be achieved by placing circular
domains in a chain. For details, see the Supplementary materials with the
corresponding Mumax script.

In the simulation with Excalibur and Spirit, we used the built-in options
allowing one to drag and invert the magnetization field interactively, as
illustrated in Supplementary Movie 1. Notably, the interactive modification of
the magnetization field with the drag tool permits for examining whether the
skyrmion has reached the optimal shape and size.

##  V Results

###  V.1 Stability diagram

Figure 1C, D, and E show elementary tailed skyrmions. The two solutions with
two and three tails can be referred to as a stick-skyrmion and a spinner-
skyrmion, respectively. We use these representative solutions to illustrate
the range of tailed skyrmion existence. In Figure 1F, the critical fields for
stick-skyrmion and spinner-skyrmion are shown with blue and red curves,
respectively. The field range corresponding to stable tailed skyrmions is
bound from above by the transition field into axially symmetric
$\pi$-skyrmion, $h_{\mathrm{t}}$, and from below by elliptic instability
field, $h_{\mathrm{e}}$. Note that Figure 1F shows the deviation of
$h_{\mathrm{t}}$ and $h_{\mathrm{e}}$ from $h^{*}$ – the critical field
corresponding to the transition between the spin-spiral state and the
ferromagnetic phase – so as to emphasize the difference between the critical
fields which are quite close to each other. The dependence $h^{*}(u)$ is the
solution of the well-known equation [5]:

$\displaystyle
2\pi-4\sqrt{h+2u}-\dfrac{\sqrt{2}h}{\sqrt{u}}\ln\dfrac{\sqrt{h+2u}+\sqrt{2u}}{\sqrt{h+2u}-\sqrt{2u}}=0,$
(3)

which is derived from the criterion that the energy of the isolated spiral
($2\pi$-domain wall) equals the energy of the saturated ferromagnetic state.
The functional dependence $h^{*}(u)$ is shown in the inset of Figure 1F.

The critical fields for any tailed skyrmion solutions always meet at two
points: i) the Bogomolnyi point $h=1$, $u=-0.5$ and ii) the phase transition
point between the spin spiral and ferromagnetic phases $h=0$, $u=\pi^{2}/8$.
Thus, the stability range of tailed skyrmions is quite significant and
includes both in-plane and strong out-of-plane anisotropy systems.

The range of the external magnetic field corresponding to stable spinner-
skyrmion (see the domain bound by the red lines in Fig. 1F) is roughly eight
times larger than that for the stick-skyrmion (see the domain bound by the
blue lines). The stability range for most of the other tailed-skyrmion
solutions considered in this work was found to coincide with high precision
with the stability range of the spinner-skyrmion. Therefore, the area bound by
the red curves in Figure 1F can be considered as a good estimate for the
stability range of all tailed-skyrmions.

###  V.2 Homotopies and minimum energy paths

Figure 1G shows the MEP connecting the $\pi$-skyrmion, stick-skyrmion,
spinner-skyrmion and four-tailed skyrmion (see Figure 1B, C, D, and E). The
calculated MEP turns out to be a homotopy as the topological charge of the
system remains the same, $Q=-1$, at every point of the MEP. Some other
examples of homotopy MEPs have already been presented in Ref. [26], reporting
the method for finding exotic three-dimensional hybrid skyrmion tubes.

Starting from the $\pi$-skyrmion solution and following the MEP we observe a
sequential increase in the number of tails in the magnetic texture. The growth
of a new tail requires overcoming an energy barrier that turns out to be
significantly larger than that for the reverse process, i.e., removal of the
tail. The latter defines the stability of the tailed states. With an
increasing external magnetic field, the barriers corresponding to the tail
contraction decrease and eventually disappear, leading to the transition of
the tailed skyrmions into the axially symmetric $\pi$-skyrmion state.

As follows from the MEP calculations, the energy of the tailed skyrmions
increases with the number of tails, meaning that tailed skyrmions represent
higher energy metastable states. The latter is consistent with the previous
study [27] where it was proven that the axially symmetric $\pi$-skyrmion
solution represents the global minimum among all solutions with $Q=-1$.

The question about the second lowest energy soliton with $Q=-1$ remains open.
For instance, $3\pi$-skyrmion and some of the skyrmion bags seem to be
promising candidates. However, comparing the energies of tailed skyrmions, we
found that the stick skyrmion may certainly compete with these solutions. In
particular, at $h=0.617$ and $u=0$, the energy of stick-skyrmion is only
$8.4\%$ higher than that of the $\pi$-skyrmion, while the energy of
$3\pi$-skyrmion for these parameters is $11\%$ higher. Moreover, the
$3\pi$-skyrmion is significantly bigger than the $\pi$-skyrmion and has a
diameter of $\sim 8L_{\mathrm{D}}$. The size of other skyrmion bags with
$Q=-1$ is even larger. This means that the stick-skyrmion is not only lower in
energy than other solutions but also has the shortest distance to the
$\pi$-skyrmion in the configuration space.

###  V.3 Thermally activated transition from $\pi$-skyrmion to stick-skyrmion

Figure 2:  Snapshots of the system during the heating-cooling protocol of
stochastic LLG simulations. (A) The initial configuration of $\pi$-skyrmion.
(B)-(D) Snapshots of the system after $10^{4}$, $5\cdot 10^{6}$ and $10^{7}$
iterations of stochastic LLG simulations with fixed temperature of
$T=0.2J/k_{\mathrm{B}}$. (E), (F) Snapshots after $10^{5}$ and $10^{6}$ LLG
time steps with $T=0$. Simulations are performed on the $8L_{\mathrm{D}}\times
8L_{\mathrm{D}}$ computational domain for $h=0.617$, $u=0$.

Since the tailed stick-skyrmion solution is the closest local minimum to the
$\pi$-skyrmion state, excitation of a $\pi$-skyrmion by external stimuli such
as thermal fluctuations can induce a nucleation of the stick-skyrmion. To
demonstrate this, we performed finite-temperature spin dynamics simulations
for $h=0.617$ and $u=0$. We used the LLG time step $\Delta
t=0.01J\gamma\mu_{\mathrm{s}}^{-1}$, and damping parameter $\alpha=0.1$. The
system was simulated at the temperature of $T=0.2J/k_{\mathrm{B}}$, which is
significantly lower than the Curie temperature, which was estimated to be
$T_{\mathrm{c}}=0.7J/k_{\mathrm{B}}$ for the given system [28].

The initial magnetic texture of axially symmetric $\pi$-skyrmion is shown in
Figure 2A. The thermal fluctuations result in skyrmion deformation, as seen in
Figure 2B-D. After turning off the temperature and performing cooling with the
LLG simulation at $T=0$, the deformed skyrmion quickly turned into a stick-
skyrmion, as shown in Figure 2E-F and in Supplementary Movie 2. A longer
heating of the system, in principle, can also give rise to the appearance of
spinner-skyrmion and other multi-tailed skyrmions. However, such events have a
much lower probability because of the large distance between corresponding
minima in the parameter space and higher energy barriers, Figure 1G.
Nevertheless, our results suggest a simple approach to the experimental
observation of tailed skyrmions.

###  V.4 Tailed skyrmions with arbitrary topological charge

So far, we have only discussed solutions with $Q=-1$. However, the tails can
in principle be added to any soliton with an arbitrary topological charge.
Adding a tail creates a new soliton, but does not change the topological
charge of the system.This section provides examples of such skyrmions. All
solutions presented in the following were obtained on a square simulation
domain of size $L_{\mathrm{x}}=L_{\mathrm{y}}=8L_{\mathrm{D}}$ for parameters
$h=0.617$ and $u=0$.

In Figure 3, we show topologically trivial solitons – skyrmionium with tails.
Note that all solitons shown in Figure 3 are homotopically equivalent despite
different number of tails and their position on the inner or outer contour of
the solitons. Noticeably, the skyrmionium with more than three tails on the
outer contour, see, e.g., Figure 3G-I, has rotational symmetry of order equal
to the number of tails.

Figure 3:  Tailed skyrmions with $Q=0$. Skyrmionium with tails on the inner
contour only (A), (B), on both contours (C), and on the outer contour only (D)
- (I). The scale bar is the same for all images.

Figure 4 provides examples of tailed skyrmion bags of various topological
charges and demonstrates diversity of tailed skyrmions. For example, Figure
4C-E shows three skyrmion configurations with the same topological charge
$Q=-1$ but different symmetry and different number of tails. Figure 4G and H
demonstrates two different skyrmion bags with the same positive topological
charge, $Q=+3$.

Figure 4:  Skyrmion bags with tails stabilized for $u=0$ and $h=0.617$.
Figure 5:  Skyrmions with tails and chiral kinks stabilized for $h=0.617$ and
$u=0$. (B) Detailed view of the spin texture near the chiral kink of skyrmion
shown in (A).  Figure 6:  (A) axially symmetric $\pi$-skyrmion. (B), (C)
stick-skyrmions of different lengths. (D) chiral droplet soliton. (E) and (F)
stick-skyrmion with one tail and one chiral kink. All solutions co-exist for
$u=1.194$ and $h=0.01$. The energy minimization was performed on a large size
domain, $L_{\mathrm{y}}=4L_{\mathrm{D}}$ and
$L_{\mathrm{x}}=16L_{\mathrm{D}}$, to avoid the influence of periodic boundary
conditions. All images have identical scale, as marked in (A). (G)
$y$-component of magnetization along the dashed symmetry line of the skyrmions
(A)-(C). Hollow circles denote the zeros of $n_{\mathrm{y}}(x)$.

Another qualitatively new feature appears in the case of skyrmions with chiral
kinks and tails. Figure 5 provides several examples of such skyrmions having
various topological charges. Note that the chiral kink on the outer contour of
the skyrmion carries the topological charge $Q=+1$. Because of this, the
topological charge of skyrmions in Figure 5 A and B equals zero and does not
depend on the number of tails. The skyrmions with antikinks of topological
charge $Q=-1$ are only stable at strong out-of-plane anisotropy [9].

Figures 5 C-E illustrate that the chiral kinks form coupled pairs in the
presence of tails. A similar phenomenon was observed for chiral kinks on the
$2\pi$-domain wall [9]. Generally, the chiral kinks on a straight domain wall
repel each other at large distances and attract each other only at small
distances [9]. Similarly, a pair of kinks can form a coupled state on tailed
skyrmions, with finite distance between the kinks on the order of $0.1L_{D}$.
With increasing field, the solutions with $Q=+1$ shown in Figures 5 C-E
collapse to an antiskyrmion [29]. Figures 5 F-I show skyrmion bags with
positive topological charges hosting chiral kinks and tails simultaneously.

###  V.5 Effect of strong perpendicular anisotropy on tailed skyrmions

So far, the diversity of tailed skyrmions has been demonstrated for the $u=0$
case. However, all tailed skyrmion solutions can also be obtained for any
value of $u$ in the range from $-1/2$ to $\pi^{2}/8$ (see the diagram in Fig.
1F). Moreover, for strong perpendicular anisotropy, $u\lesssim\pi^{2}/8$, we
found distinct co-existing solutions representing tailed skyrmions of
different lengths. In Figure 6B, C, we show several examples of stick-
skyrmions of different lengths, which coexist for a given $u$, $h$, and are
separated by small energy barriers. The energy of the stick-skyrmion increases
with its length. Note that such tailed skyrmions of different lengths are
observed only at strong anisotropy. As Figure 6G shows, the magnetic structure
of such stick-skyrmions is non-trivial. In particular, the $y$ component of
magnetization is modulated along the skyrmion symmetry axis giving an energy
gain sufficient for the skyrmion stabilization. The number of zeros of
$n_{\mathrm{y}}$ increases with the length of a stick-skyrmion.

Interestingly, the stick-skyrmion with one tail and one chiral kink (see
Figure 6D-F) also exhibits multiple minima corresponding to the solutions of
different lengths.

##  VI Conclusions

In conclusion, we presented in this work a new class of chiral skyrmions –
tailed skyrmions. We identified the range of uniaxial anisotropy and
perpendicular magnetic field where the tailed skyrmions remain statically
stable. We show a wide diversity of tailed skyrmions with arbitrary
topological charges. The presence of tails changes only the soliton’s symmetry
but does not change its topological charge. Thus, the transitions between
tailed and tail-free solutions are homotopies. We demonstrated that such
transitions can be identified using the GNEB method. The energy barriers
involved in the homotopies tend to be smaller than what needs to be overcome
to change the topological charge of the system. As a result, tailed skyrmions
can be obtained from a regular $\pi$-skyrmion by applying weak thermal
fluctuations, as demonstrated by the stochastic LLG simulations.

Our results contribute to the development of a complete picture of the
diversity of statically stable magnetic solitons in chiral magnets.

ACKNOWLEDGMENTS

We are grateful to F. Rybakov for fruitful discussions and providing us with
the Excalibur code. The authors acknowledge financial support from the
Icelandic Research Fund (Grant No. 217750), the University of Iceland Research
Fund (Grant No. 15673), the Swedish Re- search Council (Grant No. 2020-05110),
and the European Research Council (ERC) under the European Union’s Horizon
2020 research and innovation program (Grant No. 856538, project “3D MAGiC”).

## References

  * [1] I. Dzyaloshinsky, A thermodynamic theory of “weak” ferromagnetism of antiferromagnetics, J. Phys. Chem. Solids 4, 241 (1958).
  * [2] T. Moriya, Anisotropic superexchange interaction and weak ferromagnetism, Phys. Rev. 120, 91 (1960).
  * [3] A. N. Bogdanov and D. A. Yablonskii, Thermodynamically stable “vortices” in magnetically ordered crystals. The mixed state of magnets, Sov. Phys. JETP 68, 101 (1989).
  * [4] A. N. Bogdanov and A. Hubert, The Properties of Isolated Magnetic Vortices, Phys. Status Solidi B 186, 527 (1994).
  * [5] A. N. Bogdanov and A. Hubert, Thermodynamically stable magnetic vortex states in magnetic crystals, J. Magn. Magn. Mater. 138, 255 (1994).
  * [6] A. N. Bogdanov and A. Hubert, The stability of vortex-like structures in uniaxial ferromagnets, J. Magn. Magn. Mater. 195, 182-192 (1999).
  * [7] F. N. Rybakov and N. S. Kiselev, Chiral magnetic skyrmions with arbitrary topological charge, Phys. Rev. B 99, 064437 (2019).
  * [8] D. Foster, C. Kind, P. J. Ackerman, J.-S. B. Tai, M. R. Dennis, and I. I. Smalyukh, Two-dimensional skyrmion bags in liquid crystals and ferromagnets, Nat. Phys. 15, 655 (2019).
  * [9] V. M. Kuchkin, B. Barton-Singer, F. N. Rybakov, S. Blügel, B. J. Schroers, and N. S. Kiselev, Magnetic skyrmions, chiral kinks and holomorphic functions, Phys. Rev. B 102, 144422 (2020).
  * [10] B. Barton-Singer, C. Ross, C. and B. J. Schroers, Magnetic skyrmions at critical coupling, Commun. Math. Phys. 375, 2259 (2020).
  * [11] Magnetic Skyrmion Materials Y. Tokura and N. Kanazawa Chem. Rev. 121, 2857 (2021).
  * [12] J. Tang, et al. Magnetic skyrmion bundles and their current-driven dynamics. Nat. Nanotechnol. 16, 1086-1091 (2021).
  * [13] F. Zheng, N. S. Kiselev, L. Yang, V. M. Kuchkin, F. N. Rybakov, S. Blügel and R. E. Dunin-Borkowski, Skyrmion–antiskyrmion pair creation and annihilation in a cubic chiral magnet, Nat. Phys. 18, 863–868 (2022).
  * [14] E. Davoli, G. Di Fratta, D. Praetorius, and M. Ruggeri, Micromagnetics of thin films in the presence of Dzyaloshinskii-Moriya interaction, Preprint at arXiv:2010.15541.
  * [15] N. Romming, C. Hanneken, M. Menzel, J. E. Bickel, B. Wolter, K. Von Bergmann, A. Kubetzka, and R. Wiesendanger, Writing and Deleting Single Magnetic Skyrmions, Science 341, 636 (2013).
  * [16] I. Kézsmárki, S. Bordács, P. Milde, E. Neuber, L. M. Eng, J. S. White, H. M. Rønnow, C. D. Dewhurst, M. Mochizuki, K. Yanai, H. Nakamura, D. Ehlers, V. Tsurkan, and A. Loidl, Néel-type skyrmion lattice with confined orientation in the polar magnetic semiconductor GaV4S8, Nat. Mater. 14, 1116 (2015).
  * [17] N. Romming, A. Kubetzka, C. Hanneken, K. von Bergmann, and R. Wiesendanger, Field-Dependent Size and Shape of Single Magnetic Skyrmions, Phys. Rev. Lett. 114, 177203 (2015).
  * [18] F.N. Rybakov, A.B. Borisov, S. Blügel, and N.S. Kiselev, “New Type of Stable Particlelike States in Chiral Magnets,” Phys. Rev. Lett. 115, 117201 (2015); Supplementary material.
  * [19] P.F. Bessarab, V.M. Uzdin, and H. Jónsson, Comput. Phys. Commun. 196, 335 (2015).
  * [20] P.F. Bessarab, Phys. Rev. B 95, 136401 (2017).
  * [21] L. D. Landau and E. M. Lifshitz, On the theory of the dispersion of magnetic permeability in ferromagnetic bodies, Physik. Zeits. Sowjetunion 8, 153 (1935).
  * [22] J. H. Mentink, M. V. Tretyakov, A. Fasolino, M. I. Katsnelson and Th, Rasing, J. Phys.: Condens. Matter 22, 176001 (2010).
  * [23] F. N. Rybakov and E. Babaev, Excalibur software, http://quantumandclassical.com/excalibur/.
  * [24] A. Vansteenkiste, J. Leliaert, M. Dvornik, M. Helsen, F. Garcia-Sanchez, and B. Van Waeyenberge, The design and verification of MuMax3, AIP Advances 4, 10 (2014).
  * [25] G. P. Müller, M. Hoffmann, C. Dißelkamp, D. Schürhoff, S. Mavros, M. Sallermann, N. S. Kiselev, H. Jónsson, and S. Blügel, Spirit: Multifunctional framework for atomistic spin simulations, Phys. Rev. B 99, 224414 (2019).
  * [26] V. M. Kuchkin and N. S. Kiselev, Homotopy transitions and 3D magnetic solitons, APL Materials 10, 071102 (2022).
  * [27] C. Melcher, Chiral skyrmions in the plane, Proc. R. Soc. A 470, 20140394 (2014).
  * [28] V. M. Kuchkin, P. F. Bessarab, and N. S. Kiselev, Thermal generation of droplet soliton in chiral magnet, Phys. Rev. B 105, 184403 (2022)
  * [29] V. M. Kuchkin and N. S. Kiselev, Turning a chiral skyrmion inside out, Phys. Rev. B 101, 064408 (2020).


11institutetext: McGill University, Montréal, Qc, Canada
11email<EMAIL_ADDRESS>22institutetext: Trottier Space Insitute,
Montréal, Qc, Canada

# Phantom Fluid Cosmology
or
Ghosts for Gordon

James M. Cline

###### Abstract

The quanta of phantom dark energy models are negative energy particles, whose
maximum magnitude of energy $|E|$ must be less than a cutoff $\Lambda\lesssim
20$ MeV. They are produced by spontaneous decay of the vacuum into phantoms
plus normal particles. I review general cosmological constraints that have
been derived from the effects of such phantom fluid production, and a possible
application: the generation of boosted dark matter or radiation that could be
directly detected. Recent excess events from the DAMIC experiment can be well-
fit by such processes.

###### keywords:

Dark energy, dark matter

##  1 Happy Birthday, Gordon!

I never had the pleasure of collaborating with Gordon Semenoff, but in 2004 we
co-organized a conference at the Banff International Research Station, called
“New horizons in string cosmology.” The conference photo, Fig. 1, has a
peculiar feature, which I asked the audience of the present meeting whether
they could identify. Namely, the four participants of this meeting who also
attended the BIRS meeting are all grouped together in the bottom left-hand
corner. This seems to be a violation of thermodynamics, or causality.111There
was another peculiar feature which the audience equally failed to notice: I
was wearing exactly the same clothes on both occasions, down to the glasses.
It was a little difficult to see anything while giving the talk.

For the present contribution, I tried to think of something that could be
remotely of interest to Gordon, considering our rather different research
paths, subsequent to our near overlap in the area of string cosmology. What
follows is the best I was able to come up with. It is a model of a new kind of
dark energy, derived from a rather old one, phantom scalar fields. Unlike the
old version, which relies upon the condensate of a scalar field with a wrong-
sign kinetic term, our version assumes the condensate is absent, but the
corresponding negative-energy particles, that are inevitably produced by the
decay of the vacuum, are present [1]. I will start with a review of the
original phantom dark energy models by way of motivation.

Figure 1: Participants of the 2004 BIRS workshop “New horizons in string
cosmology.”

##  2 Phantom dark energy: a brief history

I was at CERN in 2003 on my first sabbatical when WMAP (Wilkinson Microwave
Anisotropy Probe) made its first data release. Amongst other things, these
data led to constraints on the equation of state of the dark energy that
drives the acceleration of the Universe [2], $w=p/\rho$, where $p$ is the
pressure and $\rho$ the energy density, with allowed regions extending to
values of $w<-1$, in violation of the dominant energy condition (DEC).
Theorists who had been analyzing various combinations of cosmological data,
including pre-WMAP data of the cosmic microwave background, had also found
such evidence, for example Refs. [3, 4].

Even before these experimental hints, R. Caldwell had pointed out a simple way
of violating the DEC in the context of scalar field theories [5]. Recall that
a cosmological constant with stress-energy tensor $T^{\mu\nu}=\Lambda\,{\rm
diag(1,-1,-1,-1)}$ saturates the DEC, while a homogeneous scalar field has
equation of state

$w={{\textstyle{1\over 2}}\dot{\phi}^{2}-V\over{\textstyle{1\over
2}}\dot{\phi}^{2}+V}\geq-1\,,$ (1)

where the inequality holds assuming that $V>0$. This interpolates between
cosmological constant, in the limit $\dot{\phi}\to 0$, and $w=1$ when $V\to
0$. Caldwell’s idea was to consider a scalar field with a wrong-sign kinetic
term, a “phantom” field, such that $\dot{\phi}^{2}\to-\dot{\phi}^{2}$ in Eq.
(1) (but $V$ remains $\geq 0$). Then $w<-1$ if $\dot{\phi}^{2}<2V$, and
diverges to $-\infty$ as $\dot{\phi}^{2}\to 2V$. This will lead to a “big rip”
singularity within a finite time, in which the expansion rate of the universe
diverges and gravitationally bound structures are ripped apart.222It is
interesting that Caldwell’s paper, although written in 1999, was not published
until 2002, when CMB data were hinting at the possibility of $w<-1$. This is
readily understood from the classical equation of motion of a homogeneous
phantom field,

$\ddot{\phi}+2H\dot{\phi}=+V^{\prime}(\phi)\,;$ (2)

the field climbs up the potential instead of rolling down.

Currently there is no compelling evidence for violation of the DEC; the
Particle Data Group gives $w=-1.028\pm 0.031$ [6]. However, this has not
stopped theorists from continuing to consider phantom dark energy. InSpire
database gives 1931 results, showing a steady rate of publication over the
last two decades. The vast majority of these papers consider the phantom as a
classical field.

But phantoms are not exempt from quantum mechanics, and one should ask what is
the quantum field theoretic consequence of a wrong sign kinetic term. The
answer depends on the choice of sign for the $i\epsilon$ prescription in the
propagator,

$P(p)={-i\over p^{2}-m^{2}\pm i\epsilon}\,.$ (3)

The positive sign leads to negative probability as in the case of Faddeev-
Popov ghosts, whereas the negative sign corresponds to positive probability,
but with negative energy particles propagating into the future. The latter
choice is the lesser of the two evils, and so one must consider phantom
particles to carry negative energy if they are to be taken seriously as
quantum fields.

Ref. [7] understood that the phase space of phantom particles must be cut off
at some scale $\Lambda$ to avoid divergent decay rates, since arbitrarily
negative energy phantoms could otherwise appear in the final state. But it was
Ref. [8] that pointed out the most stringent irreducible constraint, arising
from decay of the vacuum into two phantoms plus two photons, mediated by
gravity. The Feynman diagram is shown in Fig. 2. On dimensional grounds, it
leads to a constant decay rate per unit volume of order
$\Gamma\sim\Lambda^{8}/m_{P}^{4}$, where $m_{p}$ is the (reduced) Planck mass.
The ensuing photon spectrum roughly takes the form

${dn\over dE}\sim{\Lambda^{7}t_{0}\over m_{p}^{4}}\,\Theta(\Lambda-E)\,,$ (4)

where $t_{0}$ is the age of the Universe. It rises steeply in amplitude with
$\Lambda$, but only linearly in energy. The maximum value of $\Lambda$ such
that the prediction lies below the observed diffuse photon background was
estimated to be $\sim 3\,$MeV, which is rather low for an ultraviolet cutoff.
Above this scale, the phantom description must break down and revert to a more
normal-looking field theory.

Figure 2: Left: Feynman diagram for decay of the vacuum to two ghosts plus two
photons. Right: schematic illustration of predicted x-ray spectrum (for two
values of $\Lambda$, in red) versus cosmological diffuse photon background.

##  3 Revisiting phantom fluids

Since 2003, there have been a few developments motivating a second look at the
cosmological production of ghost particles. One is the Hubble tension: the
discrepancy between the CMB determination of the Hubble constant, and that
coming from more local measurments such as type Ia supernovae; see Ref. [9]
for a review. My naive idea was to suppose that the vacuum could decay into
massless ghosts plus dark matter. Due to the faster redshifting of the ghost
energy density, there would be a net gain in dark matter as a function of
time, which could explain why the Universe expands more rapidly at late times
than expected. But this should not be a viable resolution, as it has been
widely understood that late dark energy solutions do not work for the Hubble
tension, while early dark energy may be able to do so. The current mechanism
is late-time since the bulk of the accumulated energy is appearing today. The
problem with the naive argument is that the appearance of the CMB depends upon
the entire expansion history of the Universe, and not just that at early
times.

Dark matter produced by vacuum decay to ghosts must be light, with a mass
below the $\sim$MeV cutoff. At the same time, it will be produced with an
energy of order $\Lambda$, hence it is boosted to much higher speeds than
ordinary DM, which virializes in the galaxy with speed $\sim 10^{-3}c$. This
is reminiscent of another recent development, the interest in boosted light DM
[10]. Such light DM particles could not be directly detected if moving slowly,
being below the energy threshold of conventional detectors, but boosted DM is
able to deposit sufficient energy. Interest in such models was further
stimulated by excess events (now attributed to background) observed by the
XENON1T experiment in 2010 [11]. These were the motivations for our recent
work on the subject, Refs. [1, 12].

##  4 Revised limit on $\Lambda$

One of our goals was to better quantify the limit on $\Lambda$ from vacuum
decay into ghosts plus photons, improving on the original order-of-magnitude
estimate. First, we used the exact expression for the matrix element, which
takes the rather aesthetic form $|{\cal M}|^{2}=(tu/s)^{2}/8m_{p}^{4}$ [13].
Second, we computed the phase space integral accurately to obtain the full
photon spectrum. Third, we accounted for reshifting of the spectrum using an
appropriate Boltzmann equation, with the result shown in Fig. 3. Finally, we
used more recent observational data from the COMPTEL experiment. The new limit
on the cutoff is

$\Lambda<18\,{\rm MeV}\,,$ (5)

six times less stringent than the original estimate.

Figure 3: Left: spectrum of photons from decay of the vacuum, to
$\gamma\gamma\phi\phi$. Right: spectra of sterile neutrinos from vacuum decay
to $\nu_{s}\bar{\nu}_{s}\phi\phi$.

##  5 Contstraints on phantom fluid production

For definiteness, we consider a sterile neutrino $\nu_{s}$ as the light DM
candidate coupled to phantoms, but similar results would hold for scalar DM.
However the gravitational coupling, analogous to Fig. 2, turns out to be too
weak (by some 9 orders of magnitude) to be cosmologically relevant. Therefore
we consider the possibility that the decay to $\phi\phi{\nu_{s}}\bar{\nu_{s}}$
is mediated by a more strongly coupled intermediate particle, such as a heavy
vector or scalar boson. For definiteness I will assume the vector case here
(scalar exchange gives similar results); by integrating out the heavy
particle, it can be described by an effective interaction,

${\cal L}\quad\ni\quad{i\over
M_{v}^{2}}(\phi^{*}\\!\\!\stackrel{{\scriptstyle\leftrightarrow}}{{\partial}}_{\mu}\\!\phi)\,\bar{\nu_{s}}\gamma^{\mu}{\nu_{s}}\,,$
(6)

where the ghost is complex and the neutrino is Dirac, consistent with a U(1)
gauge interaction. We will show that $M_{v}$ as large as the GUT scale can
lead to observable effects, if $\Lambda$ saturates the bound (5).

Figure 4: Left: equation of state of the ghost fluid as a function of
redshift. Band indicates variation with respect to $m_{\nu_{s}}/\Lambda$ for
vector and scalar mediator models. Right: growth of the phantom fluid energy
density with scale factor, for several $m_{\nu_{s}}/\Lambda$ values.

For simplicity we take the phantoms to be massless. Then there are two
qualitatively different regimes: ${\nu_{s}}$ is massive or (nearly) massless.
The first case is straightforward since it doesn’t require cosmological
perturbation theory. Instead the main effects are via the time-dependent
energy densities of the ghosts and ${\nu_{s}}$s produced from vacuum decay.
They are described by coupled Boltzmann equations,

$\displaystyle{d\rho_{\phi}\over dt}+4H\rho_{\phi}$ $\displaystyle=$
$\displaystyle-\Gamma_{\rho},\quad{d\rho_{{\nu_{s}}}\over
dt}+3H(1+w_{{\nu_{s}}})\rho_{{\nu_{s}}}=+\Gamma_{\rho}\,,$ (7)

where the rate of energy density production of $\nu_{s}$ is

$\Gamma_{\rho}=\int\prod_{i}{d^{\,3}p_{i}\over(2\pi)^{3}2|E_{i}|}(2\pi)^{4}\delta^{(4)}\Big{(}\sum_{i}p_{i}^{\mu}\Big{)}|{\cal
M}|^{2}\cdot(E_{{\nu_{s}}}+E_{\bar{\nu_{s}}})$ (8)

with matrix element $|{\cal M}|^{2}=(2/M_{v}^{2})(s^{2}-(t-u)^{2})$, and the
Hubble parameter is

$H^{2}=H_{0}^{2}\times(\rho_{\rm
rad}+\rho_{m}+\rho_{\Lambda}+{\color[rgb]{1,0,0}\rho_{\phi}+\rho_{{\nu_{s}}}})/\rho_{c}\,,$
(9)

where $\rho_{c}$ is the critical density. Because of the terms in red, the
Boltzmann equations (7) are nonlinear, but in practice, their sum must remain
relatively small to avoid undue modifications of standard cosmology, and then
simple analytic approximate solutions can be found.

One can see that if ${\nu_{s}}$ is massless, $w_{{\nu_{s}}}=1/3$, then the
phantom fluid energy density $\rho_{g}\equiv\rho_{\phi}+\rho_{{\nu_{s}}}=0$ is
an exact, rather uninteresting solution. But if $m_{\nu_{s}}>0$, the two
species redshift differently, leading to net production of $\rho_{g}$. A
potential complication is that $w_{{\nu_{s}}}$ could be time-dependent, since
the ${\nu_{s}}$s are being produced at a constant rate:

$w_{{\nu_{s}}}={p_{{\nu_{s}}}\over\rho_{{\nu_{s}}}}={\int
d^{\,3}p\,f_{{\nu_{s}}}(p,t)\,p^{2}/3E\over\int
d^{\,3}p\,f_{{\nu_{s}}}(p,t)\,E}\sim\frac{1}{3}-\left(m_{{\nu_{s}}}\over\Lambda\right)^{2}+\cdots$
(10)

hence their distribution function $f_{\nu_{s}}$ is time-dependent. However we
find that this dependence cancels nearly exactly in the ratio, so that
$w_{\nu_{s}}$ is only a function of $m_{\nu_{s}}/\Lambda$, which greatly
simplifies the analysis.

In effect, $\rho_{g}$ behaves as an exotic form of dark energy, different from
that of the classical phantom condensate, with effective equation of state
$w_{\rm eff}=(p_{\phi}+p_{\nu_{s}})/(\rho_{\phi}+\rho_{\nu})$ as shown in Fig.
4 (left). It violates the DEC, $-1.5\leq w_{\rm eff}\lesssim-1.3$. Hence,
unlike matter, $\rho_{g}$ grows with time, as shown in Fig. 4 (right).

Figure 5: Left: Lower limit on the new physics scale $M_{v}$ versus $\Lambda$
from phantom fluid production, for $m_{\nu_{s}}>0$. Right: same for
$m_{\nu_{s}}=0$. Brown regions are excluded by the technical assumption that
the fluids do not become degenerate (for which our analysis is not valid).

For the massless ${\nu_{s}}$ case, no constraint arises at zeroth order in
cosmological perturbation theory, even with arbitrarily fast growth of
$\rho_{\phi}$ or $\rho_{\nu_{s}}$ individually, since they cancel exactly.
However, their perturbations do not exactly cancel, since their distribution
functions are different (Bose-Einstein versus Fermi-Dirac) and their phase
space limits differ (only ghosts are restricted to have $|\vec{p}|<\Lambda$).
Hence their first-order perturbations can grow and impact the integrated
Sachs-Wolfe effect and lead to constraints from the cosmic microwave
background (CMB).

Figure 6: Left: upper limit on $\Omega_{\nu_{s}}$ versus $\Lambda$ for the
case of massless ${\nu_{s}}$. Right: direct detection upper limit on cross
section for ${\nu_{s}}$-$e$ scattering times ${\nu_{s}}$ density versus
$\Lambda$, from XENONnT null results.

For both massive and massless $\nu_{s}$, we performed a Monte Carlo search of
the parameter space using CosmoMC [14] to derive limits from combined
cosmological data sets on the new physics parameter $M_{v}$, as well as on the
phenomenological quantity $\Omega_{g}$ (the present contribution of the ghost
fluid to the cosmological density) as functions of the cutoff scale $\Lambda$.
$\Omega_{g}$ cannot exceed $0.2$, so disappointingly it is not a viable
replacement for the cosmological constant $\Omega_{\Lambda}\gtrsim 0.5$, and
could only provide a fraction of the total dark energy. The ensuing lower
limits on $M_{v}$ versus $\Lambda$ are shown in Fig. 5 for massive
$m_{\nu_{s}}$. For $m_{\nu_{s}}=0$, the density of $\nu_{s}$ can grow to
enormous values, $\Omega_{\nu_{s}}\sim 10^{25}$ since it is canceled by the
ghosts. The limit is shown in Fig. 6 (left).

We investigated the extent to which the phantom fluid model can ameliorate the
Hubble tension, as well as the $S_{8}$ tension ($S_{8}$ being the parameter
describing the amplitude of fluctuations at the scale $8h^{-1}$ Mpc). As
expected on general grounds, this being a late-dark energy model, it does not
provide a full resolution, but it reduces the discrepancy to the $3\sigma$
level. The main source of the tension is the baryon acoustic oscillation (BAO)
data. If one omits these sets, excellent consistency can be obtained.

##  6 Direct detection

With the possibility of such a huge exotic dark radiation component, one is
led to wonder if it could have observable consequences. Suppose it was coupled
to electrons in a similar manner to Eq. (6),

${\cal L}\quad\ni\quad{i\over M^{\prime
2}}(\bar{e}\gamma_{\mu}e)(\bar{\nu_{s}}\gamma^{\mu}{\nu_{s}})\,,$ (11)

where we allow for a different interaction strength $M^{\prime}$. We choose
electrons since recent experimental searches for light dark matter have been
taking advantage of $e$-DM scattering: it is easier for a light particle to
impart energy to electrons than to nuclei. The dark radiation of massless
$\nu_{s}$ carries energy $\sim\Lambda$ which can be far greater than that of
thermal radiation, and easily above the energy thresholds $\sim 20$ eV of
detectors such as DAMIC [15]. In Ref. [12] we calculated the event rate for
such interactions in xenon and silicon detectors. A model-independent
constraint on the scattering cross section of ${\nu_{s}}$ on electrons times
the ${\nu_{s}}$ density can be derived from the null results of the XENONnT
experiment [16], which exclude the putative signal that was previously
observed by XENON1T [11]. It is shown in Fig. 6 (right).

Translating the XENONnT constraint into model parameters, we find an upper
bound on the product $M_{v}M^{\prime}$, corresponding to the product
$\sigma_{{\nu_{s}}e}n_{\nu_{s}}$. The signal turns out to be undetectably
small for cosmologically allowed values if we assume $M^{\prime}=M_{v}$; one
needs $M^{\prime}>M_{v}$ to saturate the bound. This is illustrated in Fig. 7
(left). The DAMIC experiment in SNOLAB has observed excess events that can be
interpreted as scattering on electrons, with a spectrum shown in Fig. 7
(right).

The phantom fluid model has essentially two parameters relevant for direct
detection: $\Lambda$, determining the shape of the spectrum and
$\sqrt{M_{v}M^{\prime}}$, determining its normalization; the value of
$m_{\nu_{s}}/\Lambda$ has a relatively small effect. We find an excellent fit
to the observed data by taking $\Lambda=10\,$keV and
$\sqrt{M_{v}M^{\prime}}=1.3\times 10^{5}\,$GeV, which is not far from the
XENONnT lower limit. It must be raised to the 8th power to get
$\sigma_{{\nu_{s}}e}n_{\nu_{s}}$, implying the signal is $1.3^{8}\sim 8$ times
lower than the current XENONnT sensitivity. There is a data bin at 0.05 keV
that is essentially empty. We confirm that the predicted signal falls to
negligible values in this bin, due to the binding energy of silicon preventing
excitations.

Figure 7: Left: XENONnT limit on $\sqrt{M_{v}M^{\prime}}$, compared to
previous cosmological limits assuming $M^{\prime}=M_{v}$. Right: massless
${\nu_{s}}$ fit to the DAMIC excess signal, with $\Lambda=10\,$keV and
$\sqrt{M^{\prime}M_{v}}=1.3\times 10^{5}$ GeV.

##  7 Conclusions

Phantom fluids through vacuum decay are an inevitable consequence of phantom
dark energy models, although it is necessary to couple the ghosts more
strongly than gravity to make them cosmologically relevant. If one does this,
for (massless) ghosts and light dark matter of mass $\lesssim 18$ MeV, it
gives an exotic dark energy component with equation of state $-1.5<w<-1.3$
that violates the dominant energy condition. Sadly, cosmological data prevents
it from contributing more than $\Omega_{g}\sim 0.2$ to the total energy
budget. Nevertheless, even a small amount will lead to a regularized big rip
in the future—regularized, since the theory must revert to being normal at
scales above the cutoff. $\Lambda\sim 10$ keV is sufficient to rip apart bound
objects like the sun in a future epoch. Apart from these annoying features,
phantom fluids can also do something useful by providing a new target for
detection of boosted light dark matter or radiation.

This work was supported by NSERC (Natural Sciences and Engineering Research
Council, Canada).

## References

  * [1] J. M. Cline, M. Puel, T. Toma and Q. S. Wang, “Phantom fluid cosmology: Impact of a phantom hidden sector on cosmological observables,” Phys. Rev. D 108, no.9, 095033 (2023) [arXiv:2308.12989 [hep-ph]].
  * [2] D. N. Spergel et al. [WMAP], “First year Wilkinson Microwave Anisotropy Probe (WMAP) observations: Determination of cosmological parameters,” Astrophys. J. Suppl. 148, 175-194 (2003) [arXiv:astro-ph/0302209 [astro-ph]].
  * [3] S. Hannestad and E. Mortsell, “Probing the dark side: Constraints on the dark energy equation of state from CMB, large scale structure and Type Ia supernovae,” Phys. Rev. D 66, 063508 (2002) [arXiv:astro-ph/0205096 [astro-ph]].
  * [4] A. Melchiorri, L. Mersini-Houghton, C. J. Odman and M. Trodden, “The State of the dark energy equation of state,” Phys. Rev. D 68, 043509 (2003) [arXiv:astro-ph/0211522 [astro-ph]].
  * [5] R. R. Caldwell, “A Phantom menace?,” Phys. Lett. B 545, 23-29 (2002) [arXiv:astro-ph/9908168 [astro-ph]].
  * [6] R. L. Workman et al. [Particle Data Group], “Review of Particle Physics,” PTEP 2022, 083C01 (2022)
  * [7] S. M. Carroll, M. Hoffman and M. Trodden, “Can the dark energy equation-of-state parameter $w$ be less than $-1$?,” Phys. Rev. D 68, 023509 (2003) [arXiv:astro-ph/0301273 [astro-ph]].
  * [8] J. M. Cline, S. Jeon and G. D. Moore, “The Phantom menaced: Constraints on low-energy effective ghosts,” Phys. Rev. D 70, 043543 (2004) [arXiv:hep-ph/0311312 [hep-ph]].
  * [9] E. Di Valentino, O. Mena, S. Pan, L. Visinelli, W. Yang, A. Melchiorri, D. F. Mota, A. G. Riess and J. Silk, “In the realm of the Hubble tension—a review of solutions,” Class. Quant. Grav. 38, no.15, 153001 (2021) [arXiv:2103.01183 [astro-ph.CO]].
  * [10] K. Agashe, Y. Cui, L. Necib and J. Thaler, “(In)direct Detection of Boosted Dark Matter,” JCAP 10, 062 (2014) [arXiv:1405.7370 [hep-ph]].
  * [11] E. Aprile et al. [XENON], “Excess electronic recoil events in XENON1T,” Phys. Rev. D 102, no.7, 072004 (2020) [arXiv:2006.09721 [hep-ex]].
  * [12] J. M. Cline, M. Puel and T. Toma, “Boosted dark matter from a phantom fluid,” Phys. Lett. B 848, 138377 (2024) [arXiv:2308.01333 [hep-ph]].
  * [13] N. Bernal, M. Dutra, Y. Mambrini, K. Olive, M. Peloso and M. Pierre, “Spin-2 Portal Dark Matter,” Phys. Rev. D 97, no.11, 115020 (2018) [arXiv:1803.01866 [hep-ph]].
  * [14] A. Lewis and S. Bridle, “Cosmological parameters from CMB and other data: A Monte Carlo approach,” Phys. Rev. D 66, 103511 (2002) [arXiv:astro-ph/0205436 [astro-ph]].
  * [15] A. Aguilar-Arevalo et al. [DAMIC, DAMIC-M and SENSEI], “Confirmation of the spectral excess in DAMIC at SNOLAB with skipper CCDs,” [arXiv:2306.01717 [astro-ph.CO]].
  * [16] E. Aprile et al. [XENON], “Search for New Physics in Electronic Recoil Data from XENONnT,” Phys. Rev. Lett. 129, no.16, 161805 (2022) [arXiv:2207.11330 [hep-ex]].