supplement

# How Correlated Adsorbate Dynamics on Realistic Substrates Can Give Rise to
1/$\omega$ Electric-Field Noise in Surface Ion Traps

Benjamin Foulon  Department of Chemical Engineering, Brown University,
Providence, RI 02912, USA  Keith G. Ray  Quantum Simulations Group, Lawrence
Livermore National Laboratory, Livermore, CA 94550, USA  Chang-Eun Kim
Quantum Simulations Group, Lawrence Livermore National Laboratory, Livermore,
CA 94550, USA  Yuan Liu  Center for Ultracold Atoms, Research Laboratory of
Electronics, Massachusetts Institute of Technology, Cambridge, MA 02139, USA
Brenda M. Rubenstein  Department of Chemistry, Brown University, Providence,
RI 02912, USA  Vincenzo Lordi  Quantum Simulations Group, Lawrence Livermore
National Laboratory, Livermore, CA 94550, USA

###### Abstract

Ion traps are promising architectures for implementing scalable quantum
computing, but they suffer from excessive “anomalous” heating that prevents
their full potential from being realized. This heating, which is orders of
magnitude larger than that expected from Johnson-Nyquist noise, results in ion
motion that leads to decoherence and reduced fidelity in quantum logic gates.
The exact origin of anomalous heating is an open question, but experiments
point to adsorbates on trap electrodes as a likely source. Many different
models of anomalous heating have been proposed, but these models have yet to
pinpoint the atomistic origin of the experimentally-observed $1/\omega$
electric field noise scaling observed in ion traps at frequencies between
0.1-10 MHz. In this work, we perform the first computational study of the ion
trap electric field noise produced by the motions of multiple monolayers of
adsorbates described by first principles potentials. In so doing, we show that
correlated adsorbate motions play a definitive role in producing $1/\omega$
noise and identify candidate collective adsorbate motions, including
translational and rotational motions of adsorbate patches and multilayer
exchanges, that give rise to $1/\omega$ scaling at the MHz frequencies
typically employed in ion traps. These results demonstrate that multi-
adsorbate systems, even simple ones, can give rise to a set of activated
motions that can produce the $1/\omega$ noise observed in ion traps and that
collective, rather than individual, adsorbate motions are much more likely to
give rise to low-frequency heating.

##  I Introduction

Despite their outstanding potential, fault-tolerant quantum computers capable
of scaling to many qubits remain beyond our current grasp, in large part due
to noise.Preskill (2018) One particularly pernicious and ubiquitous form of
noise that limits the coherence times of trapped ion,Wineland  _et al._
(1998); Turchette  _et al._ (2000); Brownnutt  _et al._ (2015); Brown  _et
al._ (2021) superconducting qubit, Wang  _et al._ (2015); Martinis (2009); Lee
_et al._ (2014) Rydberg atom,Saffman (2016) nitrogen-vacancy center-based,Kim
_et al._ (2015) and many other quantum architectures is surface noise. In ion
traps – in which cations that possess multiple, addressable hyperfine states
are confined near metallic electrode surfaces – surface-derived electric field
noise has long been known to limit the fidelity of quantum gates by heating
the very ion motions on which these gates depend.Brownnutt  _et al._ (2015)
The exact microscopic source of this electric field noise, however, remains
far less understood. Although intuition would reflexively point to Johnson
noiseJohnson (1928); Nyquist (1928) – the thermal noise generated by currents
that run through the trap electrodes and circuitry – as the dominant source of
noise, a number of researchers have independently shown that a different,
still unknown and therefore “anomalous” source of noise dominates the Johnson
noise in ion traps.Hite  _et al._ (2013a); Brown  _et al._ (2021) The origin
of this noise thus remains an open question in ion trap science whose
resolution would dramatically improve the performance of ion trap-based
quantum information processors,Haffner  _et al._ (2008); Benhelm  _et al._
(2008); Bruzewicz  _et al._ (2019) quantum sensors,McCormick  _et al._ (2019);
Degen  _et al._ (2017); Baumgart  _et al._ (2016); Leibfried  _et al._ (2004)
and clocks.Kozlov  _et al._ (2018); Chou  _et al._ (2010); Huntemann  _et al._
(2016)

Experiments directed at characterizing noise have revealed that the three key
parameters that control the electric field noise, $S_{E}$, are the distance at
which the ion is trapped above the electrode surface, $d$, the trapping
frequency, $\omega$, and the trap temperature, $T$, such that

$S_{E}\propto\omega^{-\alpha}d^{-\beta}f(T),$ (1)

where $f(T)$ is a function of the temperature.Brownnutt  _et al._ (2015)
Although aspects of this functional form and its exponents are still a subject
of debate,An  _et al._ (2019); Brownnutt  _et al._ (2015) most experiments
point to a frequency scaling involving an $\alpha$ of between 0 and 2 at MHz
frequencies.Turchette  _et al._ (2000); Daniilidis  _et al._ (2011); Sedlacek
_et al._ (2018a, b); Noel  _et al._ (2019); Labaziewicz  _et al._ (2008) This
scaling is in line with the practical experience that most sources of noise
have a “$1/\omega$” character with $\alpha=1$. Many of these same experiments
also suggest a distance scaling of $d^{-4}$ with $\beta=4$,Turchette  _et al._
(2000); Sedlacek  _et al._ (2018a, b); Deslauriers  _et al._ (2006);
Daniilidis  _et al._ (2011); Hite  _et al._ (2017); Labaziewicz  _et al._
(2008) which stands in stark contrast with Johnson noise, which scales as
$d^{-2}$.

Over the past decade, a number of different microscopic modelsNoel  _et al._
(2019); Kumph  _et al._ (2016); Teller  _et al._ (2021) have been advanced
that aim to both reproduce and explain these scalings, which serve as useful
constraints on the possible mechanisms that could give rise to anomalous
heating. One of the earliest models proposed was the patch potential
model,Turchette  _et al._ (2000) which espouses that local variations of
electrode potentials can induce ion motion with the $d^{-4}$ noise scaling
observed in experiments. The patch potential model, however, does not identify
the source of the local fluctuations on which it is based. In light of
experiments demonstrating that different treatments that remove surface
adsorbates, including ion milling,Hite  _et al._ (2012); Daniilidis  _et al._
(2014); McKay  _et al._ (2014) plasma treatment,McConnell  _et al._ (2015) and
laser cleaning,Allcock  _et al._ (2011) reduced trap noise by up to two orders
of magnitude, it is logical to attribute this noise to dipole fluctuations
caused by adsorbates bound to trap electrodes.

Along these lines, two particularly compelling noise models are the adatom
dipole and diffusion models. The adatom dipole model posits that atoms and
molecules that adsorb onto electrode surfaces to form layers or patches
develop induced dipole moments that locally lower the work function, and
therefore the potential, of the metal.Safavi-Naini  _et al._ (2011); Brownnutt
_et al._ (2015) The dipoles may then be caused to fluctuate by phonon-induced
transitions among different vibrational states of the adatom-surface
potential, giving rise to the requisite fluctuation spectra. Initial
investigations of this model that treated single adsorbates using plausible
values for the their masses and forms for the their binding potentials
illustrated that it can yield electric field fluctuation spectra that are of
the same magnitude as those observed in experiments.Safavi-Naini  _et al._
(2011, 2013) These spectra were flat as a function of frequency at low
frequencies (the ‘white noise’ regime), yet transitioned to decaying as
$1/\omega^{2}$ at high frequencies with a brief $1/\omega$ regime in between.
Such behavior is consistent with two-level fluctuator models prevalent in the
signal processing literature.Dutta and Horn (1981); Zanolla (2009); Weissman
(1988) More recent work employing a first principles treatment of single
hydrocarbon, water, and other adsorbate binding potentials and dipole moments
similarly exhibited a narrow $1/\omega$ region in its spectra.Ray  _et al._
(2019) That work additionally included in-plane adsorbate vibrations, which
produced the highest electric field noise and exhibited $1/\omega$ scaling at
low frequencies orders of magnitude larger than the typical 1 MHz trap
frequency. Interestingly, work released during the submission of this
manuscript has suggested that phonons that correlate the motions of multiple
adsorbates may also produce $1/\omega$ frequency scalings.Lloyd  _et al._
(2021) In the adatom diffusion model, dipole fluctuations are instead viewed
as originating from the diffusion of an adataom or adsorbate across the
surface of an electrode.Kim  _et al._ (2017) Noninteracting adatoms diffusing
over a surface covered with patches with differing work functions may be shown
to produce fluctuation spectra that scale as $1/\omega^{3/2}$.Hite  _et al._
(2013b); Gesley and Swanson (1985) Even though these models do not produce
extended regions with $1/\omega$ frequency scalings in their current forms,
the inclusion of realistic features such as adatom-adatom interactions and
surface corrugation may lead to more complicated spectral features that have
yet to be illuminated or understood. This thus raises the question of whether
a more realistic treatment of adatom dipole fluctuations that incorporates
diffusion and interparticle interactions may exhibit the so far theoretically
elusive $1/\omega$ scaling expected from experiments.

In this work, we examine the frequency-dependence of the electric field
fluctuation spectra produced by the dynamics of multiple monolayers of
interacting methane adsorbates on gold, a commonly-used electrode material in
ion traps. To study the adsorbate dynamics, we modeled them using classical
Molecular Dynamics (MD) on potential energy surfaces (PES) constructed based
on first principles Density Functional Theory (DFT) calculations. From the MD
trajectories, we obtained coverage- and temperature-dependent dipole-dipole
fluctuation spectra using DFT-derived dipole moments. This enabled us to make
a clear connection between the features of the fluctuation spectra at
different frequencies and the underlying adsorbate motions that give rise to
these features.

In doing so, we demonstrate that correlated dipole fluctuations are
overwhelmingly responsible for producing “$1/\omega$” noise and identify
several key microscopic motions that naturally lead to a “$1/\omega$”
frequency scaling in the MHz regime used in most ion traps. In particular, we
find that correlated rotational and translational motions of adsorbates within
clusters can give rise to “$1/\omega$” electric field noise at submonolayer
coverages, while interlayer particle exchanges among the first two layers of
adsorbates are the largest contributors to such noise at supermonolayer
coverages. We furthermore show how many of the MHz features of the electric
field spectra for these systems may be reproduced by a two-level fluctuator
model with these motions as rare events. Even though the methane-gold system
we study here is a model system that contains simplifications, our work
provides clear evidence for the types of adsorbate motions that can give rise
to anomalous heating and demonstrates the crucial role that collective, rather
than individual, adsorbate motions assume in the noise generation process at
experimental trap frequencies.

This paper is organized as follows. In Section II, we begin by describing how
we model the dipole-dipole fluctuation spectra of methane on a gold substrate
using our combined DFT-MD approach. We next describe our results, including
our key findings regarding the $1/\omega$-$1/\omega^{2}$ frequency-dependence
of the dipole-dipole fluctuation spectra we obtain for varying adsorbate
surface coverages and temperatures in Section III. In the same section, we
additionally present our data directly linking specific adsorbate motions with
features of the frequency spectra and demonstrate how our spectra can be
reproduced using simple two-state models. Lastly, in Section IV, we place our
findings that a realistic model of multilayer adsorbate dynamics can give rise
to $1/\omega$ noise in the context of the ongoing search for the microscopic
origins of anomalous heating, and discuss the limitations and natural
extensions of our current model. Additional information regarding our
simulations and interpretation of the data may be found in the Supplemental
Information.

##  II Methods

In this work, we model methane adsorbate dynamics on gold substrates at a
range of temperatures and surface coverages by running classical MD
trajectories on a DFT-derived substrate potential energy surface (PES). Based
upon the trajectories obtained, we then compute adsorbate dipole-dipole
correlation functions and Fourier transform them to acquire dipole-dipole
fluctuation spectra. These spectra are subsequently analyzed for their
frequency-dependent behavior, which can be used to predict ion trap noise.

###  II.1 Calculation of Heating Rates

As has been shown in previous work,Turchette  _et al._ (2000) the coupling of
electric field fluctuations with the motion of trapped ions gives rise to a
heating rate, $\dot{\hat{n}}$, that is given by:

$\dot{\hat{n}}=\frac{q^{2}}{4m_{I}\hbar\omega_{t}}S_{E}(\omega_{t}).$ (2)

Here, $q$ is the charge of the ion, $m_{I}$ is the mass of the ion,
$\omega_{t}$ is the frequency at which the ion is trapped (typically 0.1-10
MHz), and $\hbar$ is the reduced Planck’s constant. $S_{E}(\omega)$ is the
frequency spectrum of the electric field fluctuations. Patch potential models
assert that electric field noise in traps stems from local potential
fluctuations above the surfaces of the electrodes.Brownnutt  _et al._ (2015);
Dubessy  _et al._ (2009); Rossi and Opat (1992) These fluctuations can arise
from regions of the electrodes with varying crystal orientations or adsorbate
surface motion. In the latter case, which is the focus of this investigation,
fluctuations in the dipole moment of the adsorbates caused by electronic
interactions with the electrode surface give rise to electric field noise.
Specifically, for a conventional planar trap,Brownnutt  _et al._ (2015)
$S_{E}(\omega)$ may be obtained from out-of-plane adsorbate dipole
fluctuations in a surface patch by:Ray  _et al._ (2019)

$S_{E}(\omega)=\frac{3\pi\sigma
S_{\mu}(\omega)}{2(4\pi\epsilon_{0})^{2}d^{4}}$ (3)

and from in-plane fluctuations by:

$S_{E}(\omega)=\frac{\pi\sigma S_{\mu}(\omega)}{(4\pi\epsilon_{0})^{2}d^{4}}.$
(4)

In Equations (3) and (4), $\sigma$ is the average area density of the surface
patches, $d$ is the ion-electrode distance, $\epsilon_{0}$ is the permittivity
of free space, and $S_{\mu}(\omega)$ is the dipole-dipole fluctuation
spectrum. Notably, these expressions for $S_{E}(\omega)$ reflect the $d^{-4}$
scaling that has often been observed experimentally in heating rates and is
one key reason why we and others continue to examine this model.

The dipole-dipole fluctuation spectrum of a surface patch can be calculated by
taking the Fourier Transform of the dipole-dipole autocorrelation function:

$S_{\mu}(\omega)=\int_{-\infty}^{+\infty}d\tau
e^{i\omega\tau}C_{\mu,\mu}(\tau),$ (5)

where $C_{\mu,\mu}(\tau)$ represents the dipole-dipole autocorrelation
function of the total patch dipole moment at time $\tau$. This autocorrelation
function may be expressed as:

$\displaystyle C_{\mu,\mu}(\tau)$ $\displaystyle=$
$\displaystyle\langle(\mu_{z}(\tau)-\langle\mu_{z}\rangle)(\mu_{z}(0)-\langle\mu_{z}\rangle)\rangle$
(6) $\displaystyle=$
$\displaystyle\langle\Delta\mu_{z}(\tau)\Delta\mu_{z}(0)\rangle.$

In the above, $\mu_{z}(\tau)$ is the $z$-component of the total dipole moment
of a patch at time $\tau$ and $\langle\mu_{z}\rangle$ is its equilibrium
average. Although a full treatment would use $\vec{\mu}$, we instead use
$\mu_{z}$ unless otherwise specified. The values of the $x$ and $y$ components
of our adsorbate dipoles are orders of magnitude smaller than the $z$
components. There are only a small handful of cases where the $x$ or $y$
components would affect spectra generated by $z$ components alone (see
Supplemental Information).

In this work, we assume that $C_{\mu,\mu}(\tau)$ stems from how the dipole
moments of $N$ different adsorbates change as they move across a surface with
a spatially-varying binding potential. The quantity $\Delta\mu_{z}(\tau)$
represents a fluctuation of the electric dipole moment of the simulation cell,
which is a sum over the individual adsorbate dipole fluctuations. This
definition remains valid for Equations (3) and (4) as long as the simulation
cell is small compared to the ion-electrode distance, $d$, and it captures the
effect of adsorbate correlations on the noise as long as the simulation cell
is larger than the correlation length.Brownnutt  _et al._ (2015) The first
condition is satisfied as the typical ion-electrode distance is 40
$\mathrm{\mu}$m and the simulation cell dimensions are less than 100 nm. The
adsorbate-adsorbate correlation length is shown to be smaller than the
simulation size in Section III.2.

In this limit, we may express the surface patch dipole-dipole correlation
function in terms of the adsorbate dipole fluctuations as:

$\displaystyle C_{\mu,\mu}(\tau)$ $\displaystyle=$
$\displaystyle\Big{\langle}\sum_{i}^{N}\Delta\mu_{i,z}(\tau)\sum_{j}^{N}\Delta\mu_{j,z}(0)\Big{\rangle}$
$\displaystyle=$ $\displaystyle\Big{\langle}\sum_{\begin{subarray}{c}i,j\\\
i=j\end{subarray}}^{N,N}\Delta\mu_{i,z}(\tau)\Delta\mu_{j,z}(0)\Big{\rangle}$
$\displaystyle+$ $\displaystyle\Big{\langle}\sum_{\begin{subarray}{c}i,j\\\
i\neq
j\end{subarray}}^{N,N}\Delta\mu_{i,z}(\tau)\Delta\mu_{j,z}(0)\Big{\rangle},$

where $\Delta\mu_{i,z}(\tau)$ represents the fluctuation of the $z$ component
of the electric dipole moment of adsorbate $i$ at time $\tau$. In the last
expression, we have divided the sum into a sum over dipole-dipole correlation
functions of the individual ($i=j$) adsorbates and a sum over those
correlations between distinct ($i\neq j$) adsorbates. If there were no
correlation between different adsorbates, the second term would vanish. We
calculate both of these sums and demonstrate the importance of dipole-dipole
correlations between different adsorbates in Section III.2.

In practice, we evaluate Equation (6) by averaging over all of the individual
adsorbate dipole-dipole autocorrelation functions taken across the entire
simulation runtime, $\mathcal{T}$:

$\quad C_{\mu,\mu}(\tau)=\\\
\frac{1}{(\mathcal{T}-\tau)}\sum_{k=1}^{\mathcal{T}-\tau}\left[\sum_{i}^{N}\Delta\mu_{i,z}(\tau+t_{k})\sum_{j}^{N}\Delta\mu_{j,z}(t_{k})\right],$
(8)

where $\mathcal{T}-\tau$ denotes the number of timesteps taken between these
two times.

###  II.2 Generation of the Electrode Potential Energy Surface

In order to calculate the fluctuations of the dipole moments required by
Equation (8), a PES of CH4 physisorbed onto Au(111) was interpolated from DFT
CH4-Au(111) binding energies taken at a representative set of positions on the
Au(111) surface. In all of our simulations, the methane adsorbates were
treated as point particles, a reasonable simplification given their spherical
character, and we did not consider the effects of phonons.

The DFT calculations were performed with VASPKresse and Furthmüller (1996)
using the vdW-DF-cx functional Berland and Hyldgaard (2014) at a single
k-point and a plane-wave cutoff of 600 eV. The resulting potential energy
points were interpolated across a periodic surface. The interpolated PES
employed has dimensions of 2.856 Å, 4.947 Å, and 29.0 Å in the $x$, $y$, and
$z$ directions, respectively. Figure 1a depicts this PES in the $(x,y)$-plane
for one Au(111) unit cell at the minimum-energy distance from the carbon of
the adsorbate to the surface of 3.16 Å. The simulation cell contained 72 such
unit cells (twelve in the $x$ direction and six in the $y$ direction). This
unit cell tiling was chosen so as to make the dimensions of the overall
simulation cell roughly the same in both the $x$ and $y$ directions to reduce
any spurious adsorbate ordering effects from anisotropic boundary conditions.
As can be seen in Figure 1a, the hollow site marks the lowest-energy binding
position on this surface (-0.1625 eV), while the atop site marks the highest-
energy position (-0.1575 eV). Figure S3 projects the PES along the $z$
direction. Additional simulation parameters may be found in Table 1, and
additional binding energies at different positions on the surface may be found
in Table SI.

Figure 1: (a) Visualization of the CH4-Au(111) potential energy surface
generated from DFT calculations spanning one Au(111) unit cell at an
adsorbate-surface distance of 3.16 Å, which is the minimum-energy, out-of-
plane distance. The red, orange, and pink circles designate the positions of
the Au(111) hollow, atop, and bridged sites, respectively. (b) Visualization
of the magnitudes of the surface-induced methane dipole moments as a function
of methane position in the same Au(111) cell and at the same height as on the
left.

In addition to the surface-adsorbate interaction, each adsorbate was modeled
as interacting with its neighbors through the methane-methane Optimized
Potentials for Liquid Simulations (OPLS) potentialJorgensen  _et al._ (1984)
such that:

$V_{adsorbate}=V_{surface}+V_{interparticle},$ (9)

where

$V_{interparticle}=4\epsilon_{CH_{4}}\left[\left(\frac{\sigma_{CH_{4}}}{r_{ij}}\right)^{12}-\left(\frac{\sigma_{CH_{4}}}{r_{ij}}\right)^{6}\right].$
(10)

In this equation, $\epsilon_{CH_{4}}$ denotes the OPLS CH4-CH4 interaction
energy, while $\sigma_{CH_{4}}$ denotes the CH4-CH4 interaction radius.
$r_{ij}$ represents the interparticle distance between the centers of mass of
methanes $i$ and $j$. Based upon the OPLS interaction radius, the methanes
would prefer to be spaced 4.19 Å from one another. As we will see, this
spacing favored by the OPLS potential competes with the particle spacing
favored by the surface potential to produce different particle surface
configurations at different temperatures, which in turn give rise to the
different types of particle motions that are reflected in the dipole spectra.

###  II.3 Molecular Dynamics Using the Atomic Simulation Environment

The Atomic Simulation Environment (ASE) was used to perform the MD simulations
with the interpolated CH4-Au(111) and OPLS potentials. Running MD simulations
based upon this interpolated landscape is much more computationally efficient
than running fully ab initio simulations and is therefore key to reaching the
timescales required to observe events happening in the 0.1-10 MHz frequency
range at which ions are trapped.

ASE is an open-source software platform for atomistic simulations.ase (2017)
It is also Python-based, which makes it possible to seamlessly incorporate the
Python-interpolated PES directly into the MD setup. The MD scripts used can be
found in this paper’s Github repository.Git

NVT MD simulations were performed using Langevin dynamics with a friction
coefficient of 0.01 atomic units (a.u., where 1 a.u. $=4.13\times 10^{16}$
s-1). Particle forces and energies were updated every 5 fs, and positions were
recorded for analysis every 1 ps to ensure that configurations were not
artificially correlated. These parameters were verified to be computationally
efficient while providing comparable system equilibration and dynamics to
those obtained using smaller friction coefficients and step sizes (see Figures
S1 and S2). In order to isolate collective motions of interest from
thermostat-induced cluster translations and drifts, the thermostat was
configured to fix the center-of-mass position and zero the center-of-mass
momentum. Simulation runtimes ranged from 6 to 19 $\mathrm{\mu}$s for
different sets of coverages (see Supplementary Info). The MD cell employed
periodic boundary conditions in the $x$ and $y$ directions, but not in the $z$
direction.

CH4-Au(111) Potential | DFT, vdW-DF-cx
---|---
CH4-Au(111) min. energy distance | 3.16 Å
CH4-Au(111) min. energy | -0.1625 eV
OPLS min. energy distance | 4.19 Å
$\sigma_{\ch{CH4}}$ for OPLS | 3.73 Å
$\epsilon_{\ch{CH4}}$ for OPLS | 12.75 meV
MD Environment | ASE
Thermostat | Langevin
Friction Coefficient | 0.01 a.u.
MD Force Update Step | 5 fs
MD Position Recording Step | 1 ps
Table 1: Simulation details and parameters.

###  II.4 Calculation of Dipole Moments

The electric dipole moments of the methanes adsorbed over a range of positions
on the gold surface were also computed with DFT using VASP. In particular, we
integrated the charge density $\times$ position, $\rho(r)\times r$, over a
40$\times$33.6$\times$29.1 Å3 simulation cell containing an Au slab “island”
with a finite extent in the in-plane directions and a single CH4 molecule on
top. A finite slab is required so that, if the adsorbate is moved within one
Au surface unit, charge is not induced to move across the periodic boundaries
and the integral remains well-defined. The slab island consists of a
4$\times$4 arrangement of Au surface primitive cells. To map the surface-
induced dipole moments as depicted in Figure 1b, we calculated the dipole
moments for CH4 at the surface positions illustrated in one surface unit cell
located in the center of this slab and at varying heights. Electric dipole
corrections were utilized in VASP and a small linear ramp in the dipole moment
was subtracted to remove remaining finite-size effects and ensure that the
dipole moment is periodic with respect to the surface in-plane periodicity. As
in Section II.2, these calculations were also performed using the vdW-DF-cx
functional at a single k-point and with a 600 eV plane-wave cutoff.

The resulting dipole moments were interpolated using an interpolation script
to yield the $x$, $y$, and $z$ components of the dipole vector as a function
of adsorbate position. The values of the $x$ and $y$ components of the dipole
moments are negligible compared to the $z$ component values across the entire
periodic cell, and their fluctuations do not produce much noise except in a
handful of cases. Thus, we chose to make use of the $z$ components of the
dipole moments for the rest of this study (see the Supplemental Information
for supporting data). A visualization of the $z$ dipole moment surface at the
minimum-energy out-of-plane surface-adsorbate distance can be found in Figure
1b. Comparing the dipole moment surface with the PES, it can be seen that the
largest induced dipole moments occur at the locations with the largest
potential energies, and vice versa. The dipole interpolation script that
produced this surface is provided in this paper’s Github repository.Git

###  II.5 Dipole-Dipole Fluctuation Spectra

Dipole-dipole autocorrelation functions were computed using Equation (8) based
upon how the adsorbate dipole moments change in time as they traverse the
potential surface. As we modeled systems with varying surface coverages (see
Table 2), the exact value of $N$ used in this equation depended upon the
coverage studied. The first 1,000 timesteps (spanning 1 ns) were omitted to
allow for equilibration. Dipole-dipole fluctuation spectra were then
calculated by applying a Discrete Fast Fourier TransformPress and Rybicki
(1989) to these correlation functions (per Equation (5)). A standard smoothing
technique, Blackman smoothing,Blackman and Tukey (1958) was applied before
plotting.

###  II.6 Two-Level Fluctuators as a Source of $1/\omega^{2}$ and $1/\omega$
Noise

As will be discussed further below, our electric field fluctuation spectra
exhibit several key frequency scalings. Many of these frequency scalings have
well-known connections to certain types of fluctuations in the signal
processing literature that also arise in our simulations. For instance, a
quantity that fluctuates randomly about a single average value with no
correlation in time leads to featureless (white) noiseShynk (2012) in the
frequency domain. On the other hand, a quantity that transitions between two
discrete values produces a region of $1/\omega^{2}$ frequency scaling; systems
that fluctuate in this manner are known as two-level fluctuatorsDutta and Horn
(1981) (TLFs) or random telegraph signals (RTS).Zanolla (2009)

More specifically, a quantity that switches between states 0 and 1 differing
in value by $\Delta I$ and that has mean state residence times of $\tau_{0}$
and $\tau_{1}$ will yield noise spectra of the following Lorentzian
form:Machlup (1954)

$S(\omega)=\frac{4(\Delta
I)^{2}}{(\tau_{0}+\tau_{1})[(1/\tau_{0}+1/\tau_{1})^{2}+\omega^{2}]}.$ (11)

When $\tau_{0}=\tau_{1}=\tau$, this expression simplifies to:

$S(\omega)=\frac{4(\Delta I)^{2}\tau}{1+\tau^{2}\omega^{2}}.$ (12)

Equations (11) and (12) produce white noise at low frequencies followed by a
rounding-off to a $1/\omega^{2}$ descent at frequencies greater than a cutoff
frequency, given by $\omega_{c}=(1/\tau_{0}+1/\tau_{1})$. Figure S16 shows how
TLF noise spectra vary with $\tau_{0}$ if $\tau_{1}$ is held constant. When
$\tau_{1}$ is constant, the TLF noise reaches a maximum when
$\tau_{0}=\tau_{1}$. When $\tau_{0}>\tau_{1}$, the TLF noise decreases while
$\omega_{c}$ remains mostly unchanged; when $\tau_{0}<\tau_{1}$, the TLF noise
decreases and $\omega_{c}$ increases, effectively shifting the frequency curve
to the right.Zanolla (2009) In general, decreasing $\tau_{0}$ and $\tau_{1}$
will shift the $1/\omega^{2}$ region to higher frequencies.

In addition to a $1/\omega^{2}$ scaling, TLFs can give rise to $1/\omega$
regions in several ways. Every individual TLF has sections of $1/\omega$ and
near-$1/\omega$ scaling in the rounding-off region between the white noise and
$1/\omega^{2}$ regimes. Moreover, an ensemble of TLFs with different
transition frequencies can produce an aggregate $1/\omega$ region over several
decades of frequencies;Zanolla (2009) fluctuators with more than two levels
can also combine together to produce $1/\omega$ regions in a similar manner.

##  III Results and Discussion

###  III.1 Adsorbate Surface Geometries

Figure 2: Representative snapshots of the positions of the adsorbates on the
Au(111) surface at 20 K for (a) 0.36 ML and (b) 1.86 ML coverages. Note that
the intermolecular potential is stronger than the in-plane Au(111)-CH4
potential barriers when the particles are near the surface. Because the
adsorbates cannot simultaneously reside in the differing minima of their
intermolecular and surface-adsorbate potentials, they do not all occupy lowest
energy hollow sites on the surface; instead, they organize into a more
complicated, but repeating pattern. First layer particles can have a maximum
of six first layer and three second layer neighbors.

We begin by analyzing the geometries that the adsorbates assume at the surface
coverages studied. Because of the competing influences of the interparticle
and surface-adsorbate potentials, the geometries are not straightforward to
predict. Here, we visualize and discuss several emblematic coverages.

Within the simulation box employed in this work, we find that 69 adsorbates
produce the full monolayer (ML) surface coverage with the lowest energy per
adsorbate (see Figure S10). We thus define our runs in terms of $(N/69)$ ML
coverages. In Figure 2, for example, the 0.36 ML coverage on the left has 25
particles, whereas the 1.86 ML coverage on the right has 128 particles. As a
guide to our subsequent discussion, we list all of the surface coverages and
temperatures at which our simulations were performed in Table 2.

Type | Coverage ($N_{\ch{CH4}}$)
---|---
Submonolayer | 0.36 ML (25)
Submonolayer | 0.52 ML (36)
Submonolayer | 0.71 ML (49)
Monolayer | 1.0 ML (69)
Supermonolayer | 1.16 ML (80)
Supermonolayer | 1.42 ML (98)
Supermonolayer | 1.86 ML (128)

Temp. $<$ 70 K | 20, 25, 30, 35, 40, 45, 50, 55, 60
---|---
Temp. $\geq$ 70 K | 70, 80, 90, 100, 110, 125, 150, 184, 220

Table 2: Table of coverages and temperatures at which simulations were
performed in this work.

Figure 2a depicts a representative surface geometry of the 0.36 ML trajectory
at 20 K. At this coverage, the methanes cluster together and simultaneously
try to minimize their interparticle energies while residing in the hollow-site
minima of the Au(111) surface. Notably, the methanes cannot reside at every
available hollow site because such sites are separated by only 2.856 Å, which
is significantly less than $\sigma_{\ch{CH4}}$=3.73 Å, the zero-energy OPLS
distance. Since the OPLS well depth is greater than any surface site energy
differences (see Tables 1 and SI), the methanes occasionally settle for
bridged or even atop sites in order to minimize interparticle repulsions. In
so doing, they maintain a CH4-CH4 distance of $>$3.73 Å, close to the
interparticle potential minimum-energy distance of 4.19 Å. Habitually, the
methanes lying near atop sites get pulled toward hollow sites and push
neighboring hollow-site methanes into less favorable sites, and so on and so
forth throughout the clusters. This in-plane rattling occurs often – within
1-2 Å and at much higher frequencies than trapped ion frequencies – but does
not disrupt the cohesion of the clusters at low temperatures. The methanes can
also become perturbed enough to undergo more significant individual and
collective movements, and while these significantly affect the noise spectra
(see Section III.3), they do not prevent the surface clusters from remaining
intact and mostly stationary at low temperatures. At higher temperatures, the
particles have plenty of thermal energy to leave both the surface and CH4-CH4
interaction wells, resulting in a more fluid and less-defined surface
geometry.

Although most previous works have focused on adsorbate dynamics at
submonolayer coverages assuming that diffusion across the electrode surface is
the dominant source of noise,Kim  _et al._ (2017) here we additionally examine
supermonolayer adsorbate dynamics. Studying multiple monolayers is more true
to experimental reality, as several experiments have shown multiple
hydrocarbon monolayers to be present on untreated gold traps.Hite  _et al._
(2012); Kim  _et al._ (2017) Figure 2b depicts a representative 1.86 ML
bilayer coverage at 20 K. Not surprisingly, the second layer stacks on top of
the first such that each second layer particle is roughly equidistant (in the
$(x,y)$ plane) from its three nearest first layer neighbors, while typically
staying at least 3.4-3.6 Å above the first layer (6.56-6.76 Å above the
surface) in the $z$ direction. This leads to clustering in the second layer,
similar to how clusters formed on the surface at submonolayer coverages.

Second layer clusters are affected by interparticle interactions from both
first and second layer particles, but second layer particles are not as
strongly affected by the surface potential. In fact, at typical second layer
distances from the surface, the magnitude of the adsorbate-surface potential
is comparable to the minimum OPLS energy (see Figure S3). Beginning around 70
K, the second layer clusters become less coherent as interparticle
interactions become less restrictive to particle motion. This leads to more
fluid second layer particles, both parallel and perpendicular to the surface,
some of which begin to stochastically exchange layers. Such exchanges can be
reciprocal or unidirectional; in the latter case, the population of the
surface changes. Such exchanges can be reciprocal or unidirectional; in the
latter case, the population of the surface changes. Although 69-adsorbate
surfaces are slightly favored in monolayer conditions, other surface
populations can manifest in supermonolayer coverages, including 68- and
70-adsorbate surfaces. The switching between different surface populations
will become important in our interpretation of the dipole-dipole fluctuation
spectra.

Due to the dwindling influence of the surface-electrode potential at further
distances from the electrode surface, the propensity to form a cohesive third
layer is significantly reduced. In addition, the dwindling electrode-adsorbate
interaction strength makes it more difficult to keep weakly-held second layer
particles from getting pushed up or drifting to higher surface heights as the
second layer gets more crowded. For these reasons, we limit our focus to a
maximum coverage of 1.86 ML. More visualizations of and information about
different coverages can be found in the Supplementary Information.

###  III.2 Effects of Correlated Adsorbate Motion on Fluctuation Spectra

Figure 3: (a) Patch total dipole fluctuation spectrum,
($\mathcal{F}(\langle\sum\Delta\mu_{i,z}(\tau)\\\
\Delta\mu_{j,z}(0)\rangle)(\omega)$), in blue, and that spectrum decomposed
into individual, orange ($i=j$), and distinct, dashed green ($i\neq j$)
adsorbate contributions, which significantly cancel each other. (b) Dipole-
dipole fluctuation spectrum of an individual adsorbate, black line, and the
spectra of pairs of adsorbates separated by less than 6.0 Å, mostly blue, and
between 6.0 Å to 9.0 Å, mostly red, corresponding to surface adsorbate nearest
neighbors and next-nearest neighbors. The colors of the nearest-neighbor and
next-nearest neighbor pair spectra are indicative of their relative phase
compared to the individual adsorbate fluctuations, with $\cos(\phi)=1(-1)$
being dark red (blue). (c)
$\cos(\phi)\mathcal{F}(\langle\Delta\mu_{i,z}(\tau)\Delta\mu_{j,z}(0)\rangle)(\omega=\text{3.5
MHz})$ vs. distance between adsorbates $i$ and $j$. $\phi$ is the phase
difference between the distinct $i/j$ adsorbate spectrum and the individual
adsorbate spectrum at that frequency.

To demonstrate the role that dipole-dipole correlations between distinct
surface adsorbates play in the generation of the electric field noise
affecting trapped ions, we decompose the total simulation patch dipole-dipole
correlation spectrum into contributions from individual adsorbates and from
distinct adsorbate pairs. This decomposition is expressed in Equation (II.1)
and presented in Figure 3 for the case of a 0.52 ML coverage of adsorbed CH4
at 30 K. The trends shown are representative of those also observed for
different coverages and temperatures. We notice in Figure 3a that, while the
individual adsorbate dipole-dipole spectrum and the distinct adsorbate dipole-
dipole spectrum both exhibit roughly a white-noise spectrum up to 10 MHz, the
total dipole-dipole correlation spectrum scales as $1/\omega$ beyond 10 MHz.
Therefore, capturing dipole-dipole correlations between distinct adsorbates is
essential for reproducing the $1/\omega$ electric field noise measured
experimentally at frequencies of a few MHz.

Interestingly, the magnitude of the total dipole-dipole spectrum is smaller
than that of either the individual or distinct pair adsorbate contributions.
This occurs because the pair adsorbate dipole spectrum is close to the
individual adsorbate dipole spectrum in magnitude, but is roughly $\pi$ phase
shifted and thus partially cancels the noise from the individual adsorbate
dipole spectrum. In Figure 3b, we investigate this behavior in more detail by
plotting the average spectrum of an individual adsorbate along with the dipole
spectra between adsorbates that are nearest and next-nearest neighbors,
corresponding to adsorbate separations of less than 6.0 Å and 6.0-9.0 Å,
respectively. Positions were averaged over a window of 2000 timesteps prior to
spectral calculation. The colors of the plotted pair spectra indicate the
phase as $\cos(\phi)$, with dark blue signifying $\cos(\phi)=-1$, $\phi=\pi$
and dark red indicating $\cos(\phi)=1$, $\phi=0$. Both the magnitude of the
spectrum and its color are indicative of the correlation with the individual
adsorbate dipole fluctuation. Nearest-neighbor dipoles are slightly
anticorrelated on average, indicated by the dark blue color and a magnitude
that is smaller than that of the individual adsorbate spectrum, while next-
nearest neighbors have an even smaller correlation that is positive on
average. While the magnitude of the dipole-dipole spectrum of an average
nearest-neighbor pair is less than that of the individual adsorbate spectrum,
there are 5.3 nearest neighbors on average for this coverage at this
temperature, so the sum of all distinct pairs produces the dashed curve in
Figure 3a.

In Figure 3c, we plot
$\cos(\phi)\mathcal{F}(\langle\Delta\mu_{i,z}(\tau)\Delta\mu_{j,z}(0)\rangle)(\omega=\text{3.5
MHz})$ vs. the distance between adsorbates $i$ and $j$, where $\phi$ is the
phase difference between the distinct $i/j$ adsorbate spectrum and the
individual adsorbate spectrum at that frequency and separation. From this
plot, we can see how the magnitude and phase of the pair dipole-dipole
fluctuation spectra evolve with separation. As above, $\cos(\phi)=-1$
signifies anticorrelated adsorbate dipole fluctuations and $\cos(\phi)=1$
signifies correlated adsorbate dipole fluctuations, and adsorbate positions
were averaged over a window of 2000 timesteps. In Figure 3c, we again see the
anticorrelation of the nearest-neighbor adsorbates, which is strongest at 2.5
Å (see also the blue line in Figure 3b). We note, however, that most nearest-
neighbor adsorbates are close to 4.5 Å apart, a distance at which the
anticorrelation is weaker than at 2.5 Å. We can quantify the correlation
between adsorbates as close as an angstrom or two apart because we average
over many adsorbate configurations and a few such rare configurations arise
given enough sampling. The correlation of the next-nearest neighbors, first
captured by the red line in Figure 3b, is positive and most next-nearest
neighbor adsorbates are around 8 Å apart. We estimate the correlation length
to be 5.7 Å, much smaller than our simulation cell.

###  III.3 Submonolayer Dynamics and Spectra

Here, we examine the spectra produced by adsorbate dynamics at submonolayer
coverages, as depicted in Figure 4. The spectra notably contain regions that
scale as $1/\omega$ between 106 and 109 Hz, which overlaps with the range of
typical trap frequencies of 0.1-10 MHz.

The spectra can be divided into two general categories based on temperature.
For moderate to high temperatures ($>$ 50 K), the noise is independent of
frequency (white noise) and its magnitude increases with temperature. For low
temperatures (20-40 K), the high-frequency noise is also flat and greater at
higher temperatures, but the low-frequency noise behaves quite differently. At
these temperatures, the noise levels increase with decreasing frequency with a
mixture of $1/\omega$ and $1/\omega^{2}$ scalings before eventually flattening
to a second white noise region. Moreover, increasing the temperature leads to
a decrease in the low-frequency noise ceiling and an increase in the cutoff
frequency at which the ceiling is reached; this is similar to the behavior of
TLF-derived systems, as discussed in Section II.6.

Figure 4: $S_{\mu}(\omega)$ spectra obtained from dipole moment time series
data of the submonolayer coverage runs. The spectra were smoothed over a
window of 12 timesteps. Figure 5: Selected collective motions observed in
submonolayer trajectories. Arrows point from each particle’s initial position
(represented by the dimmed circles) to its final position. The Au(111) induced
dipole surface (from Figure 1b) is faintly shown beneath the particles. The
Axis Realignment (a), in which the particles reorient around a new axis, and
Cluster Rotation (b) examples are taken from the 0.36 ML trajectory at 20 K,
while the Island Shift (c) example is from the 0.71 ML trajectory at 30 K.
Such motions can take on the order of 10’s or the low 100’s of ps to occur.
Figure 6: Two-dimensional normalized histograms of the relative positions of
CH4 adsorbates on the Au(111) unit cell surface for selected submonolayer
coverages. (a) Surface facets (for reference). (b) Unit cell adsorbate
position histograms. At lower temperatures, a larger proportion of the
particles reside near less favorable atop and bridged sites due to the
strength of the adsorbate clustering. As the temperature increases, the
adsorbate clusters become more fluid, and individual particles are able to
migrate along the surface more often and end up residing in more favorable
sites.

Each of the different aspects of the low temperature spectra can be matched to
different types of movements. The high-frequency white noise stems from
ubiquitous in-place adsorbate rattling of $<$2 Å – both in-plane and out-of-
plane – which causes random and uncorrelated dipole fluctuations. The low-
frequency noise comes from movements – both individual and collective (see
Figure 5) – that are sufficiently rare and substantial to give rise to
different distinct dipole regimes between which the system fluctuates. Figure
S17 illustrates this point by comparing the featureless spectrum of a particle
with effectively no average dipole variation with the
$1/\omega^{2}$-containing spectrum of a particle that transitioned between two
distinct dipole regimes.

At low temperatures, the source – and rarity – of the motions responsible for
low-frequency noise can be traced back to the cohesive strength of the
adsorbate clusters. As Figure 6 shows, these clusters grow more cohesive with
decreasing temperature and increasing $N$, and, at their strongest, are so
cohesive that adsorbates routinely inhabit less-favorable surface facets in
order to minimize their collective CH4-CH4 interactions. In these clusters,
each adsorbate remains in a distinct location on the surface: although the
adsorbates frequently rattle in the $x,y$, and $z$ directions, their average
surface positions do not change. Crucially, this also means that their average
dipole moments do not change (although the continual rattling makes their
dipole values very noisy as a function of time). However, when perturbed
enough, the clusters can undergo reconfigurations. This can be caused by
individual particle motions – such as edge hopping – or, more often,
correlated multi-particle motions – such as collective island shifts, internal
ripples, and cluster rotations (Figure 5). Either way, these reconfigurations
change the global surface positions – and, with them, the average dipole
moments – of one or more adsorbates; dipole changes can result from changes in
the surface sites or $z$ heights the different adsorbates occupy because of
these motions. Such configuration-altering moves are rare at low temperatures
since the CH4-CH4 attraction is high and the system’s kinetic energy is low
(lower than even some of the small surface barriers – see Table SI). Because
they are rare, these reconfigurations effectively divide the trajectories of
each particle into distinct dipole fluctuation regions, each with different
average dipole values.

This collection of dipole fluctuations among the system’s particles aggregate
to produce multiple average dipole states in the summed system dipole time
series data, an example of which can be seen in Figure 7. These multi-level
fluctuators give rise to $1/\omega$ features in the corresponding dipole-
dipole fluctuation spectra in Figure 4.

Figure 7: 10,000-frame moving average of the summed system dipole from a 0.71
ML coverage trajectory run at 25 K. Multiple average dipole states emerge
across the course of the simulation. This set of fluctuations gives rise to
$1/\omega$ features in the corresponding $S_{\mu}(\omega)$ spectra.

At higher temperatures, the CH4-CH4 interactions become smaller relative to
each particle’s kinetic energy. This leads to the clusters losing their
cohesiveness (see Figure 6) as particles more commonly move across the surface
to new global positions. Instead of sampling different dipole regimes in
different adsorbate cluster configurations (as in the low-temperature case),
the particles now freely sample the whole surface with minimal hindrance
(since the adsorbate surface barriers are not significant – see Table SI).
Unlike with the rare collective movements observed at low temperatures, this
situation involves more common but substantially less correlated particle
motions. As stated earlier, trajectories that randomly sample a single
distribution produce white noise; thus, it is not surprising that these fluid
surfaces give rise to featureless spectra. The magnitude of the high-
temperature white noise increases with temperature because more of the higher-
energy $z$ positions (and thus larger dipoles) are accessed by the particles
at higher temperatures.

Temp. | Clusters | Surf. Moves | Corr. | $1/\omega$ | $1/\omega^{2}$
---|---|---|---|---|---
20-40 K | Strong | Rare | High | 106-109 Hz | 106-108 Hz
$>$ 50 K | Weak | Common | Low | Not Seen | Not Seen
Table 3: Summary of submonolayer results.

Table 3 summarizes the trends discussed for the low-temperature and high-
temperature submonolayer spectra.

Identifying correlated motions is not always straightforward, as they come in
many varieties (as can be seen in Figure 5). We used a number of methods to
identify relevant motions. Tracking the number of neighbors each particle can
identify particles migrating around or approaching cluster edges. Examining
how cluster axial angles vary can identify cluster rotations. However, the
most reliable way of identifying meaningful movements was to parse dipole
moment time-series data and look for substantial changes in the average dipole
moment, such as that shown in Figure S17. The averaging step is critical, as
the dipole data is very noisy from timestep to timestep.

It is important to emphasize that the frequency of the movements discussed may
vary from run to run. Because these events are rare, sampling sizes and
stochasticity will invariably play a role in how often they happen in a
(relatively short) simulation, and such effects will be felt more greatly at
lower temperatures. These effects will in turn influence the corresponding
feature frequencies and magnitudes in the dipole-dipole fluctuation spectra.
Nevertheless, these results illustrate how a series of collective motions that
are sufficiently rare (due to temperature, cluster cohesion, and weak-but-
existent surface barriers) can produce $1/\omega$ scaling at low frequencies
on a model trap surface.

###  III.4 Supermonolayer Dynamics and Spectra

We next examine the spectra produced by adsorbate dynamics at monolayer and
supermonolayer coverages. As Figure S14 illustrates, monolayer trajectories
produce featureless spectra. With the surface fully covered, there is little
room for significant surface movements. As with the submonolayer coverages,
some particles reside in higher-energy atop sites due to the competition
between the surface-adsorbate and OPLS interactions, but unlike the
submonolayer adsorbates, the monolayer adsorbates do not move significantly
and thus cannot significantly contribute to the noise or dipole-dipole
fluctuation spectra.

The spectra produced by supermonolayer dynamics, on the other hand, show
several sets of features, as seen in Figure 8. At low temperatures ($<$ 60 K),
spectra display scalings of $\omega^{-1}$ mixed with some $\omega^{-2}$ until
around 35 K, when the spectra become nearly flat (see also Figure S15). At
higher temperatures ($>$ 70 K), however, a new set of larger-magnitude TLF-
like features supersede the low-temperature features. These features increase
in noise magnitude until roughly 90-110 K (depending on the coverage) before
declining with a corresponding increase in cutoff frequency. These high-
temperature supermonolayer features cover a greater range of noise magnitudes
and frequencies than any of the other spectral features we found.

Figure 8: $S_{\mu}(\omega)$ spectra obtained from three different
supermonolayer coverages over a range of temperatures. The spectra were
smoothed over a window of five timesteps.

#### Low-Temperature Supermonolayer Dynamics

Figure 9: Depiction of (a) collective and (b) individual second layer particle
movements and their effect on the average dipoles of surface adsorbates. Both
examples are taken from segments of the 1.16 ML trajectory at 20 K. Arrows
point from each particle’s initial position (represented by the dimmed
circles) to its final position. Here, only first layer particles are colored
with the colormap; second layer particles, which exhibit virtually no dipole
variation (see Figure 11), are colored yellow. As second layer particles move,
they significantly impact nearby first layer particles, especially those that
gain or lose second layer neighbors. In this way, the movement of particles in
the unfilled second layer can cause surface dipole changes and position
shifts.

At low temperatures, adsorbates on filled supermonolayer surfaces are mainly
static (apart from in-place rattling), but they can be affected by second
layer adsorbates in ways that lead to individual and collective first layer
motions (see Figure 9). Even a stationary second layer affects first layer
particles: adsorbates with more neighboring second layer particles (at most
three, as can be seen in Figure 2b) are pushed down to lower surface heights
than those with no such neighbors (see Figure 10). But, second layer particles
can and do migrate around the unfilled second layer to different global
positions. Similar to the submonolayer case, such movements happen more often
with increasing temperature, but second layer particles are even more weakly
held by the surface potential than surface particles despite experiencing OPLS
interactions from first layer particles (see Figure S3) and thus migrate more
frequently at lower temperatures. As these second layer particles migrate – at
rare intervals on the ps time scale but plenty of times over the course of a
$\mathrm{\mu}$s-scale trajectory – they end up influencing different sets of
first layer particles at different points in the trajectory, thus giving rise
to different average surface height (and, with them, average dipole) regimes
for those surface particles at those times (see Figure 9). When enough
particles experience these second layer induced changes, the total system
dipole is moved into a different regime as well, leading to features in the
dipole-dipole fluctuation spectra.

Figure 10: An example of the effect of second layer neighbors on surface
adsorbates, taken from a segment of the 1.16 ML trajectory at 30 K. (a) Number
of second layer particles neighboring Particle #47. (b) 100-frame rolling
average of Particle #47’s distance from the surface. (c) 100-frame rolling
average of Particle #47’s dipole moment. When a second layer particle migrates
to Particle #47’s neighborhood, it pushes the particle to a distinctly
different average dipole regime.

Supermonolayer surfaces can also undergo collective motions as a result of
second layer particle migrations that perturb the surface in the right way, as
can be seen in the examples in Figure 9. However, these motions only seem to
consist of slight cluster shifts and ripples rather than the diversity of
submonolayer collective motions highlighted by Figure 5; this is likely due to
submonolayer surfaces having far more room to maneuver than filled-up
supermonolayer surfaces.

Motions | Coverages | Avg. Dipole Changes
---|---|---
Collective Surface Motions |  $<$1 ML |  $<$5 meÅ
Layer 2 Pushing |  $>$1 ML | 3-10 meÅ
Surface Shifts |  $>$1 ML |  $<$3 meÅ
Layer Exchanges |  $>$1 ML | 38-48 meÅ
Table 4: Summary of motions and the corresponding range of average per-
particle dipole moment magnitude changes that they cause.

Although second layer particles are responsible for first layer dipole regime
changes, their own dipole changes are minimal and do not contribute to overall
noise levels. Figure 11 shows why this is: at typical second layer distances
from the surface, adsorbate dipole moments are nearly zero and have virtually
no variation. Although this makes movement within the second layer noiseless,
it makes movement between layers significant in dipole terms, as the next
section will detail.

Figure 11: Plot of all possible particle $z$ dipole moments (across all $x$
and $y$ surface positions) versus the particle’s distance from the electrode
surface. The dark blue shading indicates how often the (dipole, distance) pair
is generated from the interpolation script at different $(x,y)$ values for
that distance; the bumps around 3-3.75 and 4-5 Å are from dipole sensitivity
to $(x,y)$ changes at those distances. The two layers are marked with mean
surface distances (dotted lines) and 95% confidence intervals (shaded regions)
calculated from supermonolayer particles at 60 K. The first layer experiences
a range of dipole moment changes from small shifts in any direction. Although
the second layer contains virtually no dipole variation, the change in dipole
moment from going between layers is an order of magnitude greater than any
dipole changes within the first layer.

#### High-Temperature Supermonolayer Dynamics

Starting at intermediate temperatures (between 70 and 90 K), the
supermonolayer spectra diverge more decisively from the submonolayer spectra.
Whereas submonolayer spectra do not contain lower frequency noise features at
higher temperatures, supermonolayer spectra demonstrate large low-frequency
shoulders orders of magnitudes larger than any of the submonolayer low-
frequency features. The emergence of such shoulders indicates that there is an
additional, larger noise source present at high temperatures that is absent at
lower coverages.

Figure 12: (a) Surface-adsorbate distances of selected adsorbates taken from
the 1.16 ML trajectory at 70 K. Adsorbate #78 undergoes a layer exchange that
alters its surface-adsorbate distance around 2.75 $\mathrm{\mu}$s, while
adsorbates #76 and #77 remain in their original layers. (b) Corresponding
single-particle dipole-dipole fluctuation spectra of the selected adsorbates.

It turns out that this noise stems from movements of adsorbates between
layers, which we refer to as layer exchanges. Unlike the submonolayer motions
in Figure 5, layer exchanges occur in the span of a few ps, and they have a
larger impact on particle and system dipoles (see Table 4).

Figure 12 demonstrates how layer exchanges translate into the larger spectral
noise observed. Of the three particles depicted, only the one that exchanged
layers produces any low frequency features. Figure 11 shows why layer
exchanges produce more noise than any surface adsorbate motions. When a
particle transitions from the first to the second layer, its dipole moment
drops from roughly 0.044 to nearly zero eÅ, which is a drop several orders of
magnitude larger than any fluctuation in the dipole moment due to first layer
surface motions alone (and compared to the motions shown in Figure 5).

With large dipole changes and long mean state residence times, layer exchanges
can produce significant low frequency TLF noise in the spectra of individual
particles (per Equation (11)). However, the spectra produced by the system’s
total dipole moment will only be affected if layer exchanges change the makeup
of the surface, where the most significant induced dipole moments are possible
(see Figure 11). Therefore, only unidirectional layer exchanges, rather than a
set of reciprocal exchanges, will contribute to the spectra seen in Figure 8.

As noted in Section III.1, our simulation cell can support 68- and
70-adsorbate surfaces in addition to the previously defined monolayer surface
of 69 adsorbates. Second-layer OPLS interactions can make less-favorable
monolayer surface populations more accessible as well. Unidirectional layer
exchanges act to transition the system between these different surface
populations. An example of this is shown in Figure 13. It is this effect of
losing or gaining a particle (and its large induced dipole) from the surface
that leads to the system state fluctuation that produces features in the
spectra in Figure 8.

As with the submonolayer surface movements, the frequency of layer exchanges
may vary from run to run. With mean residence times in the microseconds for
some runs, the spectra can be sensitive to additional exchange events or a
change in the time that they occur. For this reason, these results are not
meant to precisely recreate experimentally-observed spectra. Nevertheless,
these results illustrate how thermally-activated layer exchanges can give rise
to $1/\omega$ scaling regions at low frequencies on a model trap surface.

Figure 13: Layer exchange that results in a surface population change, taken
from a 1.16 ML trajectory segment at 70 K. The yellow circle denotes the
second layer particle before exchanging layers, while the blue circle shows
the same particle after descending to the first layer.

#### Layer Exchanges as Two-State Thermally-Activated Fluctuations

One can consider layer exchanges as a two-state activated process of the kind
discussed in Section II.6. Adsorbates can either inhabit Layer 1 (nonzero
dipole moment) or Layer 2 (near-zero dipole moment) and can stochastically
change layers abruptly. The layer exchanges in our simulations do not
perfectly adhere to this paradigm, as Layer 1 particles exhibit a range of
values with a variance of around $4.5\times 10^{-3}$ eÅ rather than one single
value. However, the Layer 2 particles can be represented by a state with a
dipole value of zero with little loss of precision (see Figure 11). The
typical Layer 1$\rightarrow$2 dipole difference of $4.4\times 10^{-2}$ eÅ is
an order of magnitude larger than the Layer 1 dipole variance, ensuring that
the magnitudes of Layer 1$\leftrightarrow$2 transitions are outside the
distribution of Layer 1 dipole values. Although layer exchanges are individual
motions, this treatment can be extended to the total system dipole moment,
where states can be defined by layer populations each differing by roughly
$4.4\times 10^{-2}$ eÅ.

To see how well this RTS model describes our simulated spectra, we calculated
analytical spectra from Equation (11) based on parameters from the
corresponding simulations. $\Delta I$ was set to $4.4\times 10^{-2}$ eÅ to
represent the magnitude of the typical dipole fluctuation from one adsorbate
vacating or joining the surface. To estimate $\tau_{0}$ and $\tau_{1}$, we
calculated the average residence times of each surface population that the
system experienced. For simplicity, here we present analysis performed on
transitions between two layer populations. In our trajectories, typical
$\tau_{0}$ and $\tau_{1}$ values ranged from the 100’s and 1000’s of ns at
90-110 K to the 10’s of ns at 150 K.

Figure 14 compares TLF analytical spectra calculated from Equation (11) with
spectra obtained from selected 1.42 ML trajectories. In addition, Figure S18
repeats this comparison using per-particle spectra. Although the simulated
spectra are noisier, they follow the same frequency and magnitude trends as
the analytical spectra. In particular, the relative separation between the
different trajectories and the ordering of the low-frequency leveling-off is
captured by the analytical spectra.

Figure 14 and Equation (11) provide a rationale for why the low-frequency
magnitudes of the supermonolayer spectra increase with decreasing temperature.
Increasing system temperature results in more frequent layer exchanges; in RTS
terms, this can be seen as decreasing both $\tau_{0}$ and $\tau_{1}$, since
the more often layer exchanges occur, the less time a particle spends in a
given layer. In the low-frequency limit of Equation (11), the noise magnitude
is inversely proportional to the $((1/\tau_{0})+(1/\tau_{1}))^{2}$ term, which
means that shorter residence times will decrease the magnitude of the noise.
This behavior is similar to that displayed in Figure S16.

Figure 14: Comparison of simulation spectra (lighter shade) and analytical
spectra calculated from the two-level fluctuator model (darker shade) for
selected 1.42 ML trajectories.

The addition of the high-frequency white noise and two-state fluctuator
regions of the spectra gives us the spectral form we see in our supermonolayer
coverages above 70 K.

##  IV Conclusions

In this work, we have analyzed the dynamics of methane adsorbates on gold
electrodes based upon first principles potentials with the aim of identifying
the adsorbate motions that give rise to anomalous heating in ion traps. Based
upon an extensive set of molecular dynamics simulations run at a wide range of
temperatures and surface coverages, we were able to correlate adsorbate
motions such as cluster reconfigurations and layer exchanges with the
different frequency-dependent features of the adsorbate dipole-dipole
fluctuation spectra that can ultimately be related to trap noise. Counter to
prevailing intuition, we find that accounting for adsorbate-adsorbate
correlations is pivotal for observing $1/\omega$ frequency scalings and that
rare, fundamentally collective adsorbate motions are responsible for the
$1/\omega$ features in these spectra at MHz frequencies. In contrast, single
adsorbate well-rattling motions or surface diffusion led to white noise at
frequencies orders of magnitude larger than those realized in ion traps. In
particular, we find that the collective rotations or translations of strongly-
interacting adsorbates within clusters give rise to $1/\omega$ scalings at
submonolayer coverages, while layer exchanges of adsorbates between the first
two adsorbate monolayers are the greatest source of these scalings at
supermonolayer coverages. Key features of our simulated spectra can be
reproduced by two-level fluctuator models parameterized with the intensities
and residence times associated with these motions, strongly corroborating our
findings. This work thus establishes rare, classical adsorbate motions as
atomistic sources of anomalous heating using some of the most realistic
molecular dynamics simulations undertaken to date and dispels the notion that
adsorbate motions are too high in frequency to be key contributors.

That said, our model involved a number of notable simplifications. First, our
methane adsorbates were treated as classical point masses, which enabled us to
treat the methane-methane interactions via an OPLS potential rather than in a
fully ab initio manner. While this is a reasonable simplification for a
spherically-symmetric molecule like methane, it would be a far more severe
approximation for larger hydrocarbons containing many more vibrational and
rotational degrees of freedom. Second, we did not take any nuclear quantum
effects into account. These could become important at the lowest temperatures
explored here and would open up the possibility of the methane adsorbates more
readily traversing the gold surface. While tunneling would not impact the
$1/\omega$ scaling that emerges due to layer exchanges at $>$70 K, it could
affect the forms of the submonolayer spectra, an issue which would have to be
investigated during future studies. Moreover, our single-crystalline gold
electrodes were unrealistically smooth. In practice, most electrodes will be
formed from either polycrystalline materials or possess a wide range of
defects. It is likely that adsorbates traversing such uneven surfaces would
also generate dipole fluctuations and, consequently, trap noise. Lastly,
throughout this work, we have relied upon the patch potential model to link
dipole fluctuations to trap noise even though it is but one approximate model
of many potential models for this noise.Kumph  _et al._ (2016) For these
reasons, our model is not meant to precisely reproduce experimental spectra,
but rather to provide insights regarding the microscopic origins of anomalous
heating that we believe will generalize to even more realistic models and
systems.

Despite these simplifications, our model still reveals a variety of regions
with $1/\omega$ scaling. The fact that it did so without involving multiple
types of adsorbates or accounting for vibrational degrees of freedom suggests
that the mere presence of even simple molecules on electrode surfaces can lead
to anomalous heating. This points to a universal molecular origin for this
phenomenon. Even so, more massive hydrocarbons involved in the trap
fabrication process and that are more likely to be circulating in the ambient
air, such as isopropanol, are much more likely than methane to be adsorbed
onto electrode surfaces. These more massive hydrocarbons would be expected to
interact more strongly with the surface and to exhibit larger, more varied
dipole moments, most likely shifting the $1/\omega$ scaling regimes we
observed to even lower frequencies. More strongly interacting hydrocarbons may
also remedy a deficiency of our results – the fact that we only observe
$1/\omega$ scalings at temperatures well below room temperature – by shifting
the spectra to the higher temperatures employed in most experiments. Molecules
with more vibrational degrees of freedom are additionally more likely to give
rise to motions that can be activated in different frequency regimes, which
would result in the emergence of multiple, overlapping two-level or multi-
level fluctuators that can combine to yield $1/\omega$ scaling over a much
larger frequency regime.Zanolla (2009) We thus look forward to future
simulations of the noise produced by a variety of interacting adsorbates with
many degrees of freedom, which may hold the key to reproducing experimental
spectra with even greater fidelity.

## Data Availability

The data that support this study, including the DFT adsorption energies and
dipole moments used to compute dipole-dipole autocorrelation functions and
fluctuation spectra, are available at https://www.github.com/blfoulon/ion-
trap-anomalous-heating-simulation/.

###### Acknowledgements.

The authors thank Dustin Hite, Jeremy Sage, Jonathon Sedlacek, John
Chiaverini, and Hartmut Häffner for helpful discussions. This work was
performed under the auspices of the U.S. Department of Energy by Lawrence
Livermore National Laboratory under Contract DE-AC52-07NA27344. We also thank
the Brown Presidential Fellowship Program for supporting Yuan Liu. Part of
this research was conducted using computational resources and services at the
Center for Computation and Visualization, Brown University.
B.R. and V.L. oversaw this work equally.

## References

  * Preskill (2018) J. Preskill, Quantum 2, 79 (2018).
  * Wineland  _et al._ (1998) D. Wineland, C. Monroe, W. M. Itano, D. Leibfried, B. King, and D. M. Meekhof, J. Res. Natl. Inst. Stand. Technol. 103, 259 (1998).
  * Turchette  _et al._ (2000) Q. A. Turchette, Kielpinski, B. E. King, D. Leibfried, D. M. Meekhof, C. J. Myatt, M. A. Rowe, C. A. Sackett, C. S. Wood, W. M. Itano, C. Monroe, and D. J. Wineland, Phys. Rev. A 61, 063418 (2000).
  * Brownnutt  _et al._ (2015) M. Brownnutt, M. Kumph, P. Rabl, and R. Blatt, Rev. Mod. Phys. 87, 1419 (2015).
  * Brown  _et al._ (2021) K. R. Brown, J. Chiaverini, J. M. Sage, and H. Haffner, Nature Reviews Materials (2021), 10.1038/s41578-021-00292-1.
  * Wang  _et al._ (2015) C. Wang, C. Axline, Y. Y. Gao, T. Brecht, Y. Chu, L. Frunzio, M. H. Devoret, and R. J. Schoelkopf, Applied Physics Letters 107, 162601 (2015), https://doi.org/10.1063/1.4934486 .
  * Martinis (2009) J. M. Martinis, Quantum Information Processing 8, 81 (2009).
  * Lee  _et al._ (2014) D. Lee, J. L. DuBois, and V. Lordi, Phys. Rev. Lett. 112, 017001 (2014).
  * Saffman (2016) M. Saffman, Journal of Physics B: Atomic, Molecular and Optical Physics 49, 202001 (2016).
  * Kim  _et al._ (2015) M. Kim, H. J. Mamin, M. H. Sherwood, K. Ohno, D. D. Awschalom, and D. Rugar, Phys. Rev. Lett. 115, 087602 (2015).
  * Johnson (1928) J. B. Johnson, Phys. Rev. 32, 97 (1928).
  * Nyquist (1928) H. Nyquist, Phys. Rev. 32, 110 (1928).
  * Hite  _et al._ (2013a) D. Hite, Y. Colombe, A. Wilson, D. Allcock, D. Leibfried, D. Wineland, and D. Pappas, MRS Bulletin 38, 826–833 (2013a).
  * Haffner  _et al._ (2008) H. Haffner, C. F. Roos, and R. Blatt, Physics Reports 469, 155 (2008).
  * Benhelm  _et al._ (2008) J. Benhelm, G. Kirchmair, C. F. Roos, and R. Blatt, Nature Physics 4, 463 (2008).
  * Bruzewicz  _et al._ (2019) C. D. Bruzewicz, J. Chiaverini, R. McConnell, and J. M. Sage, Applied Physics Reviews 6, 21314 (2019).
  * McCormick  _et al._ (2019) K. C. McCormick, J. Keller, S. C. Burd, D. J. Wineland, A. C. Wilson, and D. Leibfried, Nature 572, 86 (2019).
  * Degen  _et al._ (2017) C. L. Degen, F. Reinhard, and P. Cappellaro, Rev. Mod. Phys. 89, 035002 (2017).
  * Baumgart  _et al._ (2016) I. Baumgart, J.-M. Cai, A. Retzker, M. B. Plenio, and C. Wunderlich, Phys. Rev. Lett. 116, 240801 (2016).
  * Leibfried  _et al._ (2004) D. Leibfried, M. D. Barrett, T. Schaetz, J. Britton, J. Chiaverini, W. M. Itano, J. D. Jost, C. Langer, and D. J. Wineland, Science 304, 1476 (2004).
  * Kozlov  _et al._ (2018) M. G. Kozlov, M. S. Safronova, J. R. Crespo López-Urrutia, and P. O. Schmidt, Rev. Mod. Phys. 90, 045005 (2018).
  * Chou  _et al._ (2010) C. W. Chou, D. B. Hume, J. C. J. Koelemeij, D. J. Wineland, and T. Rosenband, Phys. Rev. Lett. 104, 070802 (2010).
  * Huntemann  _et al._ (2016) N. Huntemann, C. Sanner, B. Lipphardt, C. Tamm, and E. Peik, Phys. Rev. Lett. 116, 063001 (2016).
  * An  _et al._ (2019) D. An, C. Matthiesen, E. Urban, and H. Haffner, Phys. Rev. A 100, 063405 (2019).
  * Daniilidis  _et al._ (2011) N. Daniilidis, S. Narayanan, S. A. Moller, R. Clark, T. E. Lee, P. J. Leek, A. Wallraff, S. Schulz, F. Schmidt-Kaler, and H. Haffner, New Journal of Physics 13, 013032 (2011).
  * Sedlacek  _et al._ (2018a) J. A. Sedlacek, J. Stuart, D. H. Slichter, C. D. Bruzewicz, R. McConnell, J. M. Sage, and J. Chiaverini, Phys. Rev. A 98, 063430 (2018a).
  * Sedlacek  _et al._ (2018b) J. A. Sedlacek, A. Greene, J. Stuart, R. McConnell, C. D. Bruzewicz, J. M. Sage, and J. Chiaverini, Phys. Rev. A 97, 020302 (2018b).
  * Noel  _et al._ (2019) C. Noel, M. Berlin-Udi, C. Matthiesen, J. Yu, Y. Zhou, V. Lordi, and H. Häffner, Phys. Rev. A 99, 063427 (2019).
  * Labaziewicz  _et al._ (2008) J. Labaziewicz, Y. Ge, P. Antohi, D. Leibrandt, K. R. Brown, and I. L. Chuang, Phys. Rev. Lett. 100, 013001 (2008).
  * Deslauriers  _et al._ (2006) L. Deslauriers, S. Olmschenk, D. Stick, W. K. Hensinger, J. Sterk, and C. Monroe, Phys. Rev. Lett. 97, 103007 (2006).
  * Hite  _et al._ (2017) D. A. Hite, K. S. McKay, S. Kotler, D. Leibfried, D. J. Wineland, and D. P. Pappas, MRS Advances 2, 2189–2197 (2017).
  * Kumph  _et al._ (2016) M. Kumph, C. Henkel, P. Rabl, M. Brownnutt, and R. Blatt, New Journal of Physics 18, 023020 (2016).
  * Teller  _et al._ (2021) M. Teller, D. A. Fioretto, P. C. Holz, P. Schindler, V. Messerer, K. Schüppert, Y. Zou, R. Blatt, J. Chiaverini, J. Sage, and T. E. Northup, arXiv e-prints , arXiv:2103.13846 (2021), arXiv:2103.13846 [quant-ph] .
  * Hite  _et al._ (2012) D. A. Hite, Y. Colombe, A. C. Wilson, K. R. Brown, U. Warring, R. Jördens, J. D. Jost, K. S. McKay, D. P. Pappas, D. Leibfried, and D. J. Wineland, Phys. Rev. Lett. 109, 103001 (2012).
  * Daniilidis  _et al._ (2014) N. Daniilidis, S. Gerber, G. Bolloten, M. Ramm, A. Ransford, E. Ulin-Avila, I. Talukdar, and H. Haffner, Phys. Rev. B 89, 245435 (2014).
  * McKay  _et al._ (2014) K. S. McKay, D. A. Hite, Y. Colombe, R. Jördens, A. C. Wilson, D. H. Slichter, D. T. C. Allcock, D. Leibfried, D. J. Wineland, and D. P. Pappas, arXiv e-prints , arXiv:1406.1778 (2014), arXiv:1406.1778 [physics.atom-ph] .
  * McConnell  _et al._ (2015) R. McConnell, C. Bruzewicz, J. Chiaverini, and J. Sage, Phys. Rev. A 92, 020302 (2015).
  * Allcock  _et al._ (2011) D. T. C. Allcock, L. Guidoni, T. P. Harty, C. J. Ballance, M. G. Blain, A. M. Steane, and D. M. Lucas, New Journal of Physics 13, 123023 (2011).
  * Safavi-Naini  _et al._ (2011) A. Safavi-Naini, P. Rabl, P. F. Weck, and H. R. Sadeghpour, Phys. Rev. A 84, 023412 (2011).
  * Safavi-Naini  _et al._ (2013) A. Safavi-Naini, E. Kim, P. F. Weck, P. Rabl, and H. R. Sadeghpour, Phys. Rev. A 87, 023421 (2013).
  * Dutta and Horn (1981) P. Dutta and P. M. Horn, Rev. Mod. Phys. 53, 497 (1981).
  * Zanolla (2009) N. Zanolla,  (2009), 10.6092/unibo/amsdottorato/1333.
  * Weissman (1988) M. B. Weissman, Rev. Mod. Phys. 60, 537 (1988).
  * Ray  _et al._ (2019) K. G. Ray, B. M. Rubenstein, W. Gu, and V. Lordi, New Journal of Physics 21, 053043 (2019).
  * Lloyd  _et al._ (2021) P. N. T. Lloyd, V. Walther, and H. Sadeghpour, arXiv e-prints , arXiv:2105.00559 (2021), arXiv:2105.00559 [quant-ph] .
  * Kim  _et al._ (2017) E. Kim, A. Safavi-Naini, D. A. Hite, K. S. McKay, D. P. Pappas, P. F. Weck, and H. R. Sadeghpour, Phys. Rev. A 95, 033407 (2017).
  * Hite  _et al._ (2013b) D. Hite, Y. Colombe, A. Wilson, D. Allcock, D. Leibfried, D. Wineland, and D. Pappas, MRS Bulletin 38, 826–833 (2013b).
  * Gesley and Swanson (1985) M. A. Gesley and L. W. Swanson, Phys. Rev. B 32, 7703 (1985).
  * Dubessy  _et al._ (2009) R. Dubessy, T. Coudreau, and L. Guidoni, Phys. Rev. A 80, 031402 (2009).
  * Rossi and Opat (1992) F. Rossi and G. I. Opat, Journal of Physics D: Applied Physics 25, 1349 (1992).
  * Kresse and Furthmüller (1996) G. Kresse and J. Furthmüller, Phys. Rev. B 54, 11169 (1996).
  * Berland and Hyldgaard (2014) K. Berland and P. Hyldgaard, Phys. Rev. B 89, 035412 (2014).
  * Jorgensen  _et al._ (1984) W. L. Jorgensen, J. D. Madura, and C. J. Swenson, Journal of the American Chemical Society 106, 6638 (1984).
  * ase (2017) Journal of Physics: Condensed Matter 29, 273002 (2017).
  * (55) https://www.github.com/blfoulon/ion-trap-anomalous-heating-simulation/.
  * Press and Rybicki (1989) W. H. Press and G. B. Rybicki, Astrophys. J. 338, 277 (1989).
  * Blackman and Tukey (1958) R. B. Blackman and J. W. Tukey, Bell System Technical Journal 37, 185 (1958), https://onlinelibrary.wiley.com/doi/pdf/10.1002/j.1538-7305.1958.tb03874.x .
  * Shynk (2012) J. J. Shynk,  _Probability, Random Variables, and Random Processes: Theory and Signal Processing Applications_ (Wiley, 2012).
  * Machlup (1954) S. Machlup, Journal of Applied Physics 25, 341 (1954).