# Effects of dilution in a 2D topological magnon insulator

Miguel S. Oliveira  Centro de Física das Universidades do Minho e do Porto,
LaPMET, Departamento de Física e Astronomia, Faculdade de Ciências,
Universidade do Porto, 4169-007 Porto, Portugal  T. V. C. Antão  Laboratório
de Instrumentação e Física Experimental de Partícuals (LIP), University of
Minho, 4710-057 Braga, Portugal  Eduardo V. Castro  Centro de Física das
Universidades do Minho e do Porto, LaPMET, Departamento de Física e
Astronomia, Faculdade de Ciências, Universidade do Porto, 4169-007 Porto,
Portugal  Beijing Computational Science Research Center, Beijing 100084, China
Nuno Peres  Centro de Física das Universidades do Minho e do Porto (CF-UM-UP)
e Departamento de Física, Universidade do Minho, P-4710-057 Braga, Portugal
International Iberian Nanotechnology Laboratory (INL), Av Mestre José Veiga,
4715-330 Braga, Portugal  POLIMA—Center for Polariton-driven Light-Matter
Interactions, University of Southern Denmark, Campusvej 55, DK-5230 Odense M,
Denmark

###### Abstract

We study the effect of diluting a two-dimensional ferromagnetic insulator
hosting a topological phase in the clean limit. By considering the
ferromagnetic Heisenberg model in the honeycomb lattice with second nearest-
neighbor Dzyanshikii-Moriya interaction, and working in the linear spin-wave
approximation, we establish the topological phase diagram as a function of the
fraction $p$ of diluted magnetic atoms. The topological phase with Chern
number $C=1$ is robust up to a moderate dilution $p_{1}^{*}$, while above a
higher dilution $p_{2}^{*}>p_{1}^{*}$ the system becomes trivial.
Interestingly, both $p_{1}^{*}$ and $p_{2}^{*}$ are below the classical
percolation threshold $p_{c}$ for the honeycomb lattice, which gives physical
significance to the obtained phases. In the topological phase for
$p<p_{1}^{*}$, the magnon spectrum is gapless but the states filling the
topological, clean-limit gap region are spatially localized. For energies
above and below the region of localized states, there are windows composed of
extended states. This is at odds with standard Chern insulators, where
extended states occur only at single energies. For dilutions
$p_{1}^{*}<p<p_{2}^{*}$, the two regions of extended states merge and a
continuum of delocalized states appears around the middle of the magnon
spectrum. For this range of dilutions the Chern number seems to be ill defined
in the thermodynamic limit, and only for $p>p_{2}^{*}$, when all states become
localized, the system shows $C=0$ as expected for a trivial phase. Replacing
magnetic with non-magnetic atoms in a systems hosting a magnon Chern insulator
in the clean limit puts all the three phases within experimental reach.

##  I Introduction

Electronics based on the spin degree of freedom (Spintronics) harbor well-
known advantages from the present silicon-based technologies (Rajput  _et al._
, 2022). Magnonics (Lenk  _et al._ , 2011), a subfield of Spintronics,
addresses the use of spin waves (magnons) to transmit and process information
(Serga  _et al._ , 2004; Demidov  _et al._ , 2009; Jorzick  _et al._ , 2002;
Podbielski  _et al._ , 2006). By replacing charge currents by spin currents,
it is possible to avoid Joule heating, making these technologies of great
practical interest (Barman  _et al._ , 2021). These waves have already been
measured in Yttrium Iron Garnet thin-films (YIG) (Jungfleisch  _et al._ ,
2015). On the other hand, the study of magnons furnishes fundamental insight
about collective excitations and low-energy properties of quantum magnetic
systems.

The observation of the anomalous thermal Hall Effect for magnons – intrinsic
transverse heat transport response to a longitudinal temperature gradient – in
the insulating ferromagnet Lu_2V_2O_7 with a pyrochlore structure (Onose  _et
al._ , 2010; Hirschberger  _et al._ , 2015a), and more recently in a Kagomé
magnet (Hirschberger  _et al._ , 2015b), confirmed that magnonic systems may
host topological properties (Chisnell  _et al._ , 2015) just like other
systems with bosonic quasiparticles, such as photons (Rikken and van Tiggelen,
1996) and phonons (Zhang  _et al._ , 2010). Topology was first introduced to
condensed matter systems through the 2D electron gas in the quantum Hall
regime (Thouless  _et al._ , 1982) and through the electronic models of
Haldane (Haldane, 1988) and Kane-Mele (Kane and Mele, 2005) in the honeycomb
lattice. Soon after, topological insulators, a novel type of materials
characterized by an insulating bulk and metallic boundary states, became a hot
spot of research for their exotic properties and possible applications (Fu,
2010; Vijay  _et al._ , 2015; Dennis  _et al._ , 2002). The presence of edge
states in these materials is ensured by the bulk-edge correspondence, a
topological property which dictates that edge states properties are deeply
connected to the bulk and its symmetries. For small disorder which does not
break essential symmetries, these edge states are immune to back-scattering
making them the ideal transport states (Lado  _et al._ , 2015; Murakami,
2011). In the case of neutral quasi-particles such as magnons, the non trivial
topology may come from the spin-orbit interaction, imposed by Dzyaloshinskii-
Moriya (DM) (Moriya, 1960; McClarty, 2022) contribution to the exchange
interaction between localized magnetic moments.

A key property of topological insulators is their robustness to moderate
disorder (Wu  _et al._ , 2017; Li, 2019; Xiao  _et al._ , 2010), making them
attractive for real applications where disorder cannot be avoided. Disorder
may even be an essential ingredient for the experimental observation of a
topological response, as is the case of the quantum Hall effect (Kramer and
MacKinnon, 1993). In some cases, disorder may also benefit the appearance of
topological properties in regimes where the unperturbed model is trivial. The
topological Anderson insulator is a well known example of this class (Groth
_et al._ , 2009; Li  _et al._ , 2009; Gonçalves  _et al._ , 2018), which has
recently been realized in cold atoms (Meier  _et al._ , 2018), photonic
crystals (Stützer  _et al._ , 2018), and electrical circuits (Zhang  _et al._
, 2019). Disorder induced nontrivial topology has also been proposed (Li  _et
al._ , 2020; Yang  _et al._ , 2021a; Agarwala  _et al._ , 2020; Wang  _et al._
, 2021; Peng  _et al._ , 2022; Lóio  _et al._ , 2023) and observed (Zhang  _et
al._ , 2020) in higher-order topological insulators. For topological magnon
insulators, however, the role of disorder has been much less appreciated.

A particularly interesting type of disorder in magnetic systems is dilution,
which can be achieved by replacing a certain fraction $p$ of magnetic atoms by
non-magnetic ones. Diluted quantum magnets have been studied previously, both
diluted ferromagnets (Edwards and Jones, 1971) and antiferromagnets (Wan  _et
al._ , 1991). Since long range order is limited by the classical percolation
threshold of the underlying lattice, the question then was whether a quantum
phase transition could be induced by dilution prior to the expected classical
percolation transition. This question is even more relevant in 2D due to the
low dimensionality, but a stochastic series expansion approach for the diluted
2D antiferromagnetic Heisenberg model clearly showed that the percolating
cluster at $p_{c}$ (the critical classical dilution threshold) showcases long-
range order (Sandvik, 2001, 2002). It was later shown that critical exponents
involving dynamical correlations are different from the classical percolation
values, even though the transition is driven by the underlying classical
percolation (Vojta and Schmalian, 2005).

In this article we study the interplay between nontrivial topology and
magnetic dilution in 2D topological magnon insulators. We address in
particular the question of whether a topological transition takes place with
increasing dilution and whether or not there is a relation with the classical
percolation threshold of the underlying lattice. Since topological properties
are deeply related with localization and transport properties, we also study
the effect of dilution on localization of the magnon states. Moreover,
dilution falls in the strong disorder classification even for small
percentages of vacancies. There is no consensual answer to whether
localization properties for this type of disorder follow a one parameter
scaling theory (Abrahams  _et al._ , 1979), like other types of disorder such
as Anderson disorder (Anderson, 1958). Using the analysis of the transmission
coefficient associated to a simple nearest neighbor model in a diluted
lattice, it was shown that for increasing disorder and for a fixed finite
energy, the system undergoes a transition from delocalized states to power-law
localized states, and for large dilution the system is exponentially localized
(Islam and Nakanishi, 2008; Dillon and Nakanishi, 2014). We anticipate that
the complexity of the localization properties for this type of disorder may
introduce interesting topological behaviors.

The paper is organized as follows: In Sec. II, we introduce the model and the
methods used to characterize the system. The topological, spectral, and
localization properties are discussed in Sec. III. In Sec. IV we provide a
thorough discussion of the obtained results. The key results are summarized in
Sec. V and some conclusions are drawn. In Appendix A we compute the classical
site percolation threshold $p_{c}$ for the diluted honeycomb lattice with
first and second neighbor connections.

##  II Model and Methods

###  II.1 Magnons Model

####  II.1.1 Clean case

We consider the Heisenberg model for localized spins with DM (Moriya, 1960)
and Zeeman exchange interaction for a ferromagnetic system in the honeycomb
lattice with the Hamiltonian written as

${\cal H}=J\sum_{\left\langle ij\right\rangle}{\bf S}_{i}\cdot{\bf
S}_{j}+\sum_{\left\langle\left\langle ij\right\rangle\right\rangle}{\bf
D}_{ij}\cdot\left[{\bf S}_{i}\times{\bf S}_{j}\right]-B\sum_{i}S_{i}^{z},$ (1)

where $\left\langle\dots\right\rangle$ and
$\left\langle\left\langle\dots\right\rangle\right\rangle$ refers respectively
to the summation over the NN and NNN elements of the lattice and the coupling
${\bf D}_{ij}$ in the honeycomb lattice is given by $\nu_{ij}D\hat{z}$, with
$\nu_{ij}=2/\sqrt{3}\left({\bf d}_{1}\times{\bf d}_{2}\right)\cdot{\bf z}=\pm
1$, where ${\bf d}_{1},{\bf d}_{2}$ are two unit vectors along the bonds
connecting the NNN $\left\langle\left\langle i,j\right\rangle\right\rangle$.
In the following analysis we fix $\hbar=1$. Considering the low temperature
regime, it is possible to study this model by focusing on deviations from the
ferromagnetic ground state, and using Holstein-Primakoff formalism (Holstein
and Primakoff, 1940) to define the creation $(a_{i}^{\dagger})$ and
annihilation $(a_{i})$ bosonic operators associated to spin deviations
$a_{i}^{\dagger}a_{i}=n_{i}=S-S_{i}^{z}$, where $n_{i}$ is the number operator
of the site $i$. The spin operators may then be rewritten as

$\displaystyle S_{i}^{-}$
$\displaystyle=\sqrt{2S}a_{i}^{\dagger}\sqrt{1-\frac{a_{i}^{\dagger}a_{i}}{2S}}\approx\sqrt{2S}a_{i}^{\dagger}$
$\displaystyle S_{i}^{+}$
$\displaystyle=\sqrt{2S}\sqrt{1-\frac{a_{i}^{\dagger}a_{i}}{2S}}a_{i}\approx\sqrt{2S}a_{i},$
(2)

where we only kept the zeroth order term, disregarding interactions between
bosons, since we will focus only in the low energy regime where
$a_{i}^{\dagger}a_{i}\ll 2S$, the linear spin-wave approximation. Bosonic
statistics do not have a limit on the occupation factor, however, in this
approximation the maximum occupation number cannot exceed the number of
deviations needed to surpass the minimum spin value.

We can rewrite the Hamiltonian with respect to bosonic operators as

$\displaystyle{\cal H}=$ $\displaystyle 2JS\sum_{\left\langle
ij\right\rangle}\left(a_{i}^{\dagger}a_{i}-a_{i}^{\dagger}a_{j}\right)$
$\displaystyle+2iSD\sum_{\left\langle\left\langle
ij\right\rangle\right\rangle}\nu_{ij}a_{i}^{\dagger}a_{j}-B\sum_{i}a_{i}^{\dagger}a_{i},$
(3)

where we have omitted the constant term associated to the ferromagnetic ground
state energy. Making use of translational invariance to write ${\cal H}$ in
the momentum basis, it is possible to show that the spectrum of the model
follows the dispersion relation

$E({\bf k})=2JS\left(3\pm\sqrt{\left|f({\bf k})\right|^{2}+m({\bf
k})^{2}}\right),$ (4)

where

$\displaystyle m({\bf k})$ $\displaystyle=-2\frac{D}{J}\left(\sin({\bf
k}\cdot{\bf b}_{1})-\sin\left({\bf k}\cdot({\bf b}_{1}-{\bf
b}_{2})\right)-\sin({\bf k}\cdot{\bf b}_{2})\right)$ $\displaystyle f({\bf
k})$ $\displaystyle=1+\exp(i{\bf b}_{1}\cdot{\bf k})+\exp(i{\bf
b}_{2}\cdot{\bf k}),$ (5)

${\bf b}_{1},{\bf b}_{2}$ being the reciprocal lattice vectors.

####  II.1.2 Non-trivial topology in the clean limit

It becomes clear in Eq. (3) that in the lowest order spin-wave approximation
the system is described by a single-particle tight binding model. Apart from
the absence of a trivial mass, this bosonic model is equivalent to the
ferminonic system introduced by Haldane (Haldane, 1988), hosting topological
phases characterized by the topological index Chern number $C$. Given an
eigenstate $|\Psi_{{\bf k}}^{n}\rangle$ of our system associated to the band
$n$, the Berry curvature may be defined as

$\Omega_{xy}^{(n)}({\bf
k})=-2\text{Im}\left(\left\langle\partial_{k_{x}}\Psi_{{\bf
k}}^{n}\Big{|}\partial_{k_{y}}\Psi_{{\bf k}}^{n}\right\rangle\right).$ (6)

The Chern number corresponds to the integral of Berry curvature in Eq. (6)
over the first Brilloin zone,

$\displaystyle C_{n}$
$\displaystyle=\frac{1}{2\pi}\iint_{1ZB}\Omega_{xy}^{n}d^{2}{\bf k}.$ (7)

In the trivial phase, both bands have a null Chern number while in the
topological phase the valence and conduction band have $C=\pm 1$, keeping the
overall Chern number null.

Even though this model is a bosonic model and the notion of band filling does
not apply like in the equivalent fermionic model, the topological properties
of a non-interacting system depend exclusively on the single particle
Hamiltonian, which is independent of the statistical descriptions of the
particles which inhabit the model. In the bosonic case, the transport
properties at low temperature will not be described by the states at the
middle of the spectrum, but one can always force an excitation to have a
certain energy with a external stimulus. For example, it is possible to excite
the edge states in a magnon model, provided that the system is in a
topological phase (Shindou  _et al._ , 2013; Zhang  _et al._ , 2013a; Mook
_et al._ , 2014).

####  II.1.3 Diluted Magnons

Diluted lattices are simulated by assuming that each site is independently and
randomly occupied with a localized spin with probability $1-p$ and it is empty
with probability $p$. Whenever a site is assigned as a vacancy, it is removed
from the tight binding basis. A _cluster_ is defined as a set of occupied
neighboring sites, which in the terminology of graph theory represents a
connected graph. In the spirit of the tight binding approach, we define the
neighborhood (set of edges of one vertex) of a site as the set of occupied
sites which share non-zero hopping terms. In every disorder configuration we
choose random occupied sites and in the end retain solely the largest
connected cluster since this will be the one with physical relevance. After
collecting the cluster, we build the Hamiltonian of the system with the basis
of atomic orbitals of occupied sites only belonging to the largest connected
cluster.

Figure 1: To obtain the largest cluster of a diluted lattice, we start from a
clean honeycomb lattice (a), where the red and black colors make the bipartite
nature of the lattice explicit. Running through all the sites in the lattice,
we randomly assign the site as vacant with probability $p$ and in this case we
remove it from our graph (b). After choosing all the vacant sites, we collect
the different clusters of disconnected sites, as shown in panel (c) with
different colors. The model we will consider has next-to-nearest neighbors
hoppings, hence we consider that for 2 sites to belong to the same cluster
they should be either first neighbors or second neighbors, considering
periodic boundary conditions. From all the clusters, we retain just the one
with more sites and encode it in a tight binding Hamiltonian (d).

Periodic boundary conditions must be employed in order to calculate bulk
properties. This modification of the lattice implies occupied sites in two
distinct boundaries may still be connected and belong to the same cluster, as
exemplified in Fig. 1. A complete finite size honeycomb lattice is shown in
Fig. 1(a). In Fig. 1(b) we assign randomly the vacant sites, and in Fig. 1(c)
we cluster the disjoint sets of sites. By keeping only the largest set of
sites we obtain the largest connected cluster, shown in Fig. 1(d).

As one increases the percentage of vacancies in the simulated lattice, for
percentages larger than a critical value, the size of the largest cluster no
longer scales with the same dimensionality as the size of the clean system. It
is not possible to take the thermodynamic limit for $p>p_{c}$, as the
simulated system has no dependence on the original size of the clean system.
This critical value is the so-called classical percolation threshold $p_{c}$
and only depends on the lattice geometry (Suding and Ziff, 1999; Ziff, 1992).
A second and equivalent definition of this quantity, which highlights the
physical relevance of this property, is the following: the probability of
having a macroscopic path composed of connected sites on a simulated diluted
lattice is $0$ if $p>p_{c}$ and is $1$ if $p<p_{c}$. This implies that there
is no transport nor long range order for $p>p_{c}$. Although the classical
percolation threshold for the honeycomb lattice with nearest neighbors
connections is well known to be $p_{c}\approx 0.302957$ (Suding and Ziff,
1999), there are no results in the literature for the honeycomb with first and
second nearest neighbors connections. We determined this value to be
$p_{c}=0.640\pm 0.005$ (see Appendix A).

In the diluted regime, the Hamiltonian is rewritten as

$\displaystyle{\cal H}=$ $\displaystyle 2JS\sum_{\left\langle
ij\right\rangle}\eta_{i}\eta_{j}\left(a_{i}^{\dagger}a_{i}-a_{i}^{\dagger}a_{j}\right)$
$\displaystyle+2iSD\sum_{\left\langle\left\langle
ij\right\rangle\right\rangle}\eta_{i}\eta_{j}\nu_{ij}a_{i}^{\dagger}a_{j}-B\sum_{i}\eta_{i}a_{i}^{\dagger}a_{i},$
(8)

where dilution is introduced in the Hamiltonian with the random variables
$\eta_{i}$, which take the values $0$ or $1$ depending whether the site $i$
exists in the lattice being simulated or not. Notice that the first term of
the Hamiltonian in Eq. (8) introduces a dependency on the number of occupied
neighbors of one site to the respective on-site energy, which is not present
in the Haldane model. In this work we fixed $2JS\equiv t$, $2SD=0.1t$, $B=0$,
and use $t$ as the energy unit.

Disorder breaks translational invariance, meaning that the Hamiltonian
eigenstates will no longer be Bloch states with well defined Bloch momentum
and band index. Nevertheless, it is still possible to extend the formalism
applied in undiluted systems, and compute the Chern number with a real space
approach, as discussed next.

###  II.2 Methods

####  II.2.1 Chern number

The computation of the Chern number can be extended to systems with broken
translational invariance employing a super-cell approach with twisted boundary
conditions in real space. With Fukui’s method (Fukui  _et al._ , 2005), one
can replace the continuous integral in Eq. 7 by a discrete sum over boundary
twists $\bm{\theta}$. Additionally, with the coupling matrix method of Ref.
(Zhang  _et al._ , 2013b), we can reduce the sum to twist angles associated to
periodic boundary conditions. This method requires exact diagonalization of
two matrices with size proportional to the system’s size: the Hamiltonian, and
the coupling matrix (Zhang  _et al._ , 2013b).

To minimize finite size effects, we average the Chern number over disorder
realizations. Even though the Chern number is strictly quantized for each
disorder configuration, the averaged Chern number may not be. A finite size
scaling analysis is performed when needed in order to infer the thermodynamic
limit behavior.

####  II.2.2 Density of states

In order to build the gapless-gapped picture for the model, we used the Kernel
Polynomials Method (KPM) (Weiße  _et al._ , 2006) to compute the density of
states through a Chebyshev polynomials series expansion, equipped by the
Jackson’s Kernel. The linear scaling of this method with the system size
allowed us to explore systems of the order of $N=1024^{2}$ unit cells.

A finite cutoff of the series expansion introduces an artificial smoothing in
abrupt changing parts of the spectrum, which may lead to wrong interpretations
about the existence of a gap for small enough gaps. In these cases, an
analysis of the convergence of the DOS in the gap with respect to the cutoff
order of the series expansion was employed in hopes of minimizing the
uncertainty of the limits of the gapped-gapless phases.

####  II.2.3 Localization

The localization behavior of the system was explored using the standard tools:
level spacing statistics (Evers and Mirlin, 2008) to try to observe a
transition from a GUE (Gaussian Unitary Ensemble) distribution to a Poisson
distribution as states become localized, and Transfer Matrix Method (TMM)
(Kramer and MacKinnon, 1993) to compute the scaling of the
correlation/localization length across the spectrum for different dilution
percentages.

In order to identify the probability distribution of the level spacings, we
looked at the distribution of the ratio between consecutive level spacings,
$r_{n}=\min(s_{n},s_{n+1})/\max(s_{n},s_{n+1})$, where
$s_{n}=\epsilon_{n+1}-\epsilon_{n}$ are the level spacings. This quantity is
known to have better behavior in regions of the spectrum with low density of
states (Oganesyan and Huse, 2007; Atas  _et al._ , 2013).

In order to implement the TMM, we had to assign a high on-site energy
$V_{\text{vacancy}}$ to the vacant sites, plus a very small hopping term
$\delta t$ from adjacent sites to these vacant sites, instead of removing the
vacant sites from our tight binding basis. We used $V_{\text{vacancy}}=1000t$
and $\delta t/t=1/1000$. The limit $V_{\text{vacancy}}\to\infty,\delta t\to
0$, where both implementations of dilution should be equivalent, is discussed
in Sec. IV.

##  III Results

###  III.1 Topological phase diagram

As can be seen in Fig. 2, the topological phase is robust up to relatively
high percentages of vacancies, persisting with a well defined Chern number
$C=1$ up to $p\approx 15\%$. In the inset of Fig. 2 it is clearly seen that
for $p=15.3\%$, labeled with (a) in the inset, we have an averaged Chern
number $C\rightarrow 1$ as $1/L\to 0$. For large enough dilution, which we
overestimate as $p\gtrsim 23\%$, the system enters the trivial phase. In this
region of the phase diagram $C\rightarrow 0$ as $1/L\rightarrow 0$, as shown
for $p=23.7\%,$labeled with (d) in the inset of Fig. 2. The system becomes a
trivial insulator before reaching the classical percolation threshold for the
honeycomb lattice, $p_{c}=30\%$, and long before the threshold for the
honeycomb lattice with first and second nearest neighbors, $p_{c}=64\%$,
rendering the trivial phase physically achievable.

Figure 2: Averaged Chern number for different system sizes ($L$ is the total
number of unit cells in one spacial direction) and for different percentages
of vacancies. The standard deviation errors obtained by averaging over $500$
disorder realizations are shown for each data point. Inset: finite size
scaling analysis for the averaged Chern number at specific dilution values:
$p=15.3\%$ (a), $p=19.5\%$ (b), $p=21.6\%$ (c), and $p=23.7\%$ (d).

Between the topological and the trivial phases shown in Fig. 2, there is a
highlighted region for $15\%\lesssim p\lesssim 23\%$ where the averaged Chern
number crosses over from $C=1$ to $C=0$ as $p$ increases. Our numerical
analysis is not conclusive on whether this region shrinks to a single critical
dilution in the thermodynamic limit or remains finite. Nevertheless, a finite
size scaling analysis seems to indicate that the averaged Chern number does
not converge to quantized values as we approach the thermodynamic limit. This
is exemplified in the inset of Fig. 2 for $p$ values labeled with (b) and (c).
The presence of a crossover region which does not shrink to a single critical
dilution is also justified by the localization behavior of the system to be
discussed in Sec. III.2.2.

###  III.2 Spectral and localization properties

####  III.2.1 Spectral properties

In Fig. 3 we show the magnon DOS for different values of dilution $p$ obtained
using KPM. It is obvious that even for a small percentage of vacancies (above
$p\simeq 5\%$) the gap is completely filled. This implies that the topological
phase becomes gapless for small values of dilution. The gap remains closed as
disorder increases, even when the system becomes trivial.

The diagonal disorder seen in Eq. (8) destroys the particle-hole symmetry of
the spectrum for finite dilution, something that does not appear in analogous
fermionic systems, like the Haldane model. This asymmetry shifts the middle of
the spectrum to lower energies as dilution increases.

Figure 3: DOS for different percentages of dilution, including both the
topological and the trivial phases, as well as the crossover region. We used
systems with size $1024\times 1024$ and included $1024$ polynomials in the
KPM.

####  III.2.2 Localization properties

In order to provide a more complete analysis of the topological phase
transition, we have also characterized the localization properties of the
system for different values of $p$. The ratio between consecutive level
spacings $r_{n}$, introduced in Sec. II.2.3, allows to distinguish between
localized and extended states. For localized states $r_{n}$ should follow the
Poisson distribution, while for extended states in systems which break time
reversal symmetry we expect the Gaussian Unitary Ensemble (GUE) distribution
(Evers and Mirlin, 2008).

Figure 4: Finite size scaling for the average value $\left\langle
r\right\rangle$ of the distribution of the ratio between consecutive energy
spacings as a function of energy for increasing values of dilution percentages
$p$. We compare this average with the one expected for the GUE (horizontal
dashed blue line) and Poisson distribution (horizontal dashed red line). The
vertical dashed black line marks the middle of the spectrum, marking also the
upper limit of integration of the Berry curvature. The statistics are averaged
over a total number of $2048$ realizations of disorder for $N=32^{2}$ and
$256$ realizations of disorder for $N=64^{2}$, and the bin’s size used for the
statistics is of the order of $0.05t$.

In Fig. 4 we show the first moment of the distribution, $\left\langle
r\right\rangle$, indicating as dashed horizontal lines the expected values for
the Poisson and GUE distributions. With increasing dilution up to $p\lesssim
15\%$, there are two persistent regions, one below the middle of the spectrum
(signaled by the vertical dashed line) and the other above it, where
$\left\langle r\right\rangle$ seems to converge to GUE with increasing system
size. These two regions are separated by a region of localized states around
the middle of the spectrum, where $\left\langle r\right\rangle$ is converging
to the Poisson value. The presence of two finite energy regions where states
are delocalized is at odds with the behavior of quantum Hall systems (Kramer
and MacKinnon, 1993; Wang  _et al._ , 2014) and conventional Chern insulators
(Onoda and Nagaosa, 2003; Onoda  _et al._ , 2007; Castro  _et al._ , 2016;
Gonçalves  _et al._ , 2018), where extended states appear only at isolated
single energies. However, site dilution (or vacancies) is a particular type of
disorder where unconventional behavior is to be expected (Islam and Nakanishi,
2008; Dillon and Nakanishi, 2014). Yet another important remark to be made is
the presence of localized states around the middle of the spectrum, separating
the two regions of extended states. Since the spectrum is gapless in the
topological phase, it is the feature of having localized states filling the
gap which allows the system to retain the topological properties reminiscent
from the clean system.

If we further increase the percentage of vacancies, it can be seen in Fig. 4
that for $15\%\lesssim p\lesssim 23\%$ the states previously localized for
smaller $p$ values become extended, forming a continuum of extended states
around the middle of the spectrum. We believe this continuum to be responsible
for the crossover behavior observed in Fig. 2 regarding the topological
character of the system. The presence of this continuum prevents the usual
“levitation and pair annihilation” of the extended states carrying the
topological index (Laughlin, 1984; Onoda  _et al._ , 2007; Prodan, 2011) and
the existence of critical dilution at which the topological transition would
take place. As soon as dilution reaches the trivial phase, all the spectrum
becomes localized as is possible to observe in Fig. 2 for $p=25\%$.

Figure 5: Finite size scaling analysis of the normalized localization length
$\Lambda_{M}$ with the transverse size of the ribbon $M$ in the TMM for
different concentrations of vacancies. The longitudinal length used ensures
that the estimated error had converged to under $1\%$. The energy resolution
is of the order of $0.05t$.

In order to enrich our conclusions regarding the localization properties, we
computed the localization length with the TMM. As can be seen in Fig 5, for
small dilution up to $p\lesssim 15\%$ we observe two regions in which the
normalized localization length $\Lambda_{M}$ is scale invariant, corresponding
to delocalized states. These two regions are separated by a region where
$\Lambda_{M}$ decreases with $M$, signaling the presence of localized states.
This agrees well with the level spacing statistics results of Fig. 4. In the
crossover region $15\%\lesssim p\lesssim 23\%$, we also see the presence of a
continuum of extended states spanning a large region around the middle the
spectrum (note the collapse of $\Lambda_{M}$ for the larger $M$ values). This
confirms the unconventional localization properties of this system when
compared with conventional disordered Chern insulators (Gonçalves  _et al._ ,
2018). Upon reaching the trivial phase, all the states become localized, as is
possible to observe for $p\gtrsim 23\%$, where $\Lambda_{M}$ decreases with
the transverse size of the system, $M$, for the entire spectrum. Even though
the TMM results and those from level spacing statistics are in qualitative
agreement, there are noticeable quantitative differences which we discuss in
Sec. IV.

####  III.2.3 Spin-wave stability

As dilution increases, a finite density of states with negative energies
appears very close to zero, as can be clearly seen in Fig. 3. Negative energy
eigenstates would imply that the ferromagnetic state is not the true ground
state of the system, rendering the usage of Holstein-Primakoff transformation
not valid. This problem can be easily surpassed by making use of a finite
magnetic field, which would make a rigid shift of the spectrum. This shift
would open a energy gap between the ferromagnetic state and the lowest energy
excitations, thus stabilizing the spin-waves. The presence of single-ion
magnetic anisotropy would also give rise to a term in the Hamiltonian with
similar effects.

##  IV Discussion

The unconventional localization properties in the crossover region support the
conjecture that a topological index cannot be defined for this range of
dilutions since the system has delocalized magnons in a large window around
the middle of the spectrum (equivalent to a metallic phase in fermionic
systems). This is in agreement with the results of Fig 2 for the averaged
Chern number in the crossover region, which does not seem to converge to a
quantized value. The topological to trivial phase transition occurs when the
localized states, which separate the energy windows with extended states,
become localized. As dilution is increased, this enables a continuous flow of
Berry curvature between the two regions with delocalized states. We speculate
that the presence of a smooth crossover transition occurs because of the
appearance of this continuum near the transition.

In fermionic models belonging to the same symmetry class similar behavior was
already observed: the transition from the topological phase to the trivial one
on increasing disorder occurs through a metallic phase (Wang  _et al._ , 2015;
Qiao  _et al._ , 2016; Su  _et al._ , 2016; Yang  _et al._ , 2021b), for a
finite range of disorder values and not at a single critical disorder. The key
point to obtain the metallic behavior in these works is the presence of spin-
flip hopping disorder. Even though a clear relation between this latter type
of disorder and dilution is hard to establish, we notice that, in the clean
limit, the two sublattices of the honeycomb lattice ensure a two-component
wavefunction similar to spin-$1/2$ systems; around the middle of the spectrum
it is even possible to define a pseudo-spin from the sublattice degree of
freedom (Castro Neto  _et al._ , 2009). Pushing forward the connection,
nearest neighbor hopping disorder would mix the two components of the wave
function in a similar way as spin-flip processes mix spin-up and spin-down.
Finally, dilution may be regarded as an extreme case of hopping disorder,
where certain random hoppings are set to zero.

Another point worth mentioning is the reason for some quantitative differences
between the results from the TMM in Fig. 5 and the level spacing statistics in
Fig. 4. In the level spacing statistics we make use of exact diagonalization
of the Hamiltonian associated solely to the largest cluster for a given
configuration of vacancies, while for the TMM we cannot remove the small
isolated clusters from the computation. This means that this second method
opens the possibility of hosting for the same energy value, both localized
states in the independent clusters and localized states over the physical
cluster, which hybridize through the small hopping used to simulate vacancies
(see Sec. II.2.3). We believe this to be the reason for the absence of scaling
of $\Lambda_{M}$ with $M$ in sections of the spectrum for which the level
spacing statistics seems to tell us unambiguously that we are in the presence
of localized states (example around $E=4t$ for dilutions in the crossover
region).

##  V Conclusions

In this work we have studied the effect of magnetic dilution in a system which
hosts, in the clean limit, a topological magnon insulating phase. We have
shown that the topological phase is robust to the presence of magnetic
vacancies, surviving with a well defined Chern number $C=1$ up to a moderate
dilution $p_{1}^{*}$, which for the considered model we estimate to be around
$p_{1}^{*}\approx 15\%$. For high enough dilution $p_{2}^{*}>p_{1}^{*}$, the
system enters a trivial phase with $C=0$. For the considered model we estimate
$p_{2}^{*}\approx 23\%$. Interestingly, we found $p_{1}^{*},p_{2}^{*}<p_{c}$,
where $p_{c}$ is the classical percolation threshold, which for the honeycomb
lattice with nearest neighbor connections is $p_{c}\simeq 30\%$ (taking into
account the second neighbor connections present in our model it is
$p_{c}\simeq 64\%$).

Through an extensive characterization of the localization properties, we have
established that for $p<p_{1}^{*}$, when the system is in the topological
phase, the states that fill the clean limit gap are localized. The system
behaves effectively as a topological Chern insulator, with Berry curvature
carried by two energy regions of extended states above and below the region of
localized states. For $p>p_{2}^{*}$, when the system is in the trivial phase
with $C=0$, all states are localized. In the crossover region, for
$p_{1}^{*}<p<p_{2}^{*}$, we have found a continuum of extended states
(possibly critical) around the middle of the spectrum, in the region
previously occupied by the localized states in the topological phase for
$p<p_{1}^{*}$. Based on the localization properties and a finite scaling
analysis of the Chern number, we conjecture that, in the thermodynamic limit,
the topological index is not well defined in this region and topological
properties are not to be expected for $p_{1}^{*}<p<p_{2}^{*}$.

Replacing magnetic atoms with non-magnetic ones in a systems hosting a magnon
Chern insulator in the clean limit puts all the three phases within
experimental reach. Although most 2D ferromagnetic honeycomb materials, like
the compounds CrX_3 (X$\equiv$Cl,I and Br (Huang  _et al._ , 2017; Davis and
Narath, 1964; Narath and Davis, 1965; Samuelsen  _et al._ , 1971; Chen  _et
al._ , 2020)), do not exhibit a strong enough DM interaction to originate
nontrivial behavior, it is possible to mimic this interaction with the use of
an external time dependent electric field (Owerre, 2017) or to increase the
existing interaction via structural changes in the compound such as strain
(Koretsune  _et al._ , 2015) and the addition of more layers (Legrand  _et
al._ , 2022), allowing the clean limit topological properties to be engineered
(Antão and Peres, 2023). An interesting open question is whether dilution in
topological kagome ferromagnets leads to qualitatively similar physics (Zhuo
_et al._ , 2021, 2022).

## Acknowledgements

###### Acknowledgements.

MSO and EVC acknowledge partial support from Fundação para a Ciência e
Tecnologia (FCT-Portugal) through Grant No. UIDB/04650/2020. NMRP acknowledges
support from the Independent Research Fund Denmark (grant No. 2032-00045B),
the Danish National Research Foundation (Project No. DNRF165) and Fundação
para a Ciência e Tecnologia (FCT-Portugal) through Grant No. PTDC/FIS-
MAC/2045/2021.

##  Appendix A Classical percolation problem

In order to compute the classical site percolation threshold $p_{c}$ for the
diluted graphs where the vertices are the honeycomb lattice points and the
first and second neighbors are the edges, we computed the fraction of sites
which belong to the largest cluster of a simulated lattice for different clean
system sizes and for different percentages of vacancies. The results in the
Fig. 6 seem to indicate that for percentages below $64\%$, the fraction of
sites in the largest cluster is constant with respect to the increase of the
system size, and for percentages higher than $64\%$ there is a clear scaling
with $\propto 1/N_{\text{total}}$ showing that the largest cluster has a
constant size for these percentages. We estimate the value of $p_{c}$ to be
around $64\%$ since we have the transition from the two distinct behaviors at
this point.

Figure 6: Fraction of sites belonging to the largest cluster as a function of
the clean system size in one spacial direction for different values of $p$
averaged over $300$ configurations of dilution.

##  Appendix B Level spacing statistics sampling

For each disorder configuration, we obtained the set of energies, $e_{n}$,
with exact diagonalization. With this set, we computed the energy spacings
between consecutive levels, $s_{n}=e_{n}-e_{n}$, and used them to calculate
the ratio between consecutive energy spacings defined as
$r_{n}=\min(s_{n},s_{n+1})/\max(s_{n},s_{n+1})$. This last quantity is known
to behave better in regions where the DOS is small (Oganesyan and Huse, 2007;
Atas  _et al._ , 2013). After gathering the lists of $r_{n}$ over several
dilution configurations, we grouped them in equally spaced bins over the
spectrum to compute the average ratio $\left\langle r\right\rangle$ inside
each bin and obtain $\left\langle r\right\rangle$ as a function of energy and
compare it with the expected result for the GUE and Poisson distributions
(Fig. 4). The bin size was chosen small enough so we could probe the principal
changes in localization properties along the topological transition, but large
enough to enable $\left\langle r\right\rangle$ to have a smooth behavior.

## References

  * Rajput  _et al._ (2022) P. J. Rajput, S. U. Bhandari, and G. Wadhwa, Silicon (2022), 10.1007/s12633-021-01643-x.
  * Lenk  _et al._ (2011) B. Lenk, H. Ulrichs, F. Garbs, and M. Münzenberg, Physics Reports 507, 107 (2011).
  * Serga  _et al._ (2004) A. A. Serga, S. O. Demokritov, B. Hillebrands, and A. N. Slavin, Physical Review Letters 92, 117203 (2004).
  * Demidov  _et al._ (2009) V. E. Demidov, M. P. Kostylev, K. Rott, P. Krzysteczko, G. Reiss, and S. O. Demokritov, Applied Physics Letters 95, 112509 (2009).
  * Jorzick  _et al._ (2002) J. Jorzick, S. O. Demokritov, B. Hillebrands, M. Bailleul, C. Fermon, K. Y. Guslienko, A. N. Slavin, D. V. Berkov, and N. L. Gorn, Physical Review Letters 88, 047204 (2002).
  * Podbielski  _et al._ (2006) J. Podbielski, F. Giesen, and D. Grundler, Physical Review Letters 96, 167207 (2006).
  * Barman  _et al._ (2021) A. Barman, G. Gubbiotti, S. Ladak, A. O. Adeyeye, M. Krawczyk, J. Gräfe, C. Adelmann, S. Cotofana, A. Naeemi, V. I. Vasyuchka, B. Hillebrands, S. A. Nikitov, H. Yu, D. Grundler, A. V. Sadovnikov, A. A. Grachev, S. E. Sheshukova, J.-Y. Duquesne, M. Marangolo, G. Csaba, W. Porod, V. E. Demidov, S. Urazhdin, S. O. Demokritov, E. Albisetti, D. Petti, R. Bertacco, H. Schultheiss, V. V. Kruglyak, V. D. Poimanov, S. Sahoo, J. Sinha, H. Yang, M. Münzenberg, T. Moriyama, S. Mizukami, P. Landeros, R. A. Gallardo, G. Carlotti, J.-V. Kim, R. L. Stamps, R. E. Camley, B. Rana, Y. Otani, W. Yu, T. Yu, G. E. W. Bauer, C. Back, G. S. Uhrig, O. V. Dobrovolskiy, B. Budinska, H. Qin, S. van Dijken, A. V. Chumak, A. Khitun, D. E. Nikonov, I. A. Young, B. W. Zingsem, and M. Winklhofer, Journal of Physics: Condensed Matter 33, 413001 (2021).
  * Jungfleisch  _et al._ (2015) M. B. Jungfleisch, W. Zhang, W. Jiang, H. Chang, J. Sklenar, S. M. Wu, J. E. Pearson, A. Bhattacharya, J. B. Ketterson, M. Wu, and A. Hoffmann, Journal of Applied Physics 117, 17D128 (2015), https://doi.org/10.1063/1.4916027 .
  * Onose  _et al._ (2010) Y. Onose, T. Ideue, H. Katsura, Y. Shiomi, N. Nagaosa, and Y. Tokura, Science 329, 297 (2010), arXiv:1008.1564 [cond-mat].
  * Hirschberger  _et al._ (2015a) M. Hirschberger, J. W. Krizan, R. Cava, and N. Ong, Science (New York, N.Y.) 348, 106—109 (2015a).
  * Hirschberger  _et al._ (2015b) M. Hirschberger, R. Chisnell, Y. S. Lee, and N. P. Ong, Phys. Rev. Lett. 115, 106603 (2015b), arXiv:1502.05688 .
  * Chisnell  _et al._ (2015) R. Chisnell, J. S. Helton, D. E. Freedman, D. K. Singh, R. I. Bewley, D. G. Nocera, and Y. S. Lee, Physical Review Letters 115, 147201 (2015).
  * Rikken and van Tiggelen (1996) G. L. J. A. Rikken and B. A. van Tiggelen, Nature 381, 54 (1996).
  * Zhang  _et al._ (2010) L. Zhang, J. Ren, J.-S. Wang, and B. Li, Physical Review Letters 105, 225901 (2010).
  * Thouless  _et al._ (1982) D. J. Thouless, M. Kohmoto, M. P. Nightingale, and M. den Nijs, Physical Review Letters 49, 405 (1982), arXiv:arXiv:1011.1669v3 .
  * Haldane (1988) F. D. M. Haldane, Phys. Rev. Lett. 61, 2015 (1988).
  * Kane and Mele (2005) C. L. Kane and E. J. Mele, Physical Review Letters 95, 226801 (2005).
  * Fu (2010) L. Fu, Physical Review Letters 104, 056402 (2010).
  * Vijay  _et al._ (2015) S. Vijay, T. H. Hsieh, and L. Fu, Physical Review X 5, 041038 (2015).
  * Dennis  _et al._ (2002) E. Dennis, A. Kitaev, A. Landahl, and J. Preskill, Journal of Mathematical Physics 43, 4452 (2002), arXiv:quant-ph/0110143.
  * Lado  _et al._ (2015) J. L. Lado, N. Garcia-Martinez, and J. Fernandez-Rossier,  (2015), 10.48550/ARXIV.1502.07112.
  * Murakami (2011) S. Murakami, Journal of Physics: Conference Series 302, 012019 (2011).
  * Moriya (1960) T. Moriya, Physical Review 120, 91 (1960).
  * McClarty (2022) P. A. McClarty, Annual Review of Condensed Matter Physics 13, 171 (2022), arXiv:2106.01430 .
  * Wu  _et al._ (2017) B. Wu, J. Song, J. Zhou, and H. Jiang,  (2017), 10.48550/ARXIV.1711.10725.
  * Li (2019) Y. Li, Nature Physics 15, 4 (2019).
  * Xiao  _et al._ (2010) D. Xiao, M.-C. Chang, and Q. Niu, Reviews of Modern Physics 82, 1959 (2010).
  * Kramer and MacKinnon (1993) B. Kramer and A. MacKinnon, Rep. Prog. Phys. 56, 1469 (1993).
  * Groth  _et al._ (2009) C. W. Groth, M. Wimmer, A. R. Akhmerov, J. Tworzydło, and C. W. J. Beenakker, Physical Review Letters 103, 196805 (2009).
  * Li  _et al._ (2009) J. Li, R.-L. Chu, J. K. Jain, and S.-Q. Shen, Physical Review Letters 102, 136806 (2009).
  * Gonçalves  _et al._ (2018) M. Gonçalves, P. Ribeiro, and E. V. Castro, arXiv:1807.11247 [cond-mat, physics:quant-ph] (2018), arXiv: 1807.11247.
  * Meier  _et al._ (2018) E. J. Meier, F. A. An, A. Dauphin, M. Maffei, P. Massignan, T. L. Hughes, and B. Gadway, Science 362, 929 (2018).
  * Stützer  _et al._ (2018) S. Stützer, Y. Plotnik, Y. Lumer, P. Titum, N. H. Lindner, M. Segev, M. C. Rechtsman, and A. Szameit, Nature 560, 461 (2018).
  * Zhang  _et al._ (2019) Z.-Q. Zhang, B.-L. Wu, J. Song, and H. Jiang, Physical Review B 100, 184202 (2019), arXiv:1906.04064 .
  * Li  _et al._ (2020) C. A. Li, B. Fu, Z. A. Hu, J. Li, and S. Q. Shen, Physical Review Letters 125, 166801 (2020), 2008.00513 .
  * Yang  _et al._ (2021a) Y.-B. Yang, K. Li, L.-M. Duan, and Y. Xu, Physical Review B 103, 085408 (2021a).
  * Agarwala  _et al._ (2020) A. Agarwala, V. Juričić, and B. Roy, Physical Review Research 2, 012067 (2020).
  * Wang  _et al._ (2021) J.-H. Wang, Y.-B. Yang, N. Dai, and Y. Xu, Physical Review Letters 126, 206404 (2021).
  * Peng  _et al._ (2022) T. Peng, C.-B. Hua, R. Chen, Z.-R. Liu, H.-M. Huang, and B. Zhou, Physical Review B 106, 125310 (2022).
  * Lóio  _et al._ (2023) H. Lóio, M. Gonçalves, P. Ribeiro, and E. V. Castro,  (2023), arXiv:2305.19209 .
  * Zhang  _et al._ (2020) W. Zhang, D. Zou, Q. Pei, W. He, J. Bao, H. Sun, and X. Zhang, Physical Review Letters 126, 146802 (2020), arXiv:2008.00423 .
  * Edwards and Jones (1971) S. F. Edwards and R. C. Jones, Journal of Physics C: Solid State Physics 4, 2109 (1971).
  * Wan  _et al._ (1991) C. C. Wan, A. B. Harris, and J. Adler, Journal of Applied Physics 69, 5191 (1991).
  * Sandvik (2001) A. W. Sandvik, Physical Review Letters 86, 3209 (2001).
  * Sandvik (2002) A. W. Sandvik, Physical Review B 66, 024418 (2002).
  * Vojta and Schmalian (2005) T. Vojta and J. Schmalian, Physical Review Letters 95, 237206 (2005), arXiv:0508211 .
  * Abrahams  _et al._ (1979) E. Abrahams, P. W. Anderson, D. C. Licciardello, and T. V. Ramakrishnan, Physical Review Letters 42, 673 (1979).
  * Anderson (1958) P. W. Anderson, Physical Review 109, 1492 (1958).
  * Islam and Nakanishi (2008) M. Islam and H. Nakanishi, Physical review. E, Statistical, nonlinear, and soft matter physics 77, 061109 (2008).
  * Dillon and Nakanishi (2014) B. S. Dillon and H. Nakanishi, The European Physical Journal B 87 (2014), 10.1140/epjb/e2014-50397-4.
  * Holstein and Primakoff (1940) T. Holstein and H. Primakoff, Physical Review 58, 1098 (1940).
  * Shindou  _et al._ (2013) R. Shindou, R. Matsumoto, S. Murakami, and J.-i. Ohe, Physical Review B 87, 174427 (2013).
  * Zhang  _et al._ (2013a) L. Zhang, J. Ren, J.-S. Wang, and B. Li, Physical Review B 87, 144101 (2013a).
  * Mook  _et al._ (2014) A. Mook, J. Henk, and I. Mertig, Physical Review B 90, 024412 (2014).
  * Suding and Ziff (1999) P. N. Suding and R. M. Ziff, Physical Review E 60, 275 (1999).
  * Ziff (1992) R. M. Ziff, Physical Review Letters 69, 2670 (1992).
  * Fukui  _et al._ (2005) T. Fukui, Y. Hatsugai, and H. Suzuki, Journal of the Physical Society of Japan 74, 1674 (2005).
  * Zhang  _et al._ (2013b) Y. F. Zhang, Y. Y. Yang, Y. Ju, L. Sheng, D. N. Sheng, R. Shen, and D. Y. Xing, Chinese Physics B 22, 117312 (2013b), arXiv: 1212.6295.
  * Weiße  _et al._ (2006) A. Weiße, G. Wellein, A. Alvermann, and H. Fehske, Reviews of Modern Physics 78, 275 (2006).
  * Evers and Mirlin (2008) F. Evers and A. D. Mirlin, Reviews of Modern Physics 80, 1355 (2008).
  * Oganesyan and Huse (2007) V. Oganesyan and D. A. Huse, Physical Review B 75, 155111 (2007).
  * Atas  _et al._ (2013) Y. Y. Atas, E. Bogomolny, O. Giraud, and G. Roux, Physical Review Letters 110, 084101 (2013).
  * Wang  _et al._ (2014) C. Wang, Y. Avishai, Y. Meir, and X. R. Wang, Physical Review B 89, 45314 (2014).
  * Onoda and Nagaosa (2003) M. Onoda and N. Nagaosa, Physical Review Letters 90, 206601 (2003).
  * Onoda  _et al._ (2007) M. Onoda, Y. Avishai, and N. Nagaosa, Phys. Rev. Lett. 98 (2007), 10.1103/PhysRevLett.98.076802, 0605510 .
  * Castro  _et al._ (2016) E. V. Castro, R. de Gail, M. P. López-Sancho, and M. A. H. Vozmediano, Physical Review B 93, 245414 (2016).
  * Laughlin (1984) R. B. Laughlin, Physical Review Letters 52, 2304 (1984).
  * Prodan (2011) E. Prodan, Journal of Physics A: Mathematical and Theoretical 44, 113001 (2011), arXiv:1010.0595 .
  * Wang  _et al._ (2015) C. Wang, Y. Su, Y. Avishai, Y. Meir, and X. R. Wang, Physical Review Letters 114, 096803 (2015), arXiv:1411.4838 .
  * Qiao  _et al._ (2016) Z. Qiao, Y. Han, L. Zhang, K. Wang, X. Deng, H. Jiang, S. A. Yang, J. Wang, and Q. Niu, Phys. Rev. Lett. 117, 056802 (2016).
  * Su  _et al._ (2016) Y. Su, C. Wang, Y. Avishai, Y. Meir, and X. R. Wang, Scientific Reports 6, 33304 (2016).
  * Yang  _et al._ (2021b) H. Yang, J. Zeng, Y. Han, and Z. Qiao, Physical Review B 104, 115414 (2021b).
  * Castro Neto  _et al._ (2009) A. H. Castro Neto, F. Guinea, N. M. R. Peres, K. S. Novoselov, and A. K. Geim, Reviews of Modern Physics 81, 109 (2009), arXiv:0709.1163 .
  * Huang  _et al._ (2017) B. Huang, G. Clark, E. Navarro-Moratalla, D. R. Klein, R. Cheng, K. L. Seyler, D. Zhong, E. Schmidgall, M. A. McGuire, D. H. Cobden, W. Yao, D. Xiao, P. Jarillo-Herrero, and X. Xu, Nature 546, 270 (2017).
  * Davis and Narath (1964) H. L. Davis and A. Narath, Phys. Rev. 134, A433 (1964).
  * Narath and Davis (1965) A. Narath and H. L. Davis, Phys. Rev. 137, A163 (1965).
  * Samuelsen  _et al._ (1971) E. J. Samuelsen, R. Silberglitt, G. Shirane, and J. P. Remeika, Phys. Rev. B 3, 157 (1971).
  * Chen  _et al._ (2020) L. Chen, J.-H. Chung, T. Chen, C. Duan, A. Schneidewind, I. Radelytskyi, D. J. Voneshen, R. A. Ewings, M. B. Stone, A. I. Kolesnikov, B. Winn, S. Chi, R. A. Mole, D. H. Yu, B. Gao, and P. Dai, Phys. Rev. B 101, 134418 (2020).
  * Owerre (2017) S. A. Owerre, Journal of Physics Communications 1, 021002 (2017), arXiv:1705.04694 [cond-mat].
  * Koretsune  _et al._ (2015) T. Koretsune, N. Nagaosa, and R. Arita, Scientific Reports 5, 13302 (2015).
  * Legrand  _et al._ (2022) W. Legrand, Y. Sassi, F. Ajejas, S. Collin, L. Bocher, H. Jia, M. Hoffmann, B. Zimmermann, S. Blügel, N. Reyren, V. Cros, and A. Thiaville, Physical Review Materials 6, 024408 (2022).
  * Antão and Peres (2023) T. V. C. Antão and N. M. R. Peres, Physical Review B 107, 235410 (2023), arXiv:2303.03305 .
  * Zhuo  _et al._ (2021) F. Zhuo, H. Li, and A. Manchon, Physical Review B 104, 144422 (2021).
  * Zhuo  _et al._ (2022) F. Zhuo, H. Li, and A. Manchon, New Journal of Physics 24, 023033 (2022).


#  Quantum Scales of Galaxies from Ultralight Dark Matter

Jae-Weon Lee  Department of Electrical and Electronic Engineering, Jungwon
University, Jungwon university, 85 Munmu-ro, Goesan-eup, Goesan-gun,
Chungcheongbuk-do, 28024, Korea

###### Abstract

We propose that the ultralight dark matter (ULDM) model, in which dark matter
particles have a tiny mass of $m=O(10^{-22})eV$, has characteristic scales for
physical quantities of observed galaxies such as mass, size, acceleration,
mass flux, and angular momentum from quantum mechanics. The typical angular
momentum per dark matter particle is $\hbar$ and the typical physical
quantities are functions of specific angular momentum $\hbar/m$ and average
background density of the particles. If we use the Compton wavelength instead
for the length scale, we can obtain bounds for these physical quantities. For
example, there is an upper bound for acceleration of ULDM dominated objects,
$a_{c}={c^{3}m}/{\hbar}$. We suggest that the physical scales of galaxies
depend on the time of their formation and that these characteristic scales are
related to some mysteries of observed galaxies. Future observations from the
James Webb Space Telescope and NANOGrav can provide evidences for the presence
and evolution of these scales.

##  I Introduction

Recently, there is a growing interest in the Ultralight Dark Matter (ULDM)
model as a compelling alternative to Cold Dark Matter (CDM). In this model,
dark matter particles possess an exceptionally small mass, typically on the
order of $10^{-22}$ eV, and they exist in a Bose-Einstein condensate (BEC)
state. (For a review see Lee (2009); Suárez et al. (2014); Rindler-Daller and
Shapiro (2014); Chavanis (2011); Marsh (2016); Hui et al. (2016)). This model
is known by various other names, including Fuzzy DM, BEC DM, scalar field DM,
ultra-light axion, and wave-$\psi$ DM Baldeschi et al. (1983); Sin (1994a);
Lee and Koh (1996); Widrow and Kaiser (1993); Nucamendi et al. (2000); Arbey
et al. (2001); Goodman (2000); Hu et al. (2000); Peebles (2000); Mielke and
Schunck (2002); Mielke and Perez (2009); Sahni and Wang (2000); Alcubierre et
al. (2002a); Park et al. (2012); Sikivie and Yang (2009); Fuchs and Mielke
(2004); Matos et al. (2001); Silverman and Mallett (2001); Nucamendi et al.
(2001); Julien Lesgourgues and Salati (2002); Boehmer and Harko (2007); Eby et
al. (2015). Because of its extremely small mass, ULDM particles exhibit an
exceedingly high number density, resulting in a significant overlap of their
wave functions, and hence they are described by macroscopic wave functions.

Galaxies, the smallest celestial objects dominated by DM, exhibit a
distinctive behavior in the context of the ULDM model, where the quantum
pressure arising from the uncertainty principle acts as a counterforce,
preventing the collapse of dark matter halos. Therefore, the de Broglie wave
length $\lambda_{dB}=\hbar/mv$ of ULDM determines the typical size of
galaxies, where $v$ is the typical velocity of the halo dark matter. In this
sense one can say that quantum mechanical scale $\lambda_{dB}$ decides the
size of galaxies. It was conjectured that this scale can explain the size and
mass evolution of galaxies Lee (2009); van Dokkum et al. (2009).

Numerical simulations of galactic halos with ULDM Schive et al. (2014a)
revealed that a solitonic core exists in the halos, surrounded by granules.
Both the core and granules have a characteristic size of about $kpc$. The
characteristic length scale inherent to ULDM has been recognized for its
ability to address several of the small-scale issues associated with CDM,
including the core-cusp problem, the satellite galaxy problem and the missing
satellite problem Park et al. (2022); Salucci et al. (2003); Navarro et al.
(1996); de Blok et al. (2002); Tasitsiomi (2003). It is natural to investigate
the characteristic scales of other physical properties such as angular
momentum and acceleration in the ULDM model.

In this paper, we propose a typical scales for various physical quantities of
galaxies in this model and suggest that these scales are related to mysteries
of galaxies. We provide an easy way to calculate these scales from
parametrization of variables. In the section II the characteristic scales are
derived and briefly compared to the observed values. In the section III we
consider the time evolution of the scales. In the section IV we discuss the
results.

##  II Scales of Fuzzy dark matter

ULDM particles can be described by a macroscopic wave function $\psi$
satisfying the following Schr$\ddot{o}$dinger-Poisson equation (SPE);

$\displaystyle i\hbar\partial_{{t}}{\psi}$ $\displaystyle=$
$\displaystyle-\frac{\hbar^{2}}{2m}\nabla^{2}{\psi}+m{V}{\psi},$ (1)
$\displaystyle\nabla^{2}{V}$ $\displaystyle=$ $\displaystyle{4\pi G}\rho,$

where the the DM mass density $\rho=mN|\psi|^{2}$, and $V$ is the
gravitational potential. $N$ is the number of the DM particles in a system.
The wave function is normalized as $\int|\psi|^{2}d^{3}{\bf{x}}=1$. Here, The
SPE has a soliton solution which well approximates cores of galaxies and
represents typical structures in the ULDM models.

To derive the scales let us consider the following dimensionless variables
with hats by the parameterizations;

$\displaystyle t$ $\displaystyle\equiv
t_{c}\,\hat{t}=\frac{\hbar^{3}}{m^{3}}\frac{1}{\left(GM\right)^{2}}\hat{t},$
(2) $\displaystyle{\bf{x}}$
$\displaystyle\equiv{\bf{x}}_{c}\,\hat{{\bf{x}}}=\frac{\hbar^{2}}{m^{2}}\frac{1}{GM}\,\hat{{\bf{x}}},$
(3) $\displaystyle\psi$
$\displaystyle\equiv\psi_{c}\,\hat{\psi}=\frac{m^{3}}{\hbar^{3}}\left(GM\right)^{\frac{3}{2}}\hat{\psi},$
(4) $\displaystyle V$ $\displaystyle\equiv
V_{c}\,\hat{V}=\frac{m^{2}}{\hbar^{2}}\left({4\pi GM}\right)^{2}\hat{V},$ (5)

which lead to a dimensionless version of the SPE apt for a numerical
simulation;

$\displaystyle i\,\partial_{t}\hat{\psi}(\hat{{\bf{x}}},\hat{t})\,$
$\displaystyle=-\frac{1}{2}\nabla^{2}\hat{\psi}(\hat{{\bf{x}}},\hat{t})+\hat{V}(\hat{{\bf{x}}},\hat{t})\hat{\psi}(\hat{{\bf{x}}},\hat{t}),$
(6) $\displaystyle\nabla^{2}\hat{V}(\hat{{\bf{x}}},\hat{t})$ $\displaystyle=\
4\pi\,|\hat{\psi}|^{2}(\hat{{\bf{x}}},\hat{t}).$ (7)

The scale parameters, such as ${\bf{x}}_{c}$, represent the characteristic
scales of a system, because numerical solutions of these dimensionless
equations often have values of order unity.

For example, $t_{c}$ can represent the typical relaxation time scale of ULDM
soliton oscillations such as quasi-normal modes. Excited modes of ULDM halos
decay to the ground state over the characteristic time $t_{c}$ by
gravitational cooling which is a mechanism for relaxation by ejecting ULDM
waves gravitationalcooling . The time scale can be written as

$t_{c}=\frac{\hbar^{3}}{m^{3}}\frac{1}{\left(GM\right)^{2}}=3.7\times
10^{7}years\left(\frac{10^{-22}eV}{m}\right)^{3}\left(\frac{10^{8}M_{\odot}}{M}\right)^{2},$
(8)

where $M$ can be the mass of the core of a galaxy.

The origin of these scales can be roughly understood in the following way. The
self-gravitational force of a ULDM soliton (i.e., boson star) is balanced by
the repulsive force from the quantum pressure due to the uncertainty
principle. By equating the de Broglie wavelength of the ULDM particles to the
gravitational radius, we can obtain the equation

$\lambda_{dB}=\hbar/(mv)=\hbar/(m\sqrt{GM/\lambda_{dB}}),$ (9)

which gives

$\lambda_{dB}={\bf{x}}_{c}=\left(\frac{\hbar}{m}\right)^{2}\frac{1}{GM}=854.8~{}pc\left(\frac{10^{-22}eV}{m}\right)^{2}\frac{10^{8}M_{\odot}}{M}.$
(10)

Therefore, ${\bf{x}}_{c}$ is the gravitational Bohr radius and the de Broigle
wavelength of individual ULDM particles making the soliton. (The usual scaling
$\hat{r}=\frac{mc}{\hbar}r$ corresponds to the case when $v=c$, that is, when
we use the Compton wavelength instead of the de Broigle wavelength.) The
scales of ULDM systems are determined by two mass parameters $m$ and $M$.

Now, one may wonder how to determine the typical mass $M$ from cosmology. If
we know the average dark matter density $\bar{\rho}$ of the universe one can
argue that $M$ is of order of $\bar{\rho}{\bf{x}}_{c}^{3}$. Thus, from this
relation

$M\simeq\bar{\rho}\left(\frac{\hbar}{m}\right)^{6}\left(\frac{1}{GM}\right)^{3}$
(11)

and one can obtain the typical mass scale

$M\simeq
G^{-\frac{3}{4}}\left(\frac{\hbar}{m}\right)^{\frac{3}{2}}\bar{\rho}^{\frac{1}{4}},$
(12)

which decrease as the time flows in the universe. This mass scale turns out to
be the order of the quantum Jeans mass described in the next section. Using
this $M$, the length scale can be rewritten as a function of $m$ and
$\bar{\rho}$ in turn;

${\bf{x}}_{c}\simeq\sqrt{\frac{\hbar}{m}}(G\bar{\rho})^{-1/4},$ (13)

which implies that galaxies extend in size as the universe expands Lee (2009);
van Dokkum et al. (2009). Using Eq. (12) the time scale turns out to be

$t_{c}\simeq\frac{1}{\sqrt{G\bar{\rho}}},$ (14)

which is just the order of the Hubble time at the formation of soliton with
mass $M$. Thus, intrinsically this time scale has an cosmological origin.

From the typical scales one can easily obtain other physical scales. For
example, the typical acceleration scale is given by

$a_{c}=x_{c}/t_{c}^{2}=G^{3}m^{4}M^{3}/\hbar^{4}=1.9\times
10^{-11}~{}meter/s^{2}\left(\frac{m}{10^{-22}eV}\right)^{4}\left(\frac{M}{10^{8}M_{\odot}}\right)^{3}\simeq\sqrt{\frac{\hbar}{m}}(G\bar{\rho})^{3/4}.$
(15)

Interestingly, this scale is similar to the Modified Newtonian dynamics (MOND)
scale $a_{0}=1.2\times 10^{-10}~{}meter/s^{2}$. In Ref. Lee et al. (2019) it
is proposed that MOND is an effective phenomenon of ULDM and $a_{c}$ is
related to the observed radial acceleration relation and the baryonic Tully-
Fisher relation. Eq. (15) implies that the effective MOND acceleration scale
decreases over time, which can be verified by a redshift dependency of it in
the future observations.

The typical velocity is

$v_{c}\equiv{\bf{x}}_{c}/t_{c}=GMm/\hbar=22.4~{}km/s\left(\frac{M}{10^{8}M_{\odot}}\right)\left(\frac{m}{10^{-22}eV}\right)\simeq\sqrt{\frac{\hbar}{m}}(G\bar{\rho})^{1/4},$
(16)

which is similar to the typical velocity dispersion in a dwarf galaxy. $v_{c}$
leads to a typical angular momentum in turn

$L_{c}=M{\bf{x}}_{c}v_{c}={\hbar}\frac{M}{m}=1.1\times
10^{96}~{}\hbar~{}\left(\frac{M}{10^{8}M_{\odot}}\right)\left(\frac{10^{-22}eV}{m}\right)\simeq\frac{\left(\frac{\hbar}{m}\right)^{5/2}\bar{\rho}^{1/4}}{G^{3/4}}.$
(17)

This implies that angular momentum of galactic DM halos usually decreases over
time. Interestingly (or obviously), the typical angular momentum per dark
matter particle $L_{c}/N$ is just $\hbar$, because the typical number of the
particles $N$ is $M/m$. This result implies that each ULDM particle can have a
quantized angular momentum. Note that $v_{c},L_{c}$, and $M$ scale similarly
as $\bar{\rho}^{1/4}$, while $a_{c}$ decreases much faster.

On the other hand, the wave function scales as

$\psi_{c}=\frac{m^{3}}{\hbar^{3}}\left(GM\right)^{\frac{3}{2}}=4~{}\times
10^{-5}~{}pc^{-3/2}\left(\frac{m}{10^{-22}eV}\right)^{3}\left(\frac{M}{10^{8}M_{\odot}}\right)^{3/2}\simeq{\left(\frac{\hbar}{m}\right)^{-3/4}}{(G\bar{\rho})^{3/8}},$
(18)

which is not proportional to $\bar{\rho}^{1/2}$ as naively expected. Then, the
typical density $\rho_{c}\equiv M/{\bf{x}}_{c}^{3}$ is

$\rho_{c}=\frac{G^{3}m^{6}M^{4}}{{\hbar}^{6}}=0.16~{}M_{\odot}/pc^{3}~{}\left(\frac{m}{10^{-22}eV}\right)^{6}\left(\frac{M}{10^{8}M_{\odot}}\right)^{4}.$
(19)

Similarly, one can define the probability current

$J_{c}=\frac{\hbar}{m}Im(\psi_{c}\nabla\psi_{c}^{*})\simeq\frac{\hbar}{m}\frac{\psi_{c}^{2}}{{\bf{x}}_{c}}=\frac{G^{4}m^{7}M^{4}}{{\hbar}^{7}}\simeq\frac{Gm\bar{\rho}}{\hbar},$
(20)

which represents the flux of ULDM particles in a galaxy. Then, the mass flux
is

$MJ_{c}=3.66\times
10^{-6}M_{\odot}/pc^{2}/year\left(\frac{m}{10^{-22}eV}\right)^{7}\left(\frac{M}{10^{8}M_{\odot}}\right)^{5},$
(21)

which is a DM mass flow per $pc^{2}$ per year in a region of a galaxy.

The scale for the gravitational potential,

$V_{c}\ =\frac{m^{2}}{\hbar^{2}}\left({4\pi GM}\right)^{2}=8.8\times
10^{-7}c^{2}~{}\left(\frac{m}{10^{-22}eV}\right)^{2}\left(\frac{M}{10^{8}M_{\odot}}\right)^{2},$
(22)

becomes about $c^{2}$ when the system becomes relativistic. For example, the
soliton (a core of a galaxy) becomes relativistic when $M\simeq
10^{11}M_{\odot}$ for $m=10^{-22}eV$.

Another interesting observation is that all the typical values above, except
for $t_{c}$, are functions of a specific angular momentum $\hbar/m$ and
$\bar{\rho}$. The dependency on $\hbar/m$ means that the heavier the ULDM
particles are, the less quantum the halos are, as expected. The limit
$m\rightarrow\infty$ corresponds to the CDM case (${\bf{x}}_{c}\rightarrow
0$). The specific angular momentum is given by

$\frac{\hbar}{m}=0.019\times\left(\frac{10^{-22}eV}{m}\right)pc^{2}/year,$
(23)

and usually the typical scales in the ULDM model is a function of this
quantity Hui et al. (2016).

##  III Time evolution of the Scales

The time dependency of the physical properties of galaxies is an important
subject in astronomy. Except for ${\bf{x}}_{c}$ and $t_{c}$ all the quantities
above are decreasing functions of the time. Since the average density of the
universe decreases as the time goes, we expect we can find more massive,
compact, fast moving galaxies in the early universe than in the present
universe Lee and Lim (2010), which seems to be consistent with recent
observations Curtis-Lake et al. (2023). If future observations reveal the time
evolution of these quantities of galaxies, it will provide another support for
the ULDM model.

Let us find how the quantities evolve over time. It is well known that the SPE
has an interesting scaling property

$\\{t,{\bf{x}},\psi,\rho,V\\}\rightarrow\\{\lambda^{-2}t,\lambda^{-1}\\\ \\\
{\bf{x}},\lambda^{2}\psi,\lambda^{4}\rho,\lambda^{2}V\\},$ (24)

where $\lambda$ is a scaling parameter. This leads to a scaling law for other
quantities;

$\\{M,E,L\\}\rightarrow\\{\lambda M,\lambda^{3}E,\lambda L\\},$ (25)

where $M$ is the mass, $E$ is the energy, and $L$ is the angular momentum of a
dark matter halos. This scaling property help us to determine the time
evolution of the scales. During the matter dominated era $\bar{\rho}$ scales
as $(1+z)^{3}$, thus by setting $\lambda=(1+z)^{3/4}$ we can easily estimate
the time evolution of the galactic halos as follows,

$\displaystyle\\{t,{\bf{x}},\psi,\rho,V\\}$
$\displaystyle\rightarrow\\{(1+z)^{-3/2}t,(1+z)^{-3/4}{\bf{x}},(1+z)^{3/2}\psi,(1+z)^{3}\rho,(1+z)^{3/2}V\\},$
(26) $\displaystyle\\{M,E,L\\}$
$\displaystyle\rightarrow\\{(1+z)^{3/4}M,(1+z)^{9/4}E,(1+z)^{3/4}L\\}.$ (27)

A rescaling by $\lambda$ is equivalent to looking back to the corresponding
redshift $z$. Therefore, the ULDM model suggests that the characteristic
scales of physical quantities of galaxies depend on the time when they formed.
The evolution of these quantities can be verified by future observations by,
for example, space telescopes like James Webb space telescope. Using the above
relation the observed size evolution of very massive compact galaxies in the
early universe can be explained, if dark matter is ULDM Lee and Lim (2010).
Fig. 1 shows the time evolution of the various typical scales of galaxies as a
function of the redshift $z$.

The typical scales described above can be also justified by carefully
considering the quantum Jeans length which can be found by inserting the
Madelung relation Chavanis (2011); Suárez et al. (2014),
$\psi({\bf{x}},t)=\sqrt{\rho({\bf{x}},t)}e^{iS({\bf{x}},t)}$, into the SPE.
Then, one can get an equation for density perturbation
$\delta\rho\equiv\rho-\bar{\rho}$ Chavanis (2012); Suarez and Matos (2011),

$\frac{\partial^{2}\delta\rho}{\partial
t^{2}}+\frac{\hbar^{2}}{4m^{2}}\nabla^{2}(\nabla^{2}\delta\rho)-4\pi
G\bar{\rho}\delta\rho=0,$ (28)

which turns into an equation for density contrast
$\delta\equiv\delta\rho/\bar{\rho}=\sum_{k}\delta_{k}e^{ik\cdot{\bf{x}}}$ with
a wave vector $k$,

$\frac{d^{2}\delta_{k}}{dt^{2}}+\left[c^{2}_{q}k^{2}-4\pi
G\bar{\rho}\right]\delta_{k}=0,$ (29)

where $c_{q}=\hbar k/2m$ is a velocity for a mode with $k$. By equating two
terms in the parenthesis Hu et al. (2000) one can get the quantum Jeans length
scale at a redshift $z$,

$\lambda_{Q}(z)=\left(\frac{\pi^{3}\hbar^{2}}{Gm^{2}\bar{\rho}(z)}\right)^{1/4},$
(30)

which is of the same order of ${\bf{x}}_{c}$ in Eq. (13). The density
perturbation with a size smaller than $\lambda_{Q}$ cannot collapse to form a
galaxy. Therefore, ${\bf{x}}_{c}\simeq\lambda_{Q}$ represents the minimum
length scale of galactic halos formed at $z$ Lee (2009, 2016); Lee and Lim
(2010), which can resolve the small scale issues of the CDM Peebles (2000);
Sahni and Wang (2000). Ref. Lee and Koh, 1996 showed that repulsive self-
interaction of ULDM drastically increases the typical length scale, which is
not considered in this work because the mechanisms are very different. (See
Ref. chavanis,  for more details.) On the other hand, it is well known that
collision-less cold dark matter and self-interacting particle dark matter (not
ULDM) has a Jeans scale $v(\pi/G\bar{\rho}(z))^{1/2}$ which is different from
$\lambda_{Q}$ of ULDM. Here, $v$ is the velocity dispersion for cold dark
matter and the sound velocity for self-interacting dark matter, respectively.
Therefore, ULDM has a unique time-dependency of its typical size different
from those of other dark matter models.

It was shown Lee (2016) that the observed size evolution of the most massive
galaxies is consistent with our model, if we use the relation
$r_{*}(z)\propto\lambda_{Q}^{3/2}(z)\propto(1+z)^{-9/8}$ where $r_{*}$ is the
size of the visible part of a galaxy. Recently, a compact early galaxy
referred to as RX J2129-z95 was observed Williams et al. (2023) by James Webb
space telescope. It has stellar mass about $10^{7.63}M_{\odot}$ and its half-
light radius is just about $16.2~{}pc$ at $z=9.51$, which is much smaller than
the typical size ($O(10^{2})~{}pc$) of similar mass galaxies at present. This
observation is also roughly consistent with our model.

From $\lambda_{Q}$ the quantum Jeans mass can be derived as

$M_{J}(z)=\frac{4}{3}\pi^{\frac{13}{4}}\left(\frac{\hbar}{G^{\frac{1}{2}}m}\right)^{\frac{3}{2}}{\bar{\rho}}(z)^{\frac{1}{4}},$
(31)

which is of the same order of $M$ in Eq. (12). Because the structure formation
theory of the ULDM model naturally supports the scales ${\bf{x}}_{c}$ and $M$,
it also supports the derived scales $v_{c}$, $a_{c}$ and $L_{c}$.

Figure 1:  The evolution of the various physical quantities $Q$ of galaxies
versus redshift $z$, normalized to the current values. The dashed line
corresponds to $Q=\psi_{c}$ and $Q=V_{c}$, while the solid line represents
$Q=M$ and $Q=L_{c}$, and the dotted line represents $Q={\bf{x}}_{c}$.

Not every scale is time dependent in the ULDM model. To see this let us
consider another popular scaling relation Lee and Koh (1996) which is useful
for relativistic cases;

$\displaystyle t$ $\displaystyle\equiv
t_{c}\,\hat{t}=\frac{\hbar}{mc^{2}}\hat{t},$ (32) $\displaystyle{\bf{x}}$
$\displaystyle\equiv{\bf{x}}_{c}\,\hat{{\bf{x}}}=\frac{\hbar}{mc}\,\hat{{\bf{x}}},$
(33) $\displaystyle\psi$
$\displaystyle\equiv\psi_{c}\,\hat{\psi}=\frac{mc^{2}}{\hbar\sqrt{4\pi
G}}\hat{\psi},$ (34) $\displaystyle V$ $\displaystyle\equiv
V_{c}\,\hat{V}=c^{2}\hat{V},$ (35)

where ${\bf{x}}_{c}$ now corresponds to the Compton wavelength of the ULDM
particles, and we use a different normalization
$4\pi\int|\psi|^{2}d^{3}{\bf{x}}=M$.

Equating ${\bf{x}}_{c}$ with the gravitational radius again one can get a
typical mass scale

$M=\frac{\hbar c}{mG}=1.3\times
10^{12}M_{\odot}\left(\frac{10^{-22}eV}{m}\right),$ (36)

which is about the maximum mass a boson star (ULDM soliton) can have, beyond
which it will collapse to form a black hole Lee and Koh (1996); Khlopov
(1985). Therefore, this theory also predicts the presence of super massive
black holes with mass larger than $10^{12}M_{\odot}$ in the early universe.
Interestingly, this bound is similar to the observed mass bound of galaxies.

The typical density for this normalization is

$\rho_{c}=\frac{c^{4}m^{2}}{G{\hbar}^{2}}=5.1\times
10^{15}~{}M_{\odot}/pc^{3}~{}\left(\frac{m}{10^{-22}eV}\right)^{2}.$ (37)

With this parameterization, obviously, $v_{c}={\bf{x}}_{c}/t_{c}=c,$ and
$L_{c}=M{\bf{x}}_{c}v_{c}=N\hbar$. Note that all these scales are independent
of $\bar{\rho}$ and time independent contrary to the scales in Eq. (5).

The probability current is now

$J_{c}=\frac{c^{5}m^{2}}{4\pi G\hbar^{2}}=1.24\times
10^{14}M_{\odot}/pc^{2}/year\left(\frac{m}{10^{-22}eV}\right)^{3},$ (38)

which is actually a mass flux, due to the new normalization
$4\pi\int|\psi|^{2}d^{3}{\bf{x}}=M$.

The typical acceleration in this case

$a_{c}={c^{3}m}/{\hbar}=45.5~{}meter/s^{2}\left(\frac{m}{10^{-22}eV}\right),$
(39)

which is the maximum acceleration scale the soliton can have before the
soliton collapses to form a black hole and is different from the MOND scale.
One can see that the gravitational acceleration of ULDM halos cannot exceed
the limit without including other matter. This upper bound of acceleration is
one of the unique feature of ULDM, and any stable ULDM dominated structures
should satisfy. If we use Planck’s constant $h$ instead of $\hbar$, $a_{c}$
becomes $7.24~{}meter/s^{2}$, which is similar to the gravitational
acceleration of the earth.

##  IV Discussions

Based on simple arguments using the SPE, we showed that one can easily derive
the characteristic scale of various physical quantities for observed galaxies
from quantum mechanics, if dark matter is ULDM. The presence of the quantum
scales is a very unique feature of the ULDM model. For example, The angular
momentum per a DM particle is of order of $\hbar$, which means that the
collective wave functions describing the rotating DM halos should be one of
eigenstates of the angular momentum operators, because the ULDM particles are
bosons in a single coherent state. The maximum acceleration scale in the ULDM
model is another example which is absent in other ‘classical’ dark matter
model such as CDM or warm dark matter.

Future data from observations such as the James Webb Space Telescope and the
North American Nanohertz Observatory for Gravitational Waves (NANOGrav)
NANOGrav:2023hfp  will open new windows into the galaxy evolution, possibly
confirming these characteristic scales and their evolution. Some of these
scales have already been checked by numerical simulations, but more
sophisticated simulations with baryons and black holes are needed to compare
the theory to observations more precisely.

### acknowledgments

This work was supported by the Jungwon University Research Grant(2021-031).

## References

  * Lee (2009) J.-W. Lee, Journal of Korean Physical Society 54, 2622 (2009), eprint 0801.1442.
  * Suárez et al. (2014) A. Suárez, V. H. Robles, and T. Matos, Astrophysics and Space Science Proceedings 38, 107 (2014).
  * Rindler-Daller and Shapiro (2014) T. Rindler-Daller and P. R. Shapiro, Modern Physics Letters A 29, 1430002 (2014).
  * Chavanis (2011) P.-H. Chavanis, Phys. Rev. D 84, 043531 (2011), eprint 1103.2050.
  * Marsh (2016) D. J. E. Marsh, Phys. Rept. 643, 1 (2016).
  * Hui et al. (2016) L. Hui, J. P. Ostriker, S. Tremaine, and E. Witten (2016), eprint 1610.08297.
  * Baldeschi et al. (1983) M. R. Baldeschi, G. B. Gelmini, and R. Ruffini, Physics Letters B 122, 221 (1983).
  * Sin (1994a) S.-J. Sin, Phys. Rev. D50, 3650 (1994a).
  * Lee and Koh (1996) J.-W. Lee and I.-G. Koh, Phys. Rev. D53, 2236 (1996), eprint hep-ph/9507385.
  * Widrow and Kaiser (1993) L. M. Widrow and N. Kaiser, Astrophys. J. Lett. 416, L71 (1993).
  * Nucamendi et al. (2000) U. Nucamendi, M. Salgado, and D. Sudarsky, Phys. Rev. Lett. 84, 3037 (2000).
  * Arbey et al. (2001) A. Arbey, J. Lesgourgues, and P. Salati, Phys. Rev. D 64, 123528 (2001).
  * Goodman (2000) J. Goodman, New Astronomy Reviews 5, 103 (2000).
  * Hu et al. (2000) W. Hu, R. Barkana, and A. Gruzinov, Phys. Rev. Lett. 85, 1158 (2000).
  * Peebles (2000) P. Peebles, Astrophys. J. 534, L127 (2000).
  * Mielke and Schunck (2002) E. W. Mielke and F. E. Schunck, Phys. Rev. D 66, 023503 (2002).
  * Mielke and Perez (2009) E. W. Mielke and J. A. V. Perez, Physics Letters B 671, 174 (2009).
  * Sahni and Wang (2000) V. Sahni and L. Wang, Phys. Rev. D 62, 103517 (2000).
  * Alcubierre et al. (2002a) M. Alcubierre, et al, Class. Quant. Grav. 19, 5017 (2002a).
  * Park et al. (2012) C.-G. Park, J.-c. Hwang, and H. Noh, Phys. Rev. D 86, 083535 (2012).
  * Sikivie and Yang (2009) P. Sikivie and Q. Yang, Physical Review Letters 103, 111301 (2009).
  * Fuchs and Mielke (2004) B. Fuchs and E. W. Mielke, Mon. Not. Roy. Astron. Soc. 350, 707 (2004).
  * Matos et al. (2001) T. Matos, F. S. Guzman, L. A. Urena-Lopez, and D. Nunez, astro-ph/0102419 (2001).
  * Silverman and Mallett (2001) M. P. Silverman and R. L. Mallett, Classical and Quantum Gravity 18, L103 (2001).
  * Nucamendi et al. (2001) U. Nucamendi, M. Salgado, and D. Sudarsky, Phys. Rev. D 63, 125016 (2001).
  * Julien Lesgourgues and Salati (2002) A. A. Julien Lesgourgues and P. Salati, New Astronomy Reviews 46, 791 (2002).
  * Boehmer and Harko (2007) C. G. Boehmer and T. Harko, JCAP 0706, 025 (2007), eprint 0705.4158.
  * Eby et al. (2015) J. Eby, C. Kouvaris, N. G. Nielsen, and L. C. R. Wijewardhana, arXiv:1511.04474 (2015).
  * Lee (2009) J.-W. Lee, Phys. Lett. B681, 118 (2009).
  * van Dokkum et al. (2009) P. G. van Dokkum, M. Kriek, and M. Franx, Nature (London) 460, 717 (2009).
  * Schive et al. (2014a) H.-Y. Schive, T. Chiueh, and T. Broadhurst, Nature Phys. 10, 496 (2014a).
  * Park et al. (2022) Park, S., Bak, D., Lee, J.-W., et al. 2022, JCAP , 2022, 033.
  * Salucci et al. (2003) P. Salucci, F. Walter, and A. Borriello, Astron. Astrophys. 409, 53 (2003).
  * Navarro et al. (1996) J. F. Navarro, C. S. Frenk, and S. D. M. White, Astrophys. J. 462, 563 (1996).
  * de Blok et al. (2002) W. J. G. de Blok, A. Bosma, and S. S. McGaugh, astro-ph/0212102 (2002).
  * Tasitsiomi (2003) A. Tasitsiomi, International Journal of Modern Physics D 12, 1157 (2003).
  * (37) Seidel, E., and Suen, W. M. Phys. Rev. D 42(2), 384 (1990).
  * Lee et al. (2019) J.-W. Lee, H.-C. Kim, and J. Lee, Phys. Lett. B 795, 206 (2019), eprint 1901.00305.
  * Lee and Lim (2010) J.-W. Lee and S. Lim, JCAP 1001, 007 (2010), eprint 0812.1342.
  * Curtis-Lake et al. (2023) Curtis-Lake, E., Carniani, S., Cameron, A., et al., Nature Astronomy, 7, 622 (2023).
  * Chavanis (2012) P. H. Chavanis, Astronomy and Astrophysics 537, A127 (2012).
  * (42) P. H. Chavanis, Phys. Rev. D 107, 103503 (2023).
  * Suarez and Matos (2011) A. Suarez and T. Matos, Mon. Not. Roy. Astron. Soc. 416, 87 (2011).
  * Lee (2016) J.-W. Lee, Phys. Lett. B756, 166 (2016), eprint 1511.06611.
  * Williams et al. (2023) Williams H., et al., Science, 380, 416 (2023).
  * Khlopov (1985) M. Yu. Khlopov, B. A. Malomed, Ya. B. Zeldovich, MNRAS, 215, 575 (1985).
  * (47) G. Agazie et al. [NANOGrav], Astrophys. J. Lett. 952 no.2, L37 (2023).