# DFT+$\mu$: Density Functional Theory for Muon Site Determination

S. J. Blundell  Department of Physics, Clarendon Laboratory, Oxford
University, Parks Road, Oxford OX1 3PU, United Kingdom  T. Lancaster
Department of Physics, Centre for Materials Physics, Durham University,
Durham, DH1 3LE, United Kingdom<EMAIL_ADDRESS>and
<EMAIL_ADDRESS>
###### Abstract

The technique of muon spin rotation ($\mu$SR) has emerged in the last few
decades as one of the most powerful methods of obtaining local magnetic
information. To make the technique fully quantitative, it is necessary to have
an accurate estimate of where inside the crystal structure the muon implants.
This can be provided by density functional theory calculations using an
approach that is termed DFT+$\mu$, density functional theory with the
implanted muon included. This article reviews this approach, describes some
recent successes in particular $\mu$SR experiments, and suggests some avenues
for future exploration.

††preprint: AIP/123-QED

##  I Introduction

Very often in condensed matter physics it is necessary to have very detailed
information about the magnetic properties of materials measured at a local
level. Standard magnetic characterisation can be provided by measurements of
magnetic susceptibility $\chi$; this quantifies the magnetic response of a
sample, revealed by its magnetic moment $m$ induced by an applied field $H$,
and ideally extracted in the limit of $H\to 0$, so that $\chi=\lim_{H\to
0}M/H$, where $M$ is the magnetization $M=m/V$ and $V$ is the sample’s volume.
This has two obvious problems: (1) the limit of $H\to 0$ is hard to achieve
since the magnitude of the measurement signal is often proportional to $H$, so
measurements have to be performed in non-zero, and sometimes substantial,
applied field; (2) the measurement averages over the entire sample volume,
since $m=\int_{V}M\,{\rm d}V$, and so there is no way to distinguish between a
sample that is uniformly ordered, but with a small $M$, and one that is not
ordered, but has a small component which is ordered with a large $M$. In such
cases, what is needed is a highly-sensitive local magnetic probe, ideally one
operating at the atomic level and giving single spin detection.

Spin-polarized muons provide exactly this kind of sensitive local magnetic
probe.Cox (1987); de Réotier and Yaouanc (1997); Blundell (1999); Blundell
_et al._ (2022) In a muon-spin rotation ($\mu$SR) experiment, a beam of spin-
polarized muons is implanted in the sample to be studied. When each muon
decays, a positron and two neutrinos are emitted. A property of this decay is
that the positron is not emitted isotropically, but preferentially along the
direction of the muon spin at the moment of decay.Garwin, Lederman, and
Weinrich (1957) Detecting the direction along which the positron is emitted
allows one to infer which direction the muon spin was pointing at the moment
of its decay. Each muon lives for a different amount of time (according to the
radioactive decay law) and so each muon decay gives rise to a positron whose
detection at time $t$ at a particular angle of emission contributes to a data
point. Making many such positron detections, taken over many million muon
decays, allows one to build up a histogram $A(t)$, yielding the average muon
spin polarization $P(t)$ of the large ensemble of muons as a function of time,
allowing both static and dynamic local fields to be probed.

The technique has a number of key advantages. First of all, muon beams have
100% spin polarization (owing to parity violation in the weak interaction), in
contrast to the very weak thermal polarizations obtained in nuclei in nuclear
magnetic resonance (NMR) experiments. Secondly, muons have a larger
gyromagnetic ratio than any nucleus in NMR, resulting in high sensitivity. It
is very helpful that there is no restriction or necessity for specific nuclear
isotopes, with muons usable for measurements in any material (in contrast to
both NMR and neutron scattering). The ability to use an extended field range
for measurements, from zero field (ZF) up to 10 T, and a wide range of
temperatures (using dilution refrigerators, helium cryostats or furnaces),
means that $\mu$SR experiments are compatible with a very wide range of sample
environment. Finally, the technique allows the user to study the time
dependence of local magnetic fields and the technique gives access to
dynamical information with a very wide range of correlation times
($10^{-10}$–$10^{-5}$ s) that fit neatly between a.c. susceptibilityTopping
and Blundell (2019) on one hand and neutron scatteringBoothroyd (2019) on the
other (the longest times are limited by the muon lifetime and the shortest
times depend on the local field distribution Blundell  _et al._ (2022)).
Consequently, it has been widely used to study problems in magnetism,de
Réotier and Yaouanc (1997) superconductivity,Sonier, Brewer, and Kiefl (2000)
organic conductors,Blundell (2004) and many other systems (for an up-to-date
review of applications of $\mu$SR to a very wide range of topics, see Ref.
Blundell  _et al._ , 2022), and has become a frequently deployed method
alongside other experimental techniques, such as magnetic neutron
diffraction.Boothroyd (2019)

Despite many successes of the technique, there is a perceived drawback of the
technique that arises from the lack of knowledge of the muon’s stopping site
in materials, raising the question of which magnetic fields are actually being
probed. Consider the dipolar field $\bm{B}_{\rm dip}$ measured by the muon,
resulting from the ordered magnetic moments in the sample. The
$\alpha$-component of the dipolar ($\alpha=x$, $y$ or $z$) is given by

$B_{\rm dip}^{\alpha}=\sum_{j}\sum_{\beta}D^{\alpha\beta}_{j}m^{\beta}_{j},$
(1)

where the first sum is taken over all the magnetic moments $\bm{m}_{j}$ within
the Lorentz sphere (i.e. out to a large enough radius from the muon) and the
dipolar tensor for the $j$th moment at position $\bm{r}_{j}$ is given by

$D^{\alpha\beta}_{j}={\mu_{0}\over 4\pi
R_{j}^{3}}\left({3R_{j}^{\alpha}R_{j}^{\beta}\over
R_{j}^{2}}-\delta^{\alpha\beta}\right),$ (2)

where $\bm{R}_{j}=\bm{r}_{j}-\bm{r}_{\mu}$ and $\bm{r}_{\mu}$ is the muon
position. These expressions are often used to model ordered magnetic
arrangements and this can be a very good test of the validity of a muon site,
but notice that the expressions depend on the coordinates $\bm{R}_{j}$ which
are measured with respect to the (unknown) muon position $\bm{r}_{\mu}$. There
are other contributions to the field at the muon site in addition to
$\bm{B}_{\rm dip}$ (such as the hyperfine field and, in the case of
ferromagnets, the Lorentz and demagnetization fieldsBlundell  _et al._
(2022)), but $\bm{B}_{\rm dip}$ is the most significant and illustrates the
point that the field at the muon site is a function of the unknown quantity
$\bm{r}_{\mu}$.

Even more serious is the unknown effect that the positively charged muon has
on its local environment, raising the question of whether intrinsic behaviour
is being measured, or whether instead the experimental results are caused by
the presence of the implanted muon itself. One can suspect that, if the final
muon stopping state involves relatively little contact hyperfine coupling, the
stopped muon will resemble something like a bare particle playing the role of
an interstitial defect, albeit ‘dressed’ by interactions with the electronic
system. This might involve, for example, the muon acquiring a screening cloud
of electronic charge in a metal, or acquiring a strain field in an insulator
as it deforms atoms in its vicinity. Knowledge of the muon’s stopping site in
such cases would allow us insight into the local fields and interactions at a
known position. If the muon instead forms a bound paramagnetic state such as
muonium, or induces a large spin density at its position, knowledge of the
muon site would therefore allow us to identify the electronic state that gives
rise to the measured hyperfine coupling constants, and hence direct insight
into the local electronic structure. To address this drawback of $\mu$SR, a
density functional theory (DFT) ab-initio approach has been developed to
determine the muon site.Möller  _et al._ (2013); Bernardini  _et al._ (2013);
Möller  _et al._ (2013); Bonfà and De Renzi (2016) We call this approach
DFT$+\mu$, density functional theory with an implanted muon, and this review
will describe the method and demonstrate some of its recent successes.

##  II The site problem

The location of the muon site is determined by all the electrostatic
interactions in the system, including the distortions produced in the material
by the presence of the implanted muon, and interactions between all the
particles in the system. This is a complicated problem, but is pretty well
described by non-relativistic quantum mechanics. The system comprising the
stopped muon and its host material can therefore be described by the
Schrödinger equation, and so a computation of the wavefunction of the system
should, in principle, be possible. The system is however rather large,
comprising the muon and the large number of electrons and nuclei in the
sample. This is a many body problem and therefore extremely complex. Even if
we restrict ourselves to focussing only on the electronic part of the wave
function (which makes sense, as it is the behaviour of the more mobile
electrons that ultimately determines the muon site) the task is
formidable.Kohn (1999) In Section III we will describe how this problem can be
treated with electronic structure calculations, but in the remainder of this
section we will describe some methods that have been used previously.

###  II.1 The intuitive approach

Even without the use of sophisticated electronic structure calculations, it is
possible to think about where a muon is likely to sit using chemical
intuition. Most $\mu$SR experiments use positive muons ($\mu^{+}$) and since
the muon behaves like a light isotope of hydrogen, the $\mu^{+}$ is likely to
implant at a site that would be favoured by the H+ cation. Therefore, one can
immediately deduce some ad hoc rules of thumb for guessing muon sites. For a
start, one can conclude that positive muons favour sites close to negatively
charged anions. A good example is provided by the large number of magnetic
oxides that have been widely studied in condensed matter physics. In these it
is usually found that the muon sits around 1 Å from an O2- ion (a result which
was concluded during studies of the cuprate superconductorsBrewer  _et al._
(1990) which all consist of a collection of cations and oxygen anions). The
resulting state behaves like a hydroxyl OH group, which has a bond length
close to 1 Å.Demaison, Herman, and Lievin (2007) In fluorides, it is known
that muons form F$-\mu-$F bondsBrewer  _et al._ (1986) in which the muon sits
between the two fluoride ions, forming an analogue of the bifluoride ion
(HF${}_{2}^{-}$). We will return to this important example later. In elemental
metals, the muon would be expected to occupy a highly symmetrical site, such
as an interstitial site inside the face-centred cubic lattice of copper for
example.Luke  _et al._ (1991) These empirical rules of thumb are often useful,
and lead to the identification of highly plausible muon sites, albeit in a
limited range of materials.

###  II.2 Bayesian methods

One fruitful approach to solving (or at least avoiding) the muon site problem
is to deny that knowledge of the muon site is necessary. If the goal is to
measure a magnetic moment of a magnetic ion in a structure, is the muon site
really needed and can one include our quantified ignorance of it in the
calculation? This approach leads to a Bayesian methodBlundell  _et al._ (2012)
and is based on the observation that evaluating the dipolar field, given the
muon site and the magnetic moment, is an easy problem (just use eqns 1 and 2),
but the reverse problem is hard. The inversion between the two calculations is
accomplished using Bayes’ theorem which states that $P({\sf A}|{\sf B})=P({\sf
A})P({\sf B}|{\sf A})/P({\sf B})$. Here $P({\sf A})$ is called the prior
probability, since it is the probability of ${\sf A}$ occurring without any
knowledge as to the outcome of ${\sf B}$. The quantity which you derive is
$P({\sf A}|{\sf B})$, the posterior probability. For the muon problem, we
write Bayes’ theorem as

$P(\mu|\nu)={P(\mu)P(\nu|\mu)\over\int
P(\nu|\mu^{\prime})P(\mu^{\prime})\,{\rm d}\mu^{\prime}},$ (3)

where $\nu$ is a muon precession frequency and $\mu$ is the magnetic moment of
an ion in the crystal. Equation 3 gives the probability of a magnetic moment,
given the observed muon precession frequency (the experimental problem, which
is hard), in terms of the probability distribution of the muon precession
frequency, given the size of the magnetic moment (the computational problem,
which is easy). These are probabilities, since we do not know the muon site a
priori and so the muon site is described by some probability distribution
within the unit cell (and for the case of total ignorance of the muon site,
the probability distribution can be uniform throughout the unit cell). Thus we
can compute the distribution of dipolar fields within the unit cell under some
assumed distribution of muon sites (reflecting the level of our ignorance),
which we can express as a probability density function (pdf) $f(\nu/\mu)$,
evaluated as a function of precession frequency $\nu$ divided by magnetic
moment $\mu$ (since the precession frequency scales with the magnetic moment).
This function $f(\nu/\mu)$ allows us to evaluate $P(\nu|\mu)$ of eqn (3). Thus
$P(\nu|\mu)={1\over\mu}f(\nu/\mu)$, and since $f(\nu/\mu)$ is normalized so
that $\int f(\nu/\mu)\,{\rm d}(\nu/\mu)=1$, the factor of ${1\over\mu}$ is
needed so that $\int P(\nu|\mu)\,{\rm d}\nu=1$

Since $\nu$ is obtained from a real experiment, what we would like to know is
$g(\mu|\nu)$, the pdf of $\mu$ given the observed $\nu$. This can be obtained
from our calculated $f(\nu/\mu)$ using Bayes’ theorem in the form of Eq. (3),
which yields

$g(\mu|\nu)={{1\over\mu}f(\nu/\mu)\over\int_{0}^{\mu_{\rm
max}}{1\over\mu^{\prime}}f(\nu/\mu^{\prime})\,{\rm d}\mu^{\prime}},$ (4)

where we have assumed a prior probability [$P(\mu)$] for the magnetic moment
that is uniform between zero and $\mu_{\rm max}$, and so $P(\mu)$ is replaced
by the uniform probability density $1/\mu_{\rm max}$ [which cancels on the top
and bottom of Eq. (4)]. We choose $\mu_{\rm max}$ to take a large value,
although it is foundBlundell  _et al._ (2012) that results are insensitive to
the precise value of $\mu_{\rm max}$. When multiple frequencies $\nu_{i}$ are
present in the spectra, it is necessary to multiply their probabilities of
observation in order to obtain the chance of their simultaneous observation,
so we evaluate
$g(\mu|\\{\nu_{i}\\})\propto\prod_{i}\int_{\nu_{i}-\Delta\nu_{i}}^{\nu_{i}+\Delta\nu_{i}}f(\nu_{i}/\mu)\,{\rm
d}\nu_{i}$, where $\Delta\nu_{i}$ is the error on the fitted frequency. These
results have been applied successfully in a number of studies.Blundell  _et
al._ (2012); Steele  _et al._ (2011); Disseler (2014); Prando  _et al._ (2020)

###  II.3 Knight shift measurements

An experimental method that has been used to identify or verify a muon site
involves applying a transverse magnetic field $\bm{B}_{\rm ext}$ to the sample
and measuring the Knight shift. This occurs because the field at the muon site
$\bm{B}_{\mu}$ may not be exactly the same as the applied field $\bm{B}_{\rm
ext}$, due to contributions from the dipolar field, the hyperfine contact
interaction, the Lorentz field, and the demagnetization field. In a transverse
field experiment, one is often interested in the difference between the two,
$\bm{B}_{\mu}-\bm{B}_{\rm ext}$, but measured along $\bm{B}_{\rm ext}$. This
is because the component of $\bm{B}_{\mu}-\bm{B}_{\rm ext}$ measured
perpendicular to $\bm{B}_{\rm ext}$ makes very little difference to the
precession frequency. The Knight shift $K$ is then defined by

$K={(\bm{B}_{\mu}-\bm{B}_{\rm ext})\cdot\bm{B}_{\rm ext}\over B_{\rm
ext}^{2}},$ (5)

and hence $B_{\mu}\approx(1+K)B_{\rm ext}$. The dipolar field from ordered
moments contributes to $K$, but even in the paramagnetic state the moments can
become partially polarized in an applied field and take the value

$\bm{m}_{j}={\underline{\bm{\chi}}\bm{B}_{\rm ext}{\cal V}_{\rm
c}\over\mu_{0}},$ (6)

where ${\cal V}_{\rm c}$ is the volume per magnetic ion and
$\underline{\bm{\chi}}$ is the magnetic susceptibility tensor. Hence

$\bm{B}_{\rm dip}=\underline{\bm{{\cal D}}}\underline{\bm{\chi}}\bm{B}_{\rm
ext},$ (7)

where $\underline{\bm{{\cal D}}}={{\cal V}_{\rm
c}\over\mu_{0}}\sum_{j}\underline{\bm{D}}_{j}$ is the total dipolar tensor.
Similarly, the contact hyperfine interaction $\bm{B}_{\rm hf}$ is given by

$\bm{B}_{\rm hf}=\underline{\bm{{\cal A}}}\underline{\bm{\chi}}\bm{B}_{\rm
ext},$ (8)

where $\underline{\bm{{\cal A}}}$ is the analogous hyperfine tensor. The
contact coupling is usually independent of the field direction and so this can
be written as a scalar in the majority of cases. This provides us with all we
need to model the effect of rotating a single crystal in a constant magnetic
field and predicting the Knight shift for a particular muon site, the
anisotropy of the dipolar coupling resulting in angle dependence, while the
contact interaction gives an angle-independent contribution.

Another experimental test of the validity of a candidate muon site can be
obtained using the contribution of nearby nuclear dipoles to the decay rate of
a transverse-field precession measurement. The damping of the precession
signal arises from contributions from nearby nuclear dipoles, which are not
ordered, but sometimes add and sometimes subtract from the applied field. In a
transverse field measurement, the broadening is given at short times by a
Gaussian relaxation function $\exp(-\sigma^{2}t^{2}/2)$, where the parameter
$\sigma^{2}$ (also known as the second moment $M_{2}$ in NMR) is given byVan
Vleck (1948); Abragam (1961); Slichter (1990); Blundell  _et al._ (2022)

$\sigma^{2}=\frac{1}{3}\left({\mu_{0}\over
4\pi}\right)^{2}\hbar^{2}\gamma_{\mu}^{2}\sum_{i}\gamma_{i}^{2}I(I+1){(1-3\cos^{2}\theta_{i})^{2}\over
r_{i}^{6}}.$ (9)

This depends on the spin $I$ of the nearby nuclei, and their position $r_{i}$
and gyromagnetic ratio $\gamma_{i}$. This depolarization rate is angle-
dependent ($\theta_{i}$ is the angle between the applied field and the vector
between the muon and the nucleus) and so this broadening can be measured
experimentally. In the case of a zero-field measurement the relaxation takes
on the Kubo-Toyabe form with a value of $\Delta^{2}$ given by

$\Delta^{2}=\frac{1}{3}\left({\mu_{0}\over
4\pi}\right)^{2}\hbar^{2}\gamma_{\mu}^{2}\sum_{i}\gamma_{i}^{2}I(I+1){(5-3\cos^{2}\theta_{i})\over
r_{i}^{6}},$ (10)

which again can be checked by experiment. Moreover, in some fortunate cases
(particularly fluorides) the zero-field signal from the nuclei contains much
more structure than a simple Kubo-Toyabe relaxation (which is derived from an
assumption of a Gaussian-distributed random distribution of local field
components). In those cases, there is much more information to go on to deduce
the muon site and understand the local environment.

##  III DFT+$\mu$

The methods described in the previous section can all be useful in guessing
and experimentally checking the muon site. However, electronic structure
calculations are now proving to be extremely reliable as an ab initio
technique for deducing muon sites, as well as for calculating the resulting
distortion induced by the presence of the muon. In this section, we first
review briefly the ideas behind density functional theory (DFT) (for an
excellent introductory review, see the book by GiustinoGiustino (2014)) and
then describe how DFT can be implemented for muon site calculations.

###  III.1 Density functional theory

The many body wave function is a complex object to evaluate and store (making
demands on memory size that grow exponentially with system size, resulting in
one of the pioneers of the subject even doubting whether the many body wave
function is a legitimate scientific concept for systems of more than a handful
of atoms!Kohn (1999)). The approach taken in DFTGiustino (2014); Kohn (1999)
is to dispense with the many body wave function and instead deal with the
functional $E[n(\bm{r})]$ whose output is the ground state energy of the
$N$-electron system and whose input is the electron density function
$n(\bm{r})$. The electron density is not a function of all of the coordinates
of all of the electrons; it is simply a function of the three position
coordinates $\bm{r}=(x,y,z)$; thus we treat the electron density as a fluid,
and forget that the electrons are really individual quantum objects obeying
antisymmetry conditions. The functional is therefore a quantity averaged over
$3N-3$ of the degrees of freedom of the $N$-particle wavefunction. It seems
surprising that we can make such a severe simplification and still compute the
energy of a system, but that this is possible is the result of the two
Hohenberg-Kohn (HK) theorems.Hohenberg and Kohn (1964)

The first HK theorem establishes a one-to-one correspondence between an
external potential $V(\bm{r})$ and the ground state electronic density
$n(\bm{r})$. It can be proved that the ground state density $n(\bm{r})$
uniquely determines the potential, and hence all of the properties of the
system. So, for example, the ground state energy $E$, can be written as a
unique functional of the electron density $n(\bm{r})$ via an expression of the
form

$\displaystyle E[n(\bm{r})]=$ $\displaystyle
T[n(\bm{r})]+E_{\mathrm{ee}}[n(\bm{r})]+E_{\mathrm{en}}[n(\bm{r})]$
$\displaystyle=$ $\displaystyle
F[n(\bm{r})]+\int\mathrm{d}^{3}r\,n(\bm{r})V_{\mathrm{en}}(\bm{r}),$ (11)

where, in the final line of eqn 11, we have separated out the electron-nuclear
part of the energy $E_{\mathrm{en}}$ and have written the rest as
$F[n(\bm{r})]=T[n(\bm{r})]+E_{\mathrm{ee}}[n(\bm{r})]$, which is a functional
encoding all of the purely electronic contributions to the energy. If we have
access to $E[n(\bm{r})]$ then we might be able to solve the problem of finding
the electronic configuration corresponding to it by minimizing $E$ with
respect to the electron density $n$. That is, we continually adjust the
function $n(\bm{r})$ until we find a form that gives us the lowest-energy
electron density $E$. Is this a safe strategy? The second HK theorem tells us
that it is.

The second HK theorem states that $E[n(\bm{r})]$ gives the lowest energy if,
and only if, the input density $n(\bm{r})$ is the true ground state electron
density. If we do not use the true density, then we obtain an upper bound on
the ground state energy. Thus the idea is to minimize the functional
$E[n(\bm{r})]$ with trial functions $n(\bm{r})$ and if we find the actual
minimum then we will have the true electron density. There is, however, a
problem, in that we do not actually know the exact form of the electron
functional $F[n(\bm{r})]$. So, how do we find it?

A major simplification of the calculation of the electron density is made if
we use the Kohn-Sham formulation,Kohn and Sham (1965) the idea of which is to
reformulate the problem by replacing the interacting system of many electrons
by a non-interacting system of many electrons constrained to have the same
electron density. Thus, the idea is to approximate the interacting kinetic
energy $T[n(\bm{r})]$ of the $N$ interacting particles by the kinetic energy
$T_{\mathrm{s}}[n(\bm{r})]$ of $N$ non-interacting particles with the same
density $n(\bm{r})$. We also extract from $F[n(\bm{r})]$ the Hartree component
$U[n(\bm{r})]$. The Hartree component is simply the classical Coulomb
interaction between regions of charge density, and can be easily written in
terms of $n(\bm{r})$. We then write

$F[n(\bm{r})]=T_{\mathrm{s}}[n(\bm{r})]+U[n(\bm{r})]+E_{\mathrm{xc}}.$ (12)

We have therefore bundled up our remaining ignorance into the functional
$E_{\mathrm{xc}}$, known as the exchange-correlation energy. To summarize, the
HK theorems tell us we can find the ground state energy by minimizing a
functional
$E[n(\bm{r})]=F[n(\bm{r})]+\int\mathrm{d}^{3}r\,V_{\mathrm{en}}(\bm{r})n(\bm{r})$
with respect to the density $n(\bm{r})$. The KS formulation, with its
replacement of the interacting system with an equivalent non-interacting one,
allows us to say that this minimization will be equivalent to a simpler
procedure where we minimize the total energy of a non-interacting system
subject to an effective potential $V_{\mathrm{s}}$, whose form is

$V_{\mathrm{s}}(\bm{r})=V_{\mathrm{en}}(\bm{r})+V_{\mathrm{H}}(\bm{r})+V_{\mathrm{xc}}(\bm{r}).$
(13)

The first term in the effective potential represents the electron-nuclear
potential, the second the Hartree interaction, and the third the exchange-
correlation potential. The end result of the KS approach is that we can solve
the ground state density for the interacting system by first solving the
Schrödinger equation for a single particle in an effective potential, using
the set of KS equations defined as:

$\left[\frac{\hat{p}^{2}}{2m}+V_{\mathrm{s}}(\bm{r})\right]\phi_{i}(\bm{r})=\varepsilon_{i}\phi_{i}(\bm{r}).$
(14)

The solutions of these KS equations give a set of KS wavefunctions and energy
levels labelled by $i$. These states are filled up by the set of available
electrons and the density can be constructed from these non-interacting
wavefunctions as $n(\bm{r})=\sum|\phi_{i}(\bm{r})|^{2}$, where the sum runs
over the occupied states.

Although we now have a set of single-particle problems to solve, we don’t
actually know the potential. Our ignorance has two aspects: the first is the
functional form of $V_{\mathrm{xc}}$, the second is that the Hartree potential
$V_{\mathrm{H}}$ and the exchange correlation potential $V_{\mathrm{xc}}$ both
depend on $n(\bm{r})$, the solution to the problem, which in turn is
constructed from $\phi_{i}(\bm{r})$. The first problem can be addressed by
approximating $V_{xc}$, using for example the local density approximation
(LDA)Kohn and Sham (1965) or some more sophisticated approach.Becke (1993);
Perdew, Burke, and Ernzerhof (1996); Burke (2012) We then solve the equations
iteratively: (i) an initial guess is first made of $n(\bm{r})$ and the
potential $V_{\mathrm{s}}$ is computed; (ii) the single-particle wavefunctions
are computed; (iii) the set of occupied $\phi_{i}$ is used to compute
$n(\bm{r})$, and then we start again by returning to step (i) and
recalculating the potential. This procedure is repeated until some suitable
convergence criterion is reached.

Central to the use of DFT is the Born-Oppenheimer (BO) approximation. This
approximation makes use of the fact that atomic nuclei are much heavier than
electrons, and therefore that the energy scale for the electronic part of the
wavefunction is significantly larger than that of the nuclear part of the
wavefunction, allowing these to be separated. The approximation therefore
allows us to treat the atomic positions as parameters of an electron-only
Hamiltonian that obeys a Schrödinger equation. The Born-Oppenheimer
approximation basically clamps the nuclei at given fixed positions and the
appropriate electron wavefunction is the one corresponding to this clamped
nuclear configuration. However, nothing tells us whether the configuration of
nuclei we specify is an equilibrium one. Since one of our main purposes in
using DFT is finding the details of the nuclear structure that accommodates
the muon, we must repeat our calculation of the electronic density for
different nuclear structures until a global minimum in the ground state energy
has been found. To do this, we allow the nuclei to move under the effect of
forces, until these forces become small enough that we judge the system is at
equilibrium. This process is called a geometry optimization (or a relaxation
of the structure).

The BO approximation will be important in applying DFT techniques to a crystal
containing a muon impurity. The muon is treated as a nucleus (a reduced-mass
hydrogen) which, like all nuclei within BO, is clamped into position. The fact
that the muon is 1/9 the mass of a proton, and therefore somewhere between a
nucleus and an electron, ultimately limits the applicability of the technique,
as we discuss in Section IV.6.

###  III.2 Adding the muon

The DFT+$\mu$ technique is based around using DFT to optimize the geometry of
a material along with a muon impurity, using methods developed for first-
principle calculations of defect states.Freysoldt  _et al._ (2014) In practice
we can start by randomly assigning an initial muon site in a unit cell of the
target material and relaxing the structure by allowing both the muon and the
atoms in the crystal to adjust their positions via repeated iteration until
the forces on them are reduced below the threshold. We then evaluate the total
energy of the final relaxed configuration. For crystalline materials, we
generally use DFT codes that assume periodic boundary conditions. Since the
muons in a $\mu$SR experiment are implanted in the ultra-dilute limit, to
ensure the boundary conditions never result in muon-muon interactions, we
usually specify the structure via a supercell, containing a number of unit
cells and a single muon. The charge state of the muon is determined by the
overall charge of the defect [$+1$ for diamagnetic and neutral (zero) for
paramagnetic states] which is fixed initially; note that for the diamagnetic
case a uniform, smeared-out Hartree-like compensating charge is added to the
supercell ensure overall neutrality,Leslie and Gillan (1985); Makov and Payne
(1995); Van de Walle (1994) and the final charged state of the impurity is to
be determined by the final relaxed configuration.

The geometry optimization is repeated for several other randomly-chosen
initial muon sites, with the lowest total energy of the relaxed structure
yielding the most likely muon site, as well as providing an estimate for the
structural distortion the muon introduces. Often this procedure will lead to a
range of candidate muon sites being produced. It is then necessary to form
these into clusters of roughly-equivalent sites, which is most efficiently
done by making use of the underlying symmetry of the unperturbed material’s
structure, along with some assumptions about the nature of the distortion and
the energy ranges of sites that are likely to be realised. The method has
proved extremely successful, providing results which are in quantitative
agreement for cases in which the muon site can be independently verified (for
example, ionic fluorides provide a valuable test-bed for evaluating the method
because of the quantum coherent oscillations produced in F–$\mu$–F
states,Brewer  _et al._ (1986) see Section IV.1).

Since this procedure is rather well defined for most systems, it is possible
to create simple software solutions that allow it to be carried out
automatically. The first calculationsMöller  _et al._ (2013); Bernardini  _et
al._ (2013) utilised bespoke programs calling the open-source electronic
structure code QUANTUM ESPRESSO.Giannozzi  _et al._ (2009) Recently some
examples have started to become available for general use. MuFinderHuddart
_et al._ (2022) is a tool that enables users to carry out muon-site
calculations via a simple graphical user interface (GUI). The procedure for
calculating muon sites, by generating initial muon positions, relaxing the
structures, and then clustering and analysing the resulting candidate sites,
can be done entirely within the GUI. The software was originally configured to
make use of the plane wave electronic structure code Castep.Clark  _et al._
(2005) Once candidate muon sites have been determined, the local magnetic
field at the muon site can then also be computed within the program, making
use of the Magnetic structure and mUon Embedding Site Refinement (MuESR)
software,Bonfà, Onuorah, and Renzi (2018) allowing the connection between the
muon sites obtained and experiment to be made. Alternative software tools are
also available, including one based on a combination of ab initio random
structure searching (AIRSS) and machine learning. Liborio, Sturniolo, and
Jochym (2018) Here, AIRSS is a general scheme for randomly generating possible
structures a system can adopt and then introduce biases based on chemical,
experimental and/or symmetry grounds.Pickard and Needs (2011)

The method described so far can be quite computationally costly, owing to the
need for many iterations of the optimization procedure to relax the structure.
In view of this, an alternative scheme is the unperturbed electrostatic
potential (UEP) method, in which the charge density of the host material
obtained from a DFT calculation is used, unperturbed by the presence of the
muon. This provides an estimate the Coulomb force acting on each of the
$\mu^{+}$ in the initial muon-containing structures. This potential is a
direct output of the DFT calculation, being the opposite of the Coulomb
potential experienced by electrons, and its minimum represents the candidate
site where the muon feels zero force. Although approximate, the UEP method is
extremely fast and works reasonably well in materials in which atomic
displacements don’t play a major role in stabilising a diamagnetic muon
site.Sturniolo and Liborio (2020)

To further reduce computational cost, semi-empirical methods can be employed,
which rely on simplifying approximations and parameterizations to allow more
efficient computations. Hartree-Fock (HF) theory Lancaster and Blundell (2014)
was much-used before DFT became established. For speeding up HF calculations,
a useful series of semi-empirical methods were developed in which a set of
empirical parameters is derived for each chemical element treated by the
method. This is done by optimizing against a training set of molecular data,
which includes properties such as geometry, ionization energy, and electric
dipole moment, and is generally limited to light atoms and molecular
materials.Blundell  _et al._ (2022) In view of the success of the efficient HF
semi-empirical approach, attempts were made to apply similar semi-empirical
ideas to the DFT calculation framework. This led to the method known as
density functional tight binding (DFTB), which focuses on parametrizing the
interactions between pairs of atoms.Sturniolo, Liborio, and Jackson (2019) The
DFTB method is relatively fast, since each stage in the geometry relaxation
only requires a single diagonalization of the energy matrix, rather than the
iterative self-consistent-field loop of a usual DFT calculation. Although the
semi-empirical methods are useful for rapidly obtaining geometries for a
subset of problems, such as muoniated molecular radicals, they are not able to
produce accurate spin structures. Thus, in order to produce a reliable spin
distribution and corresponding hyperfine parameters, it is necessary to follow
on with a single-point DFT calculation after the semi-empirical geometry
optimization.Blundell  _et al._ (2022)

##  IV Examples

The location of the muon site in a solid is chiefly determined by the
electrostatic interactions of the positively-charged muon and its
surroundings. In metals, we expect to find bare muons (often called
diamagnetic muons), uncoupled from unpaired electron density. In insulators,
there is also the possibility of forming muonium (Mu): a bound state of an
electron and muon with a binding energy of 13.54 eV. However, our interest
below is mostly in bare muons. We expect a stopped bare muon to be self
trapped by its own local distortion to the lattice,Blundell  _et al._ (2022)
such that (classically) it sits at the bottom of an electrostatic potential
well. The energy scales chiefly at play determining muon sites are therefore
electrostatic ones of order 0.1–10 eV. This large range of scales, which is
often reflected in the energy differences between candidate muon-containing
structures found using DFT, follows from that found in the cohesive energies
of different types of solid, from molecular crystals (of order 0.1
eV/molecule), through metals (1 eV/atom), to insulators (10 eV/atom).Ziman
(1972) Such scales are well described by electronic-structure methods such as
DFT. Its worth noting that magnetic interactions, whose energy scales are
typically determined by interactions of order tens of millivolts and below,
are rather smaller still. So while the electronic state of a material, and the
consequent distribution of electron density, will likely affect the muon-
stopping state, the magnetic state of a material (and the consequent
distribution of electron spin density) is unlikely to determine the muon site.

###  IV.1 Fluorides

Ionic fluorides were one of the first systemsMöller  _et al._ (2013);
Bernardini  _et al._ (2013) to be studied using DFT+$\mu$. This is because the
muon site in fluoride systems can be identified using the characteristic
experimental signature of an F–$\mu$–F state.Brewer  _et al._ (1986) Even
though there is no electronic magnetism in many of these compounds, a clear
and rather complex precession signal was observed. has some rather special
properties: (i) fluorine is extremely electronegative and thus attractive to
the positively-charged muon, (ii) F- has a very small ionic radius (so that
the muon sits very close), (iii) the fluorine nucleus has spin one-half (so
that any relaxation retains a simple, yet characteristic, time dependence) and
(iv) the fluorine nucleus has a large nuclear moment, associated with an
isotope which occurs with 100% abundance (so that all states behave the same).
It was recognised that the effect could be explained if a muon stops between
two fluorine ions to form what is known as an F–$\mu$–F state. This state can
be described by the Hamiltonian ${\cal H}$ given by

${\cal H}={\cal H}_{\mu{\rm F}}+{\cal H}_{{\rm FF}},$ (15)

the sum of two terms expressing the dipolar interaction between the muon and
the fluorine nuclei

${\cal H}_{\mu{\rm F}}={\mu_{0}\over
4\pi}\sum_{i=1}^{2}{\hbar^{2}\gamma_{\mu}\gamma_{\rm F}\over
r_{i}^{3}}\left[\bm{S}_{\mu}\cdot\bm{S}_{{\rm
F}i}-{3(\bm{S}_{\mu}\cdot\bm{r}_{i})(\bm{S}_{{\rm F}i}\cdot\bm{r}_{i})\over
r_{i}^{2}}\right],$ (16)

and the dipolar interaction between the fluorine nuclei

${\cal H}_{{\rm FF}}={\mu_{0}\over 4\pi}{\hbar^{2}\gamma_{\rm F}^{2}\over
r_{\rm FF}^{3}}\left[\bm{S}_{{\rm F}1}\cdot\bm{S}_{{\rm F}2}-{3(\bm{S}_{{\rm
F}1}\cdot\bm{r}_{\rm FF})(\bm{S}_{{\rm F}2}\cdot\bm{r}_{\rm FF})\over r_{\rm
FF}^{2}}\right].$ (17)

If we ignore the interaction between the two fluorine moments, then
eigenvalues are $-1$, $\frac{1-\sqrt{3}}{2}$, $0$, and $\frac{1+\sqrt{3}}{2}$,
all in units of $\hbar\omega_{\rm d}=\mu_{0}\hbar\gamma_{\mu}\gamma_{\rm
F}/(4\pi r^{3})$ and twice-repeated, so the energy levels are four doublets.
The polycrystalline average polarization is then given by an analytical form

$\langle
P_{z}(t)\rangle=\frac{1}{2}+\frac{1}{6}\sum_{i=1}^{3}a_{i}\cos(\alpha_{i}\omega_{\rm
d}t),$ (18)

with $a_{1}=1$, $a_{2}=1-\frac{1}{\sqrt{3}}$, $a_{3}=1+\frac{1}{\sqrt{3}}$,
$\alpha_{1}=1$, $\alpha_{2}=\frac{3-\sqrt{3}}{2}$, and
$\alpha_{3}=\frac{3+\sqrt{3}}{2}$. If the interaction between the two fluorine
moments is also included then the factors in this expression change a little
bit (and the expression doesn’t look as pretty). Nevertheless, the form of the
oscillatory signal provides information concerning the nature of the
interaction (whether with one fluorine, or two, and whether the F–$\mu$–F bond
is straight or bent) and the timescale of the oscillatory signal is set by
$\omega_{\rm d}$, and hence by $r$, so that the muon-fluorine bond length can
be accurately measured. This provides crucial quantitative information on the
muon site and its local distortion, allowing a rigorous test of the accuracy
of the DFT$+\mu$ calculations.

The oscillations in the muon polarization described in eqn 18 have been
observed experimentally in many inorganic fluorides,Brewer  _et al._ (1986);
Noakes  _et al._ (1993) fluoropolymers,Pratt  _et al._ (2003); Nishiyama,
Nishiyama, and Higemoto (2003); Lancaster  _et al._ (2009) and fluorine-
containing molecular magnets Lancaster  _et al._ (2007). DFT$+\mu$
calculationsMöller  _et al._ (2013) in a variety of inorganic fluorides (LiF,
NaF, CaF2, BaF2, and CoF2) reveal a linear F–$\mu$–F state with a F–F
internuclear distance between 2.31 and 2.36 Å, whose range includes 2.36 Å,
the measured distance of the (FHF)- molecular bifluoride ion in
vacuum.Kawaguchi and Hirota (1987) The bifluoride ion has various vibrational
modes (e.g. a bending mode of frequency around 1290 cm-1), but the calculated
modes for F–$\mu$–F are around a factor of three largerMöller  _et al._ (2013)
(due to the muon mass being approximately one ninth of the proton mass). This
results in the muon in a F–$\mu$–F state having a zero-point energy of
$\approx 0.8$ eV, considerably larger than the value of $\approx 0.3$ eV for
the proton in the bifluoride ion.Möller  _et al._ (2013). We will return to
this point in Section IV.6.

The rigidity of the F–$\mu$–F bond will introduce a distortion into the
crystal. For example, in CaF2, the F–$\mu$ distance is found by experiment to
be 1.172 Å (with the DFT$+\mu$ calculation giving 1.134 Å), but this value is
much lower than $a/4=1.362$ Å, the distance between the muon site and the
fluorine ion in its undistorted position, demonstrating a significant
distortion.Wilkinson and Blundell (2020) The same effect in NaF is illustrated
in Fig. 1. In the undistorted structure, each ion sits 1.64 Å from the site of
the muon; following implantation, the distance from the muon to the nearest-
neighbour fluorine (sodiuum) ion becomes 1.20 Å (2.31 Å). The results of these
DFT$+\mu$ calculations are in excellent agreement with experimental results
for which the complex precession signals due to the dipolar couplings with all
nearby nuclei can now be calculated and high-statistics measurements have now
been performed.Wilkinson and Blundell (2020)

Figure 1:  (a) The crystal structure in NaF. (b) After implantation of the
muon, two nearby fluorine ions are pulled in and two nearby sodium ions are
pushed out. [Adapted from J. Wilkinson and S. J. Blundell, Phys. Rev. Lett.
125, 087201 (2020).Wilkinson and Blundell (2020)] Copyright 2020 American
Physical Society.

This approach has now been extended to other fluorides, including YF3 (see
Fig. 2) and $\alpha$-PbF2; in the latter compound the $\mu$SR data and
analysis provide evidence for a Mu${}^{-}=\mu^{+}\mbox{e}^{-}\mbox{e}^{-}$
stopping site in an anion Frenkel defect.Wilkinson  _et al._ (2021) Although
fluorine is the most ideal ion for exploring these effects, it is not the only
one. $\mu$SR data obtained on V3Si shows evidence for quantum coherent
oscillations which can be related to a V–$\mu$–V state.Bonfà  _et al._ (2021)
These experiments highlight the extreme sensitivity of the entangled states to
the local structural and electronic environments which, in the case of these
A15 compounds containing V (and also Nb), emerges through the quadrupolar
interaction with the electric field gradient. This demonstrates that positive
muons, usually thought of as a purely magnetic probe, can also be deployed as
quantum sensors to measure structural and charge-related phenomena.Bonfà  _et
al._ (2021)

Figure 2: $\mu$SR experimental data on YF3 with the fitted muon polarization
using the muon site obtained by DFT$+\mu$. This site is displayed on the YF3
crystal structure in the inset. [Adapted from J. Wilkinson et al. Phys. Rev. B
104, L220409 (2021).Wilkinson  _et al._ (2021)] Copyright 2021 American
Physical Society.

###  IV.2 Magnetism

Magnetism remains a key research area for $\mu$SR and many applications of the
DFT+$\mu$ techniques have concentrated on systems that are magnetically
ordered. Here a knowledge of the muon site tells us which magnetic fields are
causing the time-evolution of the muon’s spin, potentially giving us access to
more information on the local magnetic field distribution and its dynamics.

A simple, topical example is found in materials adopting the double-perovskite
structure exemplified by antiferromagnetic Sr2FeOsO6.Williams  _et al._ (2016)
In this system the muon site is found from DFT to be approximately 1 Å from an
oxygen atom in the basal plan of the oxygen octahedra surrounding the Os ions,
validating the rule of thumb that such sites close to O2- ions are most likely
in oxides. In this case the site is consistent with the dipole-field map
computed for the known magnetic structure, giving confidence that the DFT-
derived site is likely the one realised. Other ionic systems in which
DFT$+\mu$ calculations show that the muon occupies a site close to the anion
include the honeycomb systemLang  _et al._ (2016) $\alpha$-RuCl3, and the
spin-Jahn-Teller antiferromagnetsKirschner  _et al._ (2019); Lang  _et al._
(2019) CoTi2O5 and FeTi2O5. A similar approach applies in the case of the
ferromagnet Nd2Fe14B, where the muon site is identified [the 8i site
(0.6745,0.8838,0)] near the square base of a NdFe3B pyramid, leading to a
quantitative measurement (using the observed $\mu$SR precession frequencies)
of the moment on the Nd and Fe atoms.Sugiyama  _et al._ (2019)

More stringent tests of DFT-computed sites are found in materials with more
complicated, often noncollinear, magnetic structures and/or hyperfine
contributions. The helimagnet MnSi has been the subject of several
illuminating investigations of its properties using muons Amato  _et al._
(2014); Bonfá , Sartori, and De Renzi (2015) and provides a good example of
site-determination in a more complicated magnetic system, In this material the
muon site computed using DFT methods (shown in Fig. 3) is consistent with the
symmetry properties derived from Knight-shift measurements, and also with the
very restrictive details of the complicated magnetic structure of the
material.Amato  _et al._ (2014) This information proved very useful in
providing a complete analysis of the measured $\mu$SR spectra in terms of a
magnetic structure and hyperfine interaction,Amato  _et al._ (2014) that did
not require the invocation of a more exotic muon stopping state that had been
previously suggested.Storchak  _et al._ (2011)

Figure 3:  The muon site in MnSi is found experimentally to be at
$(0.532,0.532,0.532)$ (indicated by small circles at crystallographically
equivalent sites). These positions are enclosed by the regions where the
unperturbed electrostatic potential computed by DFT takes a minimum. [Adapted
from A. Amato et al., Phys. Rev. B 89, 184425 (2014).Amato  _et al._ (2014)]
Copyright 2014 American Physical Society.

Figure 4:  Left: the four classes of muon stopping site determined for GaV4S8.
Right: local geometry around the muon for each of the four classes of muon
stopping site in GaV4S8. The sites are numbered in order of increasing energy.
[Adapted from K. J. A. Franke et al., Phys. Rev. B 98, 054428 (2018).Franke
_et al._ (2018)] Copyright 2018 American Physical Society.

Noncollinear magnetism magnetic structures are especially relevant to the
study of skyrmion systems. Skyrmions Lancaster (2019) are vortex-like magnetic
excitations found in a growing number of magnetic materials. These include
MnSi, which hosts Bloch skyrmions. Another example skyrmion compound, this one
hosting Néel skyrmions, is GaV${}_{4}X_{8}$ ($X$=S, Se). Structural
relaxations of a supercell of GaV4S8Franke  _et al._ (2018) reveal four
distinct muon sites [Fig. 4]. Three of these (labeled I-III) involve the muon
sitting close to a single S atom. A fourth site (site IV) has the muon closer
to V atoms and is the highest energy site. In the lowest energy site (site I)
the muon sits between two S atoms, in the plane defined by three S atoms
within V4S4 units. The two $\mu^{+}$–S distances are unequal (1.4 Å and 2.0 Å)
with greater electron density found between the muon and the nearest S atom.
This site is therefore best described in terms of the muon forming a
$\mu^{+}$–S bond (rather than an S–$\mu^{+}$–S state), though the presence of
a second nearby S atom does seem to stabilize this geometry. Two further sites
involve the muon sitting close to a single S atom.

Perhaps the most pressing problem with the DFT+$\mu$ method is determining
which sites are actually realized given an energy-ordered list of candidate
states. This problem is particularly relevant for the results of analogous
calculations for GaV4Se8. These give stopping sites (labeled $1$ \- $4$ in
order of ascending energy) that are similar to those calculated for GaV4S8,
with three of the four sites involving the muon sitting close to a Se atom
(sites $2$ to $4$) and a site in which the muon sits above a face of a V4Se4
unit (site 1). However, the ordering of sites is inverted in the Se case,
compared to the S-containing series. In particular, the cube face site (site
1), which corresponds to the highest energy stopping site for GaV4S8, is the
lowest energy site for GaV4Se8. Although the difference in energy ordering is
interesting, without access to a method to compute the capture cross section
for each site, the interpretation of this information remains an open
question.

So far we have seen examples of the commonly encountered cases of sites near
fluorine and near oxygen. An interesting case of a system where both sites
near oxygens and fluorines are predicted is barlowite, a frustrated magnet
kagome antiferromagnet with formula [Cu4(OH)6FBr].Tustain  _et al._ (2020) In
this system, DFT suggests two distinct classes of muon stopping site. The
first localises $\approx 1.0$ Å away from the oxygen atoms in the hydroxide
groups that connect the Cu2+ ions within the kagome layers. This forms a
triangular $\mu$–OH complex, with the muon-proton distance of 1.54 Å in the
lowest-energy sites. In the second class, muons localise near the fluoride
anions in between the kagome layers, with a $\mu$–F separation of 1.1 Å. These
latter sites lie at substantially higher energies above the lowest energy
sites in the calculations ($\approx 1$ eV). Although this suggests, on purely
energetic grounds, that the formation of $\mu$–F complexes in barlowite is
unlikely compared to $\mu$–OH. However, the measured spectra in the
paramagnetic regime show clear dipole-dipole oscillations consistent with a
$\mu$–F complex, implying that during the stopping process muons are indeed
captured in these latter potential minima.

The study of H-like defects in semiconductors and insulators using first-
principles methods is a mature field,Van de Walle (1994, 2000); Freysoldt  _et
al._ (2014) and understanding the analogous muon-like states is clearly
related.Cox (2009); Hiraishi  _et al._ (2022) A recent development has come
from the suggestion that in certain oxides, such as Cr2O3, the muon can form a
charge-neutral complex composed of a muon and an electron polaron.Dehn  _et
al._ (2020) Cr2O3 contains Cr3+ (3d3) ions, but the idea is that the $\mu^{+}$
binds to an oxygen anion but an excess electron localizes on an adjacent Cr
ion, changing its charge state to Cr2+ (3d4), which is Jahn-Teller active,
resulting in a Jahn-Teller polaron. The existence of this state is supported
by first-principles calculations, and an analogous effect has been
proposedDehn  _et al._ (2021) for Fe2O3.

###  IV.3 Superconductivity

An important use of muon in superconducting materials is the investigation of
the vortex lattice in type-II superconductors, where the muon can be used to
measure the penetration depth. In such measurements, the muon site is often of
secondary importance since the length scale of the vortex state is large on
the scale of the underlying unit cell, so the muons effectively probe the
whole of the flux lattice, whatever their stopping site in the crystal.
However, the physics of unconventional superconductors is inextricably bound
up with magnetic phenomena (and sometimes questions of phase separation and
coexistence), and so in $\mu^{+}$SR studies of materials showing
unconventional superconductivity, it is usually absolutely necessary to have
knowledge of the muon stopping state.

In the cuprate superconductors, it has been thought for a long time that the
muon site will be about 1 Å away from an oxygen anion.Brewer  _et al._ (1990)
A recent study of muon sites in La2CuO4 has concluded that most muons stop
near the apical oxygen, but there are two additional sites that also receive a
tiny fraction of muons.Ramadhan  _et al._ (2022) Iron-based superconductors
provide another informative example of an unconventional superconducting
series of materials.De Renzi  _et al._ (2012); Cheung  _et al._ (2018) The DFT
calculation allows a straightforward identification of the muon site using the
UEP approach. In the case of the 1111 structure this allows the identification
of three sets of potential minima shown in Fig. 5. Similar results have been
obtained forSundar  _et al._ (2023) UTe2.

Figure 5:  Muon sites identified in the LaFeAsO unit cell. Gold-shaded areas
represent the volume of zero point displacement for the muon. [Taken from R.
De Renzi et al., Supercond. Sci. Technol. 25, 084009 (2012).De Renzi  _et al._
(2012)] Copyright 2012 Institute of Physics Publishing.

One area where $\mu$SR has proven important is in attempts to determine
instances of time-reversal symmetry breaking (TRSB) in superconducting
systems, since this can provide a tight constraint on the symmetry of the
superconducting gap.Huddart  _et al._ (2021a) The assignment of TRSB is made
through the appearance of spontaneous magnetic fields found using $\mu$SR
measurements. (Importantly, this effect is absent for most superconductors.)
The effects observed are frequently small, and sometimes only seen in muon
measurements and so it is natural to ask whether a muon-induced perturbation
is a contributing factor to the observations. A recent DFT+$\mu$ study of
these materialsHuddart  _et al._ (2021a) suggest that the muon is an innocent
probe of these materials. We shall concentrate on one candidate TRSB material:
the layered perovskite superconductor Sr2RuO4.

Figure 6:  The lowest-energy muon site in Sr2RuO4. (a) Local geometry of the
muon site. (b) Radial displacements of atoms as a function of their distances
from the muon site. PDOS for the Ru atom closest to the muon site for the
structures (c) without a muon and (d) with a muon. Energies are given with
respect to the Fermi energy. [Adapted from B. M. Huddart et al., Phys. Rev.
Lett. 127, 237002 (2021).Huddart  _et al._ (2021a)] Copyright 2021 American
Institute of Physics.

In the lowest-energy muon site found for Sr2RuO4, the muon is bonded to an
oxygen (O2) with bond distance 0.973 Å, consistent with muon sites in other
oxides including high-temperature superconducting cuprates and pyrochlores.
The small muon-induced displacement vanishes rapidly as a function of distance
from the muon site, such that significant distortions are observed only for
atoms within 6 Å of the muon site. We can therefore conclude that there is
little structural distortion, but what about a distortion to the physics at
the next largest energy scale: the electronic structure? It is found that the
dominant contribution to the electronic density of states (DOS) close to the
Fermi energy is that from the Ru atoms. The effect of muon implantation on the
projected DOS of the Ru atom closest to the muon site is shown in Figs. 6.
There is an increase in the DOS at around 1 eV below the Fermi energy caused
by small changes in the splitting of the Ru 4dzy and 4dzx states at the Fermi
level, which are not observed for Ru atoms further away from muon. However,
after summing the d-state contributions from all of the Ru ions in the
supercell, the small state splitting is no longer resolvable. The partial
density of states corresponding to the muon itself lies around 8 eV below the
Fermi energy. Since this is very large interval in energy, it is highly
unlikely the muon has any effect on the electrons near the Fermi energy and so
these results suggest that the implanted muon does not have a significant
effect on Sr2RuO4. A similar conclusion was reached for other TRSB systems
investigated in this study.Huddart  _et al._ (2021a)

###  IV.4 Molecular systems

Figure 7:  Low-energy muon sites in Cu-PM. Lighter spheres represent the ionic
positions in the unit cell without the muon. H atoms have been omitted for
clarity where appropriate. (a) The nitrate site. (b) The H2O site. (c) The
N(pym) site. [Taken from B. M. Huddart et al., Phys. Rev. B 103, L060405
(2021).Huddart  _et al._ (2021b)] Copyright 2021 American Institute of
Physics.

Muons are routinely used to elucidate the properties of molecular materials.
These systems are attractive owing to their tunability, but their chemical
complexity regularly leads to questions about the nature of the muon stopping
state. These materials are rather soft in comparison to their inorganic
counterparts, and so we might expect that the muon-induced distortion could
present a relatively large perturbation to the electronic or magnetic
properties of the system.

One class of these materials are coordination polymer magnets, where a
magnetic ion is linked via molecular ligands to form a low-dimensional magnet.
These are of interest since the low-dimensional character often suppresses
long-range magnetic order, reducing the size of the response of conventional
measurement techniques owing to the small moment size and small change in
entropy upon ordering. In these systems muons are often able to detect
transitions that prove invisible to magnetic susceptibility or specific heat.

The staggered molecular spin chain [pym-Cu(NO3)2(H2O)2] (pym = pyrimidine),
known as Cu-PM is a good example.Huddart  _et al._ (2021b) Here magnetic order
was detected using muons at $T=0.23$ K. Three distinct classes of muon
stopping site were determined [Fig. 7]. Sites where the muon sits around 1 Å
from an O atom in a nitrate group [Fig. 7(a)] or H2O ligand [Fig. 7(b)] are
the lowest and second lowest energy classes of sites, respectively. We also
find candidate sites where the muon sits 1.0 Å from an N atom in a pym ligand
[Fig. 7(c)], which are substantially higher in energy and result in larger
local distortions to the crystal structure. These sites can be mapped to
features in the ZF spectra by considering the dipolar fields resulting from
candidate antiferromagnetic structures. From dipolar field calculations one
obtains fields of 9–40 mT/$\mu_{\mathrm{Cu}}$ for the nitrate site, 57–63
mT/$\mu_{\mathrm{Cu}}$ for the H2O site, and 93–99 mT/$\mu_{\mathrm{Cu}}$ for
the N(pym) site, where $\mu_{\mathrm{Cu}}$ is the ordered moment of the Cu2+
ions in Bohr magnetons. The relative size of the calculated fields for the H2O
and nitrate sites is consistent with the ratio between observed frequencies
measured in the ordered regime. This assignment gives an estimate
$\mu_{\mathrm{Cu}}$ $\approx$ 0.38$\mu_{\mathrm{B}}$ for the ordered
moment.Huddart  _et al._ (2021b)

There have been claims,Eggert and Affleck (1992) supported by some
experimental evidence,Chakhalian  _et al._ (2003) that a muon implanted in
some spin chains could form an unusual spin-singlet state, not unlike that
found in the Kondo effect. However, results on materials such as Cu-PM and the
linear chain Cu(pyz)(NO3)2, Xiao  _et al._ (2015) do not provide evidence for
the realization of these states. An example where a substantial muon-induced
perturbation seems to be present is the molecular spin-ladder material
(Hpip)2CuBr4,Lancaster  _et al._ (2018) where the muon forms states based on
Br– $\mu$–Br bonds. This state causes a sizeable local distortion to the atoms
neighbouring the muon and consequently to the electronic structure (Fig. 8).
Although it might be feared that this would prevent the muon from faithfully
measuring the properties of the material, this seems not to be the case as the
magnetic phase diagram is determined by physics on a length scale that is very
long compared to the muon and its distortion. As a result the transitions in
the material are observed with $\mu$SR at the fields and temperatures found
using other techniques.

Figure 8:  Example muon site in the molecular spin ladder (Hpip)2CuBr4. [Taken
from T. Lancaster et al., New. J. Phys. 20, 103002 (2018).Lancaster  _et al._
(2018)] Copyright 2018 Institute of Physics Publishing.

An example of a rather different structure is the organic radical
2-(4,5,6,7-tetrafluorobenzimidazol-2-yl)-4,4,5,5-tetramethyl-4,5-dihydro-1$H$-imidazole-3-oxide-1-oxyl
(F4BImNN)Blundell  _et al._ (2013) which forms hydrogen-bonded chains and
exhibits one-dimensional ferromagnetic exchange. Here three energetically
similar sites are predicted (i) where the muon forms a covalent bond with the
oxygen in the nitronyl nitroxide group (see Fig. 4), (ii) where the muon forms
a covalent bond with the vacant nitrogen and hydrogen-bonds to one oxygen in
the nitronyl nitroxide group (this is the lowest energy site), and (iii) where
the muon bonds to the other oxygen in the nitronyl nitroxide group.
Interestingly, despite the presence of fluorine in this sytem, an F–$\mu$–F
state is not predicted.

###  IV.5 Muon-induced distortions

This example of a spin-ladder material from the previous section is a case
where the muon-induced distortion likely leads to a probe state comprising
both the muon and its surrounding perturbation that is nevertheless sensitive
to the intrinsic properties of the material. However, it is natural to ask
about other cases where the muon has a strong effect on its surroundings.

Figure 9: Muon stopping site (white sphere) and atomic positions (oxygen ions
are red, Pr ions are yellow, and Sn is suppressed) in the quantum spin ice
Pr2Sn2O7 as calculated by DFT$+\mu$. Bonds are shown as solid black lines,
with the bonds of the unperturbed lattice (yellow lines) included for
comparison. Pr ions are numbered in order of separation from $\mu$ with Pr1
being the closest. The pair Pr2 and Pr 2’ is equidistant from the muon, as are
Pr4 and Pr 4’. [Adapted from F. R. Foronda et al. Phys. Rev. B 92, 134517
(2015).Foronda  _et al._ (2015)] Copyright 2015 American Physical Society.

An important example of this was found in a family of quantum spin ices, a
muon-induced change of the local crystal field dominates the measured response
in a way that could be fully quantified,Foronda  _et al._ (2015) demonstrating
that the concerns raised above can be of critical importance. Spin ice
materials are based on the pyrochlore structure and the physics of such
pyrochlore materials, with general formula A2B2O7 is strongly dependent on the
crystal field surrounding the lanthanide cation A. Thus the introduction of a
muon into the structure could alter that crystal field and change the nature
of the ground state.

This effect is realised in the material Pr2Sn2O7, which exhibits quantum spin
ice behaviour.Zhou  _et al._ (2008) Without the muon, the Pr3+ ion is
surrounded by eight oxygen anions and adopts a doublet crystal field ground
state, resulting in an effective spin-$\frac{1}{2}$ moment.Princep  _et al._
(2013) However, DFT$+\mu$ calculations imply that the presence of the muon
[the site is identified as $(-0.013,0.047,0.203)$] distorts the local symmetry
around nearby Pr ions. The largest effect is found to be an anisotropic
distortion, with one Pr-–O bond strongly rotated (bent) and another
significantly extended (see Fig. 9). The net result is that the doublet ground
state is split on each of the neighbouring Pr cations has a singlet ground
state. Then, via a hyperfine enhancement mechanism,Bleaney (1973); Blundell
(2023) a distribution of static magnetic moments is produced in which the
average moment size grows on cooling. The model (involving DFT$+\mu$
calculations of the distortion and crystal field calculations of the resulting
Pr environments) produced quantitative agreement with the temperature
dependence of the observed muon response, demonstrating the validity of this
approach.Foronda  _et al._ (2015) The Pr3+ ion has a non-Kramers ground state
and so is particularly vulnerable to electrostatic perturbations. For the case
of Dy2Ti2O7, the Kramers doublet ground state of Dy3+ is unaffected by the
muon, even though the muon site and resulting anisotropic distortion are
pretty much identical to that shown in Fig. 9 for Pr2Sn2O7.

###  IV.6 Quantum effects

The DFT treatments described so far have assumed that the muon is a classical
particle sitting statically at its equilibrium position. This is a consequence
of the Born-Oppenheimer approximation where we make a clear distinction
between the classical nuclei (including the positive muon, which is treated as
a light proton) and the quantum electrons. However, this approach neglects
quantum mechanical effects of muons and nuclei, such as zero-point motion,
tunnelling, and large zero-point energies (ZPEs). These can change the
energetic ordering of muon sites, destabilize certain classical muon sites or
cause several sites to merge into one if ZPE overcomes the potential-energy
barriers between them. For nuceli, all of these effects are mass dependent and
therefore especially pronounced for light particles such as the muon. The
result is that, in reality, the muon wavefunction is expected to be spread
over an appreciable volume in the material, rather than concentrated at a
single point. All observable quantities such as muon coupling constants have,
in principle, to be suitably averaged over the whole extended muon
wavefunction, rather than sampled at a single classical point as in the Born-
Oppenheimer approximation. The significant spatial extent of the wavefunction
of the implanted muon is known to be relevant in many systems. Well known
examples are the quantum diffusion of muons in metallic systems such as copper
or of muonium in insulating systems such as solid nitrogen.

Dealing with quantum effects is, of course, a difficult problem (DFT being a
scheme invented to avoid this!) and so several approximations have been
proposed and attempted. Perhaps the simplest class of approximation is a
harmonic one. This starts with a DFT calculation of the ($\Gamma$-point)
phonon spectrum for the classical muon site geometry. An additional assumption
often made is that, since muons are lighter than nuclei, muon zero-point
motion is adiabatically decoupled from the lattice over some relatively short
length scale.

An example is, once again, provided by the fluorine states described in
Section IV.1. In vacuum the linear F–$\mu$–F anion has four vibrational modes:
a symmetric stretch, bending (two-fold degenerate), and an asymmetric stretch.
In a solid the twofold degeneracy of the bending mode is broken due to the
symmetry of the site: in LiF, NaF, CaF2, and BaF2 there are two fundamentally
inequivalent directions of bending; one is towards a neighboring cation and is
shifted up in frequency while the other direction is into a gap in the crystal
structure and is shifted down in frequency. From the frequencies of these
decoupled vibrational modes we can estimate the zero-point energy (ZPE) of the
muon. The strong bond in this case, combined with the small muon mass leads to
a large ZPE for the F–$\mu$–F center of 0.80 eV in vacuum. This is larger than
the ZPE of any natural triatomic molecule (the ZPEs of H2O and H${}^{+}_{3}$
are 0.56 and 0.54 eV, respectively), demonstrating the importance of quantum
effects in muon localization.Möller  _et al._ (2013)

More generally, the quantum muon problem can be addressed using a range of
single-particle approximations where the muon sits in the potential of its
surroundings and the Schrödinger equation is solved for the muon only. These
schemes can be most easily addressed in two limits: weakly- and strongly-bound
muons. In both of these, an effective single-particle potential is constructed
from total DFT energy. To do this, the muon is displaced from its classical
site, while: (i) keeping the nuclei fixed at the positions corresponding to
the unperturbed muon site (in the weakly-bound case), or (ii) letting the
nuclei relax to new lowest-energy positions (in the strongly-bound case),
while keeping the center of mass of the system fixed.

An example of an attempt to apply these approaches is the case of
$\alpha$-N2.Gomilšek  _et al._ (2022) Here the muon forms an extended
electric-dipole polaron around a central, quantum-entangled [N2–$\mu$–N2]+
complex. The ZPE in the approximation schemes is shown in Fig. 10(a) and (b).
The muon wavefunction is estimated to be significant on the length scale of
the underlying crystal, making the quantum calculation necessary.
Unfortunately, the two single-particle schemes give rather different estimates
of the ZPE. The phonon calculation [shown in Fig. 10(c) as a function of muon
mass] shows the problem here: a strong hybridization of the longitudinal muon
normal mode with intra-molecular vibrations of both N2 in the [N2–$\mu$–N2]+
complex, which implies significant muon–nuclear entanglement. As a result,
both the single particle and harmonic approximations provide only a very rough
guide to the quantum effects and do not agree between themselves. In cases
such as that of N2, an exact approach is therefore needed. This is carried out
using path-integral molecular dynamics (PIMD), in which observables are
calculated from arbitrary muon–nuclear zero-point motion. This scheme, based
on discretizing an imaginary-time path integral, is very intensive on computer
power, but can yield reliable numerical estimates. In the case of N2 it
confirms that the muon and its surroundings are strongly entangled, but leads
to a very precise estimate of the 14N nuclear quadrupolar coupling
constant.Gomilšek  _et al._ (2022)

Figure 10:  (a,b) Effective muon potential in $\alpha$–N2 on a square-root
linear scale (i.e. harmonic potentials appear as straight lines) along
different directions for different approximations (left axis) with the
corresponding muon probability density (right axis), directional
contributions, mean muon displacements, wavefunction widths and inflection
points. (c) $\Gamma$-point phonon spectra for muonated $\alpha$–N2 vs. muon
mass. Shading gives the muon amplitude squared in the normal modes. [Taken
from M. Gomilšek et al., arXiv:2202.05859.Gomilšek  _et al._ (2022)]

Continuing the theme of quantum mechanical effects, transition-state searches
provide a route to evaluating the possibility of muon diffusion between
candidate sites (an early example of this was applied to the study of muon
diffusion in copperBonfá , Sartori, and De Renzi (2015)). In a general
geometry optimization process, the coordinates of the atoms are adjusted so
that the energy of the structure is brought to a stationary point, in which
the forces on the atoms are small. A transition state (TS) is a stationary
point that is an energy maximum in one direction and an energy minimum in all
other directions. A procedure that automates a search for TSs in codes such as
Castep Clark  _et al._ (2005) is particularly useful for predicting barriers
to chemical reactions and determining reaction pathways. In a chemical
reaction, for example, starting from reactants, energy increases to a maximum
and then decreases to the energy of the products. The maximum energy along the
reaction pathway is the activation energy and the structure corresponding to
this energy is the transition state. These ideas can be used to find diffusion
barriers, which is where its use in $\mu$SR arises.

Figure 11:  Muon stopping sites in TaS2. Ta and S are represented by blue and
yellow spheres respectively. Sites 1 and 2 occupy the interlayer region. Site
3 is associated with the S layer and site 4 with the Ta layer. [Taken from F.
L. Pratt et al., npj Quantum Materials 6, 69 (2021).Mañas-Valero  _et al._
(2021)]

A recent example of the use of transition-state analysis to evaluate muon
diffusion in a solid is the case of TaS2.Mañas-Valero  _et al._ (2021) Here
the muon sites form four distinct groups (Fig. 11). Muons in site 1 bond to
one S atom in each of two adjacent TaS2 layers forming a linear S–$\mu$–S
state with unequal bond lengths. In site 2 the muon is bonded to only a single
S atom (with a bond length of 1.4 Å), with these sites being around 0.03 eV
higher in energy than site 1. In site 3, the muon stops inside the S layer,
slightly displaced from the centre of the triangle defined by three S atoms.
For site 4 the muon site is displaced by around 0.7 Å along the $c$-axis above
the centre of the triangle defined by three Ta atoms. To investigate the
possibility of muon diffusion, TS searches between each pair of distinct muon
stopping sites associated with a single TaS2 layer were carried out. The
energy barrier between sites 1 and 2 is very close to the difference in energy
between the two. This, coupled with the large estimated ZPEs of both sites
(0.35 eV and 0.51 eV respectively), means that the muon is likely to be
delocalised between these two geometries, rather than there being two distinct
stopping sites. Similarly, the barrier between sites 2 and 3 is less than
their ZPE. Thus, sites 1, 2 and 3 are expected to form a single quantum
delocalised state.

##  V Outlook and conclusion

DFT+$\mu$ is computationally costly, particularly for complex materials, and
each new material investigated has to be studied from scratch in the same way.
Many of the most exciting materials that are currently being discovered
(including spin liquids) tend to be chemically complex, and this puts severe
demands on the feasibility of electronic structure calculations. Although we
have concentrated on calculations of the muon site and associated local
distortions, DFT calculations can also be used to estimate muon contact
hyperfine fields in metals.Onuorah, Bonfà, and De Renzi (2018); Onuorah  _et
al._ (2019)

The $\mu$SR technique concentrates on using the positive muon, but experiments
with negative muons can also be performed.Blundell  _et al._ (2022) There is
no problem with identifying the $\mu^{-}$ site because it is captured by the
nucleus. This does result in an atom with atomic number $Z$ behaving as a
muonic atom with effective nuclear charge $+(Z-1)e$, the $\mu^{-}$ partially
screening the full nuclear charge. The effect on this state on the local
crystal structure can then be calculated using DFT.Gill, Pratt, and Blundell
(2023) (This method can be termed DFT$+\mu^{-}$ to distinguish it from the
usual DFT$+\mu^{+}$ technique.)

A notable feature of the method used in (positive muon) DFT$+\mu$ of randomly-
initializing muons in sites and then relaxing the structure, is that it is
rather unlike the real-life situation of the muon finding its stopping site by
impinging on the sample with an energy of 4 MeV and then shedding energy
through interactions with the solid, before eventually coming to rest.Blundell
_et al._ (2022) This means that it is not possible, given an energy-ordered
list of calculated muon sites, to tell which is realised, nor to reliably
estimate the relative occupancy of several multiply occupied sites. Ultimately
a capture cross section for each site is needed, whose value will reflect not
only the energy, but also some details of the final stages of the muon
stopping process. This is likely to be an avenue of future research whose
solution will allow further information to be squeezed from results obtained
using the $\mu$SR technique. Nevertheless, the progress obtained within the
last few years using DFT$+\mu$ has already transformed the way in which muon
data are both interpreted and understood. Though this method is in its
infancy, it has become an indispensable part of how muon experiments are now
conducted.

###### Acknowledgements.

We are grateful to the following students, postdocs, and colleagues, who have
collaborated with us in our research in this area and made important
contributions to this field: Pietro Bonfà, Davide Ceresoli, Stewart Clark,
Roberto De Renzi, Marina Filip, Francesca Foronda, George Gill, Matjaž
Gomilšek, Zachary Hawkhead, Alberto Hernandez Melían, Thomas Hicken, Benjamin
Huddart, Dominik Jochym, Mohammad Maikudi Isah, Franziska Kirschner, Franz
Lang, Leandro Liborio, Nicola Marzari, Johannes Möller, Ifeanyi John Onuorah,
Francis Pratt, Samuele Sanna, Simone Sturniolo, Johnny Wilkinson, and Hank Wu.
We acknowledge funding from EPSRC(UK). SJB acknowledges funding from UK
Research and Innovation (UKRI) under the UK government’s Horizon Europe
funding guarantee [grant number EP/X025861/1].

## References

  * Cox (1987) S. F. J. Cox, “Implanted muon studies in condensed matter science,” J. Phys. C 20, 3187–3319 (1987).
  * de Réotier and Yaouanc (1997) P. D. de Réotier and A. Yaouanc, “Muon spin rotation and relaxation in magnetic materials,” J. Phys.: Condens. Matter 9, 9113–9166 (1997).
  * Blundell (1999) S. J. Blundell, “Spin-polarized muons in condensed matter physics,” Contemp. Phys. 40, 175–192 (1999).
  * Blundell  _et al._ (2022) S. J. Blundell, R. D. Renzi, T. Lancaster, and F. L. Pratt, eds.,  _Muon Spectroscopy - An Introduction_ (Oxford University Press, Oxford, 2022).
  * Garwin, Lederman, and Weinrich (1957) R. L. Garwin, L. M. Lederman, and M. Weinrich, “Observations of the failure of conservation of parity and charge conjugation in meson decays: the magnetic moment of the free muon,” Phys. Rev. 105, 1415–1417 (1957).
  * Topping and Blundell (2019) C. V. Topping and S. J. Blundell, “AC susceptibility as a probe of low-frequency magnetic dynamics,” J. Phys.: Condens. Matter 31, 013001 (2019).
  * Boothroyd (2019) A. T. Boothroyd,  _Principles of Neutron Scattering from Condensed Matter_ (Oxford University Press, Oxford, 2019).
  * Sonier, Brewer, and Kiefl (2000) J. E. Sonier, J. H. Brewer, and R. F. Kiefl, “$\mu$SR studies of the vortex state in type-II superconductors,” Rev. Mod. Phys. 72, 769–811 (2000).
  * Blundell (2004) S. J. Blundell, “Muon-spin rotation studies of electronic properties of molecular conductors and superconductors,” Chem. Rev. 104, 5717–5735 (2004).
  * Möller  _et al._ (2013) J. S. Möller, D. Ceresoli, T. Lancaster, N. Marzari, and S. J. Blundell, “Quantum states of muons in fluorides,” Phys. Rev. B 87, 121108 (2013).
  * Bernardini  _et al._ (2013) F. Bernardini, P. Bonfà, S. Massidda, and R. De Renzi, “Ab initio strategy for muon site assignment in wide band gap fluorides,” Phys. Rev. B 87, 115148 (2013).
  * Möller  _et al._ (2013) J. S. Möller, P. Bonfà, D. Ceresoli, F. Bernardini, S. J. Blundell, T. Lancaster, R. D. Renzi, N. Marzari, I. Watanabe, S. Sulaiman, and M. I. Mohamed-Ibrahim, “Playing quantum hide-and-seek with the muon: localizing muon stopping sites,” Physica Scripta 88, 068510 (2013).
  * Bonfà and De Renzi (2016) P. Bonfà and R. De Renzi, “Toward the computational prediction of muon sites and interaction parameters,” J. Phys. Soc. Jpn. 85, 091014 (2016).
  * Kohn (1999) W. Kohn, “Nobel lecture: Electronic structure of matter—wave functions and density functionals,” Rev. Mod. Phys. 71, 1253–1266 (1999).
  * Brewer  _et al._ (1990) J. H. Brewer, R. F. Kiefl, J. F. Carolan, P. Dosanjh, W. N. Hardy, S. R. Kreitzman, Q. Li, T. M. Riseman, P. Schleger, H. Zhou, E. J. Ansaldo, D. R. Noakes, L. P. Le, G. Luke, Y. J. Uemura, K. Hepburn-Wiley, and C. E. Stronach, “Site of the positive muon in YBa2Cu3O7,” Hyp. Int. 63, 177–182 (1990).
  * Demaison, Herman, and Lievin (2007) J. Demaison, M. Herman, and J. Lievin, “The equilibrium OH bond length,” International Reviews in Physical Chemistry 26, 391–420 (2007).
  * Brewer  _et al._ (1986) J. H. Brewer, S. R. Kreitzman, D. R. Noakes, E. J. Ansaldo, D. R. Harshman, and R. Keitel, “Observation of muon-fluorine “hydrogen bonding" in ionic crystals,” Phys. Rev. B 33, 7813–7816 (1986).
  * Luke  _et al._ (1991) G. M. Luke, J. H. Brewer, S. R. Kreitzman, D. R. Noakes, M. Celio, R. Kadono, and E. J. Ansaldo, “Muon diffusion and spin dynamics in copper,” Phys. Rev. B 43, 3284–3297 (1991).
  * Blundell  _et al._ (2012) S. J. Blundell, A. J. Steele, T. Lancaster, J. D. Wright, and F. L. Pratt, “A Bayesian Approach to Magnetic Moment Determination Using $\mu$SR,” Physics Procedia 30, 113–116 (2012).
  * Steele  _et al._ (2011) A. J. Steele, P. J. Baker, T. Lancaster, F. L. Pratt, I. Franke, S. Ghannadzadeh, P. A. Goddard, W. Hayes, D. Prabhakaran, and S. J. Blundell, “Low-moment magnetism in the double perovskites Ba2MOsO6 (M = Li, Na),” Phys. Rev. B 84, 144416 (2011).
  * Disseler (2014) S. M. Disseler, “Direct evidence for the all-in/all-out magnetic structure in the pyrochlore iridates from muon spin relaxation,” Phys. Rev. B 89 (2014).
  * Prando  _et al._ (2020) G. Prando, P. Telang, S. D. Wilson, M. J. Graf, and S. Singh, “Monopole-limited nucleation of magnetism in Eu2Ir2O7,” Phys. Rev. B 101 (2020).
  * Van Vleck (1948) J. H. Van Vleck, “The dipolar broadening of magnetic resonance lines in crystals,” Phys. Rev. 74, 1168–1183 (1948).
  * Abragam (1961) A. Abragam,  _Principles of Nuclear Magnetism_ (Oxford University Press, Oxford, 1961).
  * Slichter (1990) C. P. Slichter,  _Principles of Magnetic Resonance_ , 3rd ed. (Springer, Berlin, 1990).
  * Giustino (2014) F. Giustino,  _Materials Modelling using Density Functional Theory: Properties and Predictions_ (Oxford University Press, Oxford, 2014).
  * Hohenberg and Kohn (1964) P. Hohenberg and W. Kohn, “Inhomogeneous electron gas,” Phys. Rev. 136, B864–B871 (1964).
  * Kohn and Sham (1965) W. Kohn and L. J. Sham, “Self-consistent equations including exchange and correlation effects,” Phys. Rev. 140, A1133–A1138 (1965).
  * Becke (1993) A. D. Becke, “Density-functional thermochemistry. III. The role of exact exchange,” J. Chem. Phys 98, 5648–5652 (1993).
  * Perdew, Burke, and Ernzerhof (1996) J. P. Perdew, K. Burke, and M. Ernzerhof, “Generalized gradient approximation made simple,” Phys. Rev. Lett. 77, 3865–3868 (1996).
  * Burke (2012) K. Burke, “Perspective on density functional theory,” J. Chem. Phys. 136, 150901 (2012).
  * Freysoldt  _et al._ (2014) C. Freysoldt, B. Grabowski, T. Hickel, J. Neugebauer, G. Kresse, A. Janotti, and C. G. Van de Walle, “First-principles calculations for point defects in solids,” Rev. Mod. Phys. 86, 253–305 (2014).
  * Leslie and Gillan (1985) M. Leslie and N. J. Gillan, “The energy and elastic dipole tensor of defects in ionic crystals calculated by the supercell method,” J. Phys. C 18, 973 (1985).
  * Makov and Payne (1995) G. Makov and M. C. Payne, “Periodic boundary conditions in ab initio calculations,” Phys. Rev. B 51, 4014–4022 (1995).
  * Van de Walle (1994) C. G. Van de Walle, “Energies of various configurations of hydrogen in silicon,” Phys. Rev. B 49, 4579–4585 (1994).
  * Giannozzi  _et al._ (2009) P. Giannozzi, S. Baroni, N. Bonini, M. Calandra, R. Car, C. Cavazzoni, D. Ceresoli, G. L. Chiarotti, M. Cococcioni, I. Dabo, A. D. Corso, S. de Gironcoli, S. Fabris, G. Fratesi, R. Gebauer, U. Gerstmann, C. Gougoussis, A. Kokalj, M. Lazzeri, L. Martin-Samos, N. Marzari, F. Mauri, R. Mazzarello, S. Paolini, A. Pasquarello, L. Paulatto, C. Sbraccia, S. Scandolo, G. Sclauzero, A. P. Seitsonen, A. Smogunov, P. Umari, and R. M. Wentzcovitch, “QUANTUM ESPRESSO: a modular and open-source software project for quantum simulations of materials,” J. Phys.: Condens. Matter 21, 395502 (2009).
  * Huddart  _et al._ (2022) B. Huddart, A. Hernández-Melián, T. Hicken, M. Gomilšek, Z. Hawkhead, S. Clark, F. Pratt, and T. Lancaster, “MuFinder: A program to determine and analyse muon stopping sites,” Computer Physics Communications 280, 108488 (2022).
  * Clark  _et al._ (2005) S. J. Clark, M. D. Segall, C. J. Pickard, P. J. Hasnip, M. I. J. Probert, K. Refson, and M. C. Payne, “First principles methods using CASTEP,” Zeitschrift für Kristallographie - Crystalline Materials 220, 567–570 (2005).
  * Bonfà, Onuorah, and Renzi (2018) P. Bonfà, I. J. Onuorah, and R. D. Renzi, “Introduction and a quick look at muesr, the magnetic structure and muon embedding site refinement suite,” JPS Conf. Proc. 21, 011052 (2018).
  * Liborio, Sturniolo, and Jochym (2018) L. Liborio, S. Sturniolo, and D. Jochym, “Computational prediction of muon stopping sites using ab initio random structure searching (AIRSS),” J. Chem. Phys. 148, 134114 (2018).
  * Pickard and Needs (2011) C. J. Pickard and R. J. Needs, “Ab initio random structure searching,” Journal of Physics: Condensed Matter 23, 053201 (2011).
  * Sturniolo and Liborio (2020) S. Sturniolo and L. Liborio, “Computational prediction of muon stopping sites: A novel take on the unperturbed electrostatic potential method,” J. Chem. Phys. 153, 044111 (2020).
  * Lancaster and Blundell (2014) T. Lancaster and S. Blundell,  _Quantum Field Theory for the Gifted Amateur_ (OUP Oxford, 2014).
  * Sturniolo, Liborio, and Jackson (2019) S. Sturniolo, L. Liborio, and S. Jackson, “Comparison between density functional theory and density functional tight binding approaches for finding the muon stopping site in organic molecular crystals,” J. Chem. Phys. 150, 154301 (2019).
  * Ziman (1972) J. M. Ziman,  _Principles of the theory of solids, 2nd edition_ (Cambridge University Press, Cambridge, 1972).
  * Noakes  _et al._ (1993) D. R. Noakes, E. J. Ansaldo, S. R. Kreitzman, and G. M. Luke, “The (F$\mu$F)- ion in solid fluorides,” J. Phys. Chem. Solids 54, 785 (1993).
  * Pratt  _et al._ (2003) F. L. Pratt, S. J. Blundell, I. M. Marshall, T. Lancaster, A. Husmann, C. Steer, W. Hayes, C. Fischmeister, R. E. Martin, and A. B. Holmes, “$\mu$SR in polymers,” Physica B 326, 34 (2003).
  * Nishiyama, Nishiyama, and Higemoto (2003) K. Nishiyama, S. Nishiyama, and W. Higemoto, “Asymmetric F–$\mu$–F interaction of the muon in polyfluorocarbons,” Physica B 326, 41 (2003).
  * Lancaster  _et al._ (2009) T. Lancaster, F. L. Pratt, S. J. Blundell, I. McKenzie, and H. E. Assender, “Muon-fluorine entanglement in fluoropolymers,” J. Phys.: Condens. Matter 21, 346004 (2009).
  * Lancaster  _et al._ (2007) T. Lancaster, S. J. Blundell, P. J. Baker, M. L. Brooks, W. Hayes, F. L. Pratt, J. L. Manson, M. M. Conner, and J. A. Schlueter, “Muon-fluorine entangled states in molecular magnets,” Phys. Rev. Lett. 99, 267601 (2007).
  * Kawaguchi and Hirota (1987) K. Kawaguchi and E. Hirota, “Diode laser spectroscopy of the $\nu_{3}$ and $\nu_{2}$ bands of FHF-1 in 1300 cm-1 region,” J. Chem. Phys. 87, 6838–6841 (1987).
  * Wilkinson and Blundell (2020) J. M. Wilkinson and S. J. Blundell, “Information and decoherence in a muon-fluorine coupled system,” Phys. Rev. Lett. 125, 087201 (2020).
  * Wilkinson  _et al._ (2021) J. M. Wilkinson, F. L. Pratt, T. Lancaster, P. J. Baker, and S. J. Blundell, “Muon sites in ${\mathrm{PbF}}_{2}$ and ${\mathrm{YF}}_{3}$: Decohering environments and the role of anion Frenkel defects,” Phys. Rev. B 104, L220409 (2021).
  * Bonfà  _et al._ (2021) P. Bonfà, M. M. Isah, B. A. Frandsen, E. J. Gibson, E. Brück, I. J. Onuorah, R. De Renzi, and G. Allodi, “Ab initio modeling and experimental investigation of ${\mathrm{Fe}}_{2}\mathrm{P}$ by DFT and spin spectroscopies,” Phys. Rev. Materials 5, 044411 (2021).
  * Williams  _et al._ (2016) R. C. Williams, F. Xiao, I. O. Thomas, S. J. Clark, T. Lancaster, G. A. Cornish, S. J. Blundell, W. Hayes, A. K. Paul, C. Felser, and M. Jansen, “Muon-spin relaxation study of the double perovskite insulators Sr2BOsO6 (B=Fe, Y, In),” Journal of Physics: Condensed Matter 28, 076001 (2016).
  * Lang  _et al._ (2016) F. Lang, P. J. Baker, A. A. Haghighirad, Y. Li, D. Prabhakaran, R. Valentí, and S. J. Blundell, “Unconventional magnetism on a honeycomb lattice in $\alpha$-${\mathrm{RuCl}}_{3}$ studied by muon spin rotation,” Phys. Rev. B 94, 020407 (2016).
  * Kirschner  _et al._ (2019) F. K. K. Kirschner, R. D. Johnson, F. Lang, D. D. Khalyavin, P. Manuel, T. Lancaster, D. Prabhakaran, and S. J. Blundell, “Spin Jahn-Teller antiferromagnetism in ${\mathrm{CoTi}}_{2}{\mathrm{O}}_{5}$,” Phys. Rev. B 99, 064403 (2019).
  * Lang  _et al._ (2019) F. Lang, L. Jowitt, D. Prabhakaran, R. D. Johnson, and S. J. Blundell, “${\mathrm{FeTi}}_{2}{\mathrm{O}}_{5}$: A spin Jahn-Teller transition enhanced by cation substitution,” Phys. Rev. B 100, 094401 (2019).
  * Sugiyama  _et al._ (2019) J. Sugiyama, K. Miwa, H. Nozaki, Y. Kaneko, B. Hitti, D. Arseneau, G. D. Morris, E. J. Ansaldo, and J. H. Brewer, “Magnetic moment of rare-earth elements in ${R}_{2}{\mathrm{Fe}}_{14}\mathrm{B}$ estimated with ${\mu}^{+}\mathrm{SR}$,” Phys. Rev. Mater. 3, 064402 (2019).
  * Amato  _et al._ (2014) A. Amato, P. Dalmas de Réotier, D. Andreica, A. Yaouanc, A. Suter, G. Lapertot, I. M. Pop, E. Morenzoni, P. Bonfà, F. Bernardini, and R. De Renzi, “Understanding the $\mu$SR spectra of MnSi without magnetic polarons,” Phys. Rev. B 89, 184425 (2014).
  * Bonfá , Sartori, and De Renzi (2015) P. Bonfá , F. Sartori, and R. De Renzi, “Efficient and reliable strategy for identifying muon sites based on the double adiabatic approximation,” J. Phys. Chem. C 119, 4278–4285 (2015).
  * Storchak  _et al._ (2011) V. G. Storchak, J. H. Brewer, R. L. Lichti, T. A. Lograsso, and D. L. Schlagel, “Electron localization into spin-polaron state in mnsi,” Phys. Rev. B 83, 140404 (2011).
  * Franke  _et al._ (2018) K. J. A. Franke, B. M. Huddart, T. J. Hicken, F. Xiao, S. J. Blundell, F. L. Pratt, M. Crisanti, J. A. T. Barker, S. J. Clark, A. Štefančič, M. C. Hatnean, G. Balakrishnan, and T. Lancaster, “Magnetic phases of skyrmion-hosting ${\mathrm{GaV}}_{4}{\mathrm{S}}_{8-y}{\mathrm{Se}}_{y}\phantom{\rule{4.0pt}{0.0pt}}(y=0,2,4,8)$ probed with muon spectroscopy,” Phys. Rev. B 98, 054428 (2018).
  * Lancaster (2019) T. Lancaster, “Skyrmions in magnetic materials,” Contemporary Physics 60, 246–261 (2019).
  * Tustain  _et al._ (2020) K. Tustain, B. Ward-O’Brien, F. Bert, T. Han, H. Luetkens, T. Lancaster, B. M. Huddart, P. J. Baker, and L. Clark, “From magnetic order to quantum disorder in the Zn-barlowite series of S = 1/2 kagomé antiferromagnets,” npj Quantum Materials 5, 74 (2020).
  * Van de Walle (2000) C. G. Van de Walle, “Hydrogen as a cause of doping in zinc oxide,” Phys. Rev. Lett. 85, 1012–1015 (2000).
  * Cox (2009) S. F. J. Cox, “Muonium as a model for interstitial hydrogen in the semiconducting and semimetallic elements,” Rep. Prog. Phys. 72, 116501 (2009).
  * Hiraishi  _et al._ (2022) M. Hiraishi, H. Okabe, A. Koda, R. Kadono, and H. Hosono, “Ambipolarity of diluted hydrogen in wide-gap oxides revealed by muon study,” J. Appl. Phys. 132, 105701 (2022).
  * Dehn  _et al._ (2020) M. H. Dehn, J. K. Shenton, S. Holenstein, Q. N. Meier, D. J. Arseneau, D. L. Cortie, B. Hitti, A. C. Y. Fang, W. A. MacFarlane, R. M. L. McFadden, G. D. Morris, Z. Salman, H. Luetkens, N. A. Spaldin, M. Fechner, and R. F. Kiefl, “Observation of a charge-neutral muon-polaron complex in antiferromagnetic ${\mathrm{cr}}_{2}{\mathrm{o}}_{3}$,” Phys. Rev. X 10, 011036 (2020).
  * Dehn  _et al._ (2021) M. H. Dehn, J. K. Shenton, D. J. Arseneau, W. A. MacFarlane, G. D. Morris, A. Maigné, N. A. Spaldin, and R. F. Kiefl, “Local electronic structure and dynamics of muon-polaron complexes in ${\mathrm{fe}}_{2}{\mathrm{o}}_{3}$,” Phys. Rev. Lett. 126, 037202 (2021).
  * Ramadhan  _et al._ (2022) M. R. Ramadhan, B. Adiperdana, I. Ramli, D. P. Sari, A. E. Putri, U. Widyaiswari, H. b. Rozak, W. N. Zaharim, A. Manaf, B. Kurniawan, M. I. Mohamed-Ibrahim, S. Sulaiman, T. Kawamata, T. Adachi, Y. Koike, and I. Watanabe, “Estimation of the on-site Coulomb potential and covalent state in ${\mathrm{La}}_{2}{\mathrm{CuO}}_{4}$ by muon spin rotation and density functional theory calculations,” Phys. Rev. Res. 4, 033044 (2022).
  * De Renzi  _et al._ (2012) R. De Renzi, P. Bonfà, M. Mazzani, S. Sanna, G. Prando, P. Carretta, R. Khasanov, A. Amato, H. Luetkens, M. Bendele, F. Bernardini, S. Massidda, A. Palenzona, M. Tropeano, and M. Vignolo, “Effect of external pressure on the magnetic properties of LnFeAsO (Ln = La, Ce, Pr, Sm),” Superconductor Science and Technology 25, 084009 (2012).
  * Cheung  _et al._ (2018) S. C. Cheung, Z. Guguchia, B. A. Frandsen, Z. Gong, K. Yamakawa, D. E. Almeida, I. J. Onuorah, P. Bonfá, E. Miranda, W. Wang, D. W. Tam, Y. Song, C. Cao, Y. Cai, A. M. Hallas, M. N. Wilson, T. J. S. Munsie, G. Luke, B. Chen, G. Dai, C. Jin, S. Guo, F. Ning, R. M. Fernandes, R. De Renzi, P. Dai, and Y. J. Uemura, “Disentangling superconducting and magnetic orders in ${\mathrm{NaFe}}_{1-x}{\mathrm{Ni}}_{x}\mathrm{As}$ using muon spin rotation,” Phys. Rev. B 97, 224508 (2018).
  * Sundar  _et al._ (2023) S. Sundar, N. Azari, M. R. Goeks, S. Gheidi, M. Abedi, M. Yakovlev, S. R. Dunsiger, J. M. Wilkinson, S. J. Blundell, T. E. Metz, I. M. Hayes, S. R. Saha, S. Lee, A. J. Woods, R. Movshovich, S. M. Thomas, N. P. Butch, P. F. S. Rosa, J. Paglione, and J. E. Sonier, “Ubiquitous spin freezing in the superconducting state of UTe2,” Comm. Phys. 6, 24 (2023).
  * Huddart  _et al._ (2021a) B. M. Huddart, I. J. Onuorah, M. M. Isah, P. Bonfà, S. J. Blundell, S. J. Clark, R. De Renzi, and T. Lancaster, “Intrinsic nature of spontaneous magnetic fields in superconductors with time-reversal symmetry breaking,” Phys. Rev. Lett. 127, 237002 (2021a).
  * Huddart  _et al._ (2021b) B. M. Huddart, M. Gomilšek, T. J. Hicken, F. L. Pratt, S. J. Blundell, P. A. Goddard, S. J. Kaech, J. L. Manson, and T. Lancaster, “Magnetic order and ballistic spin transport in a sine-Gordon spin chain,” Phys. Rev. B 103, L060405 (2021b).
  * Eggert and Affleck (1992) S. Eggert and I. Affleck, “Magnetic impurities in half-integer-spin Heisenberg antiferromagnetic chains,” Phys. Rev. B 46, 10866–10883 (1992).
  * Chakhalian  _et al._ (2003) J. A. Chakhalian, R. F. Kiefl, R. Miller, J. Brewer, S. R. Dunsiger, G. Morris, W. A. MacFarlane, J. E. Sonier, S. Eggert, I. Affleck, A. Keren, and M. Verdaguer, “Local Magnetic Susceptibility of the Positive Muon in the Quasi-One-Dimensional $S=1/2$ Antiferromagnet Dichlorobis (Pyridine) Copper (II),” Phys. Rev. Lett. 91, 027202 (2003).
  * Xiao  _et al._ (2015) F. Xiao, J. S. Möller, T. Lancaster, R. C. Williams, F. L. Pratt, S. J. Blundell, D. Ceresoli, A. M. Barton, and J. L. Manson, “Spin diffusion in the low-dimensional molecular quantum heisenberg antiferromagnet $\mathrm{Cu}(\mathrm{pyz}){({\mathrm{NO}}_{3})}_{2}$ detected with implanted muons,” Phys. Rev. B 91, 144417 (2015).
  * Lancaster  _et al._ (2018) T. Lancaster, F. Xiao, B. M. Huddart, R. C. Williams, F. L. Pratt, S. J. Blundell, S. J. Clark, R. Scheuermann, T. Goko, S. Ward, J. L. Manson, C. Rüegg, and K. W. Krämer, “Quantum magnetism in molecular spin ladders probed with muon-spin spectroscopy,” New Journal of Physics 20, 103002 (2018).
  * Blundell  _et al._ (2013) S. J. Blundell, J. S. Möller, T. Lancaster, P. J. Baker, F. L. Pratt, G. Seber, and P. M. Lahti, “$\mu$SR study of magnetic order in the organic quasi-one-dimensional ferromagnet F4BImNN,” Phys. Rev. B 88, 064423 (2013).
  * Foronda  _et al._ (2015) F. R. Foronda, F. Lang, J. S. Möller, T. Lancaster, A. T. Boothroyd, F. L. Pratt, S. R. Giblin, D. Prabhakaran, and S. J. Blundell, “Anisotropic local modification of crystal field levels in Pr-based pyrochlores: A muon-induced effect modeled using density functional theory,” Physical Review Letters 114, 5–9 (2015).
  * Zhou  _et al._ (2008) H. D. Zhou, C. R. Wiebe, J. A. Janik, L. Balicas, Y. J. Yo, Y. Qiu, J. R. D. Copley, and J. S. Gardner, “Dynamic Spin Ice: ${\mathrm{Pr}}_{2}{\mathrm{Sn}}_{2}{\mathrm{O}}_{7}$,” Phys. Rev. Lett. 101, 227204 (2008).
  * Princep  _et al._ (2013) A. J. Princep, D. Prabhakaran, A. T. Boothroyd, and D. T. Adroja, “Crystal-field states of Pr3+ in the candidate quantum spin ice Pr2Sn2O7,” Phys. Rev. B 88, 104421 (2013).
  * Bleaney (1973) B. Bleaney, “Enhanced nuclear magnetism,” Physica 69, 317–329 (1973).
  * Blundell (2023) S. J. Blundell, “The quantum muon,” J. Phys.: Conf. Ser. 2462, 012001 (2023).
  * Gomilšek  _et al._ (2022) M. Gomilšek, F. L. Pratt, S. P. Cottrell, S. J. Clark, and T. Lancaster, “Many-body quantum muon effects and quadrupolar coupling in solid nitrogen,”  (2022).
  * Mañas-Valero  _et al._ (2021) S. Mañas-Valero, B. M. Huddart, T. Lancaster, E. Coronado, and F. L. Pratt, “Quantum phases and spin liquid properties of 1T-TaS2,” npj Quantum Materials 6, 69 (2021).
  * Onuorah, Bonfà, and De Renzi (2018) I. J. Onuorah, P. Bonfà, and R. De Renzi, “Muon contact hyperfine field in metals: A DFT calculation,” Phys. Rev. B 97, 174414 (2018).
  * Onuorah  _et al._ (2019) I. J. Onuorah, P. Bonfà, R. De Renzi, L. Monacelli, F. Mauri, M. Calandra, and I. Errea, “Quantum effects in muon spin spectroscopy within the stochastic self-consistent harmonic approximation,” Phys. Rev. Mater. 3, 073804 (2019).
  * Gill, Pratt, and Blundell (2023) G. J. W. Gill, F. L. Pratt, and S. J. Blundell, “Extending WiMDA for the data analysis of $\mu^{-}$SR experiments,” J. Phys.: Conf. Ser. 2462, 012010 (2023).