Surface-induced decoherence and heating of charged particles

Levitating charged particles in ultra-high vacuum provides a preeminent platform for quantum information processing, for quantum-enhanced force and torque sensing, for probing physics beyond the standard model, and for high-mass tests of the quantum superposition principle. Existing setups range from single atomic ions, to ion chains and crystals, to charged molecules and nanoparticles. Future technological applications of such quantum systems will be crucially affected by fluctuating electric fields emanating from nearby electrodes, which interact with the levitated particles' monopole and higher charge moments. In this article, we provide a theoretical toolbox for describing how the rotational and translational quantum dynamics of charged nano- to microscale objects is affected by near metallic and dielectric surfaces, as characterized by their macroscopic dielectric response. The resulting quantum master equations describe the coherent surface-particle interaction, due to image charges and Casimir-Polder potentials, as well as surface-induced decoherence and heating, with the experimentally observed frequency and distance scaling. We explicitly evaluate the master equations for typical charge distributions and types of motion, thereby providing the tools required for describing and mitigating surface-induced decoherence in a variety of experiments with charged objects.

Controlling atomic, molecular, and nanoscale ions with electrodes and superconducting circuits [59][60][61][62] brings them close to metallic and dielectric surfaces. Thermally fluctuating charges and currents in the surface then unavoidably induce decoherence and motional heating of the trapped particle [10,[63][64][65][66]. Given that building a quantum sensing device or a quantum computer will likely require chip-based miniaturization [65], the role of surface induced heating and decoherence can be expected to become worse in future devices, presenting a potential roadblock. This calls for a universal theoretical framework for describing charge-induced decoherence and heating close to metallic and dielectric surfaces, which allows incorporating microscopic and heuristic models, and which points the way towards noise mitigation techniques.

B. Surface-induced heating and decoherence
The experimental characterization of surface noise and decoherence relies on measuring heating rates and coherence times of trapped ions near the surface [2,[65][66][67][68][69][70][71] and of interfering free electron beams travelling above a plane [72][73][74]. While different setups can exhibit distinct noise spectral densities, most experiments seem to agree on a characteristic inverse frequency scaling and on a proportionality to the fourth power of inverse distance [66]. That said, the microscopic origin of surface noise is still debated [2,65,66,71]. Theoretical models range from the damping of conduction electrons [75], the dissipative formation of image charges [76], the excitation of surface plasmons [77], to the fluctuations of surface-bound [78,79] or diffusing [80] adatoms, to patch potentials [81][82][83][84], and to thermally activated two-level adsorbates [85]. In an approach that does not require specifying a microscopic mechanism for surface-induced heating and decoherence, one may characterize the surface material solely by means of its dielectric response function [86][87][88][89][90].
This mechanism-independent approach can be formalized by using the framework of macroscopic quantum electrodynamics [91,92], which describes the quantum and thermal noise of the quantized electromagnetic field in terms of fluctuating bound charges and currents. This framework expresses the observable forces on the particle through empirically accessible response functions, while providing a versatile, i.e. non-material specific, description. Macroscopic quantum electrodynamics can reproduce the classical field fluctuations found above lossy dielectric layers [89,90,93], including the experimentally observed magnitude, distance, and frequency scaling [66]. This supports that the noise originates from a thin surface layer. In addition, macroscopic quantum electrodynamics has been employed to quantitatively describe the decoherence of fast electrons [74,88], of spins [94], and of the center-of-mass motion of magnetic [95] and polarizable [96] spheres above surfaces.
Future quantum technologies, which utilize molecules and nanoparticles levitating and rotating in micro-traps near electrodes and micro-structured surfaces, will require a comprehensive toolbox to predict surface-induced heating and decoherence for arbitrary charge distributions and surface configurations. This article provides this toolbox in the form of Lindblad master equations for the ro-translational quantum dynamics of arbitrarily charged particles in front of general dielectrics.

C. This work
This article derives quantum master equations for charged rigid bodies above general metallic and dielectric surfaces. These master equations contain a coherent particle-surface potential, accounting for the formation of image charges and the impact of Casimir-Polder interactions, and a dissipator, describing surface-induced decoherence and heating. We express the resulting translational and rotational heating and decoherence rates in terms of the particle charge distribution, the surface geometry and layer structure, and the surface dielectric response and surface temperature. Based on the general theory, we provide specific formulas for typical charge distributions, surface dielectric functions, and types of motion. This includes the oscillation of monopoles, the libration of dipoles, and the rotation of quadrupoles above metallic and superconducting surface covered by thin dielectric layers, whose noise spectrum can decohere the particle motion either resonantly or off-resonantly.
The chosen framework facilitates a pragmatic description of surface-induced decoherence and heating in terms of macroscopic dielectric response functions. The latter are accessible both by in situ and ex situ measurements as well as by microscopic ab-initio calculations [97,98] and molecular-dynamics simulations [99]. By disentangling the effects of the surface geometry from the influence of the surface material, the derived mas-Figure 1. A charged particle (blue) with charge distribution 0 evolving in a spatio-orientational superposition state |Ψ(R, Ω) will suffer from motional decoherence due to electric fields (red lines) originating from a fluctuating polarization (black arrows) in nearby metallic or dielectric surfaces. The fluctuation strength of the polarization at a position r depends on the dielectric material response ε0εr(r, ω) and the temperature T . For small distances, the electric field propagates quasistatically from each fluctuation to the particle, as can be described by a scalar potential φ(r, t). A layered dielectric consisting of a surface layer with ε0εs(ω) and a bulk dielectric with ε0ε b (ω) can explain the noise levels and scalings measured in ion traps [89].
ter equations can be used for the design of surfaces with favourable noise properties given a specific charge distribution and ro-translational motion, providing a bridge to device engineering and materials design [66]. We expect the here presented toolbox to be relevant for a variety of state-of-the art setups ranging from trapped ion quantum computers [66], to hybrid quantum devices with charged atoms or molecules [10], to nanoparticles in electric and optical traps for fundamental tests of physics [38].
The article is structured as follows: Section II summarizes the theoretical framework of surface-induced decoherence and presents the resulting master equations in their most generic form. The master equations are then specified in Sects. III and IV for typical charge distributions, for metallic, dielectric, and superconducting surfaces, and for planar and layered surface geometries. In Sect. V we discuss criteria how the presented toolbox can be utilized, and we illustrate this with two examples. Finally, Sect. VI discusses possible future research directions and provides our conclusions. The technical details of the derivation are presented in Apps. A to C. Appendix D discusses the ro-translational decoherence due to elastic photon scattering (Thomson scattering).

A. Particle-surface coupling
We consider the motion of a particle of mass m and with moments of inertia I 1 , I 2 , I 3 moving in vacuum at center-of-mass position R. The particle orientation Ω is specified by the rotation tensor R Ω , relating the spacefixed axes e i to the body-fixed principal axes n i = R Ω e i , so that the inertia tensor is I Ω = i I i n i ⊗ n i . The body carries a rigid charge distribution (r) = 0 [R T Ω (r − R)], where 0 is the charge distribution at the reference position R = 0 and orientation R Ω = 1.
The particle interacts with the fluctuating electromagnetic field emanating from a nearby surface, see Fig. 1. This field is due to thermal polarization currents in the surface material and can be expressed via the scalar and vector gauge potential φ(r, t) and A(r, t). The resulting particle dynamics are described by the minimal coupling Hamiltonian Here, P and J are the canonical linear and angular momentum vectors. The canonical angular momentum vector is understood in terms of the canonical momentum coordinates p Ω , conjugate to the orientation coordinates Ω, see App. A 1. (It is not identical to the kinetic angular momentum due to gauge coupling to the transverse field.) Promoting the gauge fields to Heisenberg-picture operators, their dynamics can be described within the framework of macroscopic quantum electrodynamics [91,100]. This theory relates the gauge fields at position r to the quantized polarization currents −iωP N (r , ω) at position r and frequency ω, which fluctuate in the dielectric with complex-valued permittivity ε r (r , ω) (App. A 2). The polarization noise field is modelled by with bosonic annihilation operators f e (r, ω), from now on dubbed as polarization oscillators. They satisfy the standard (dyadic) canonical commutation relations [91], e.g., In the following, we assume that the oscillators are always in a thermal state characterized by the surface temperature T , so that Ω Ω (R , R ) − 2h Im g(r, r , ωrot) r=r =Rcm with the Bose-Einstein occupation factor n(ω) The polarization noise field (2) induces an electromagnetic field at the particle position, as given by Maxwell's equations. In most experimental setups, this field is well approximated by its longitudinal part, and therefore described by a quasistatic scalar Green function (App. C) The Green function g(r, r , ω) determines the electrostatic potential at r of a charge oscillating in magnitude with frequency ω at position r in presence of a dielectric. This Green function fulfills the relation which will be useful in the following. App. C derives this relation together with further properties of g(r, r , ω).
In the quasistatic approximation (C5) the vector potential vanishes in Coulomb gauge so that the particlesurface interaction is determined by the operator-valued scalar potential where Im[ε r (r, ω)] ≥ 0 for positive frequencies [91]. The hermitian Heisenberg-picture operator (7) acts on the Hilbert space of the polarization oscillators. Having specified how the particle-field coupling Hamiltonian relates to the surface noise, one can obtain the reduced quantum dynamics of the particle by tracing out the field dynamics [101]. For realistic experimental situations, the weak-coupling Born-Markov approximation is justified (Sec. V B). The resulting quantum master equation depends on the relation of the involved timescales. We will start with the slow-particle limit before discussing the case of resonant motion.

B. Slow particle motion
In the slow-particle limit, the surface-free particle dynamics due to the Hamiltonian H 0 and the surfaceinduced particle dynamics are much slower than the surface-fluctuation correlation time.
As shown in App. B 2, the resulting master equation takes Lindblad form with Lindblad operators The limit is taken from above because Im[ε r (r, ω)] can exhibit a discontinuity at ω = 0 due to the Schwarz reflection principle. The surface-interaction Hamiltonian H si describes the coherent interaction with image charges in the dielectric, Equations (8)-(10) describe ro-translational decoherence and heating in terms of (i) the dielectric properties of the surface material, (ii) the surface geometry, and (iii) the particle charge distribution. As one expects, only the dielectric regions of space, where Im[ε r (r, ω)] > 0, contribute to the integral in (9). Note that the form of Eq. (10) implies that the Lindblad operators correspond to the (negative) displacement field induced by the particle charge distribution located at position R and orientation Ω.
Using relation (6), one can disentangle the contributions of charge distribution and surface properties by showing that (App. B 2) Here, the surface geometry and response function enter through the surface-fluctuation kernel The particle charge distribution takes the role of the Lindblad operators, leading to a localization in position and orientation. The decoherence rate in this basis, defined via Table II. The resonant surface-fluctuation kernel for particles in front of a homogeneous planar surface (top) and a planar surface covered by a dielectric layer (bottom). For the latter case we assume the particle to be at a large distance from the surface, r · e3, r · e3 ds, as described in Sec. IV. In both expressions r and r are evaluated in vacuum, and the mirror tensor M = 1 − 2e3 ⊗ e3 serves to reflect a vector at the surface plane.
Illustration Surface-fluctuation kernel Im g(r, r , ω) = − 2 4πε0 where we used g(r, r , ω) = g(r , r, ω) (see App. C 2). The decoherence rate is non-negative, as follows from (6), it vanishes for R = R and Ω = Ω , and remains bounded for large lateral superposition sizes. It thus describes how superposition states of the charged particle gradually lose coherence due to their interactions with surface charge fluctuations.
The Lindblad operators are specified in Sec. III for common charge distributions, while the form of the surface-fluctuation kernel (13) is discussed in Sec. IV for typical surface geometries and generic material responses. Table I summarizes the master equations for particles characterized by their first multipole moments and for typical types of particle motion. Table II gives the surface-fluctuation kernel for a half-space and a layered surface geometry.

C. Resonant motion
The resonant limit applies if both the free particle motion and the decay of surface-fluctuation correlations are much faster than the surface-induced particle dynamics. The resulting master equation can be stated by diagonal-izing the surface-free particle dynamics, including a possible external trapping potential V ext . This allows expanding the particle charge density by means of the operators with coefficients where ω labels the level spacings of H 0 . A rotating wave approximation of the interaction Hamiltonian with respect to the particle dynamics yields a dissipator of Lindblad form (App. B 3), with L e (r) = √ 2 d 3 s c (s)ε r (r, ω ) ∂ ∂r g(r, s, ω ).
The Lindblad operators thus describe the displacement field smeared out over the field-free motion of the charged particle.
Only the dielectric regions contribute to the integral in the master equation (19). The latter describes resonant transitions between the particle motion (16) and the dielectric, as weighted by the energy loss function Im[ε(r, ω )]/|ε(r, ω )| 2 . The first term in (19) accounts for spontaneous and stimulated de-excitation of the particle motion; the second term accounts for surface-induced excitations of the particle dynamics. This leads to heating, damping, and thermalization of the particle motion with the surface temperature.
The associated surface-interaction Hamiltonian (Lamb shift) reads in the zero-temperature limit It contains two contributions: The first accounts for the interaction of the particle charge distribution with the particle-induced surface polarization -i.e. the imagecharge potential (see App. B 3). The second contribution describes modifications of the particle dynamics induced by zero-point fluctuations in the surface, similar to the Casimir-Polder potential. One can again disentangle the contribution of particle and surface properties by using Eq. (6). This yields where the surface properties now enter through the imaginary part of the Green function in (5), evaluated on resonance with the particle motion. The Lindblad operators c (r) are fully specified by the particle charge distribution 0 and the field-free Hamiltonian H 0 , implying that the integrands in (22) are localized within the smeared-out charge distribution.
The form of these Lindblad operators is discussed in Sec. III for typical charge distributions. The resonant surface-fluctuation kernel Im[g(r, r , ω )] of typical surface geometries and generic material responses is then specified in Sec. IV, and Sec. V discusses the resulting heating and decoherence rates for state-of-the-art setups. Master equations for typical charge distributions and types of motion are illustrated in Tab. I. The reso-nant surface-fluctuation kernel is given in Tab. II for a half-space and a layered surface geometry.

III. LINDBLAD OPERATORS FOR CHARGE MULTIPOLES
In the master equations (12) and (22) the particle properties, entering through the Lindblad opertors, are separated from the dielectric surface, entering through the surface-fluctuation kernel. In the slow particle limit, only the instantaneous charge distribution appears, while the resonant limit blurs the charge over the fast particle motion, as described by the operators c (s). This section gives the Lindblad operators for typical experimental setups. The resulting master equations are provided in Tab. I.

A. Slow particle motion
a. Slow monopole-The Lindblad operator of a point monopole with total charge q at position r is independent of the particle orientation, b. Slow dipole-A point particle with a body-fixed dipole moment p 0 , but no other multipoles, is characterized by 0 (r) = −p 0 · ∂ r δ(r), yielding the Lindblad operator which involves the orientation state-dependent dipole vector p Ω = R Ω p 0 . After inserting (24) into the master equation (12), the spatial derivatives of the δdistributions can be shifted to h(r, r ) by partial integration. c. Slow quadrupole-A point particle with a body fixed quadrupole tensor Q 0 , but no other multipoles, is characterized by 0 (r) = ∂ r ·Q 0 ∂ r δ(r)/6, leading to the Lindblad operators Similar to the case of a slow dipole, the Liouvillan (12) depends on the orientation state-dependent quadrupole tensor Q Ω = R Ω Q 0 R T Ω and the Hessian of h(r, r ) at the particle position.
General charge distribution-The dissipator for a particle with a charge distribution characterized by several multipole moments is not simply obtained by adding the multipole dissipators. The Lindblad operator in (12) is given by the operator-valued charge distribution, implying that the master equation contains cross-terms involving different multipole moments.

B. Resonant motion
d. Oscillating monopole-The quantum dynamics of an oscillating monopole at position R = R eq can be decomposed into three bosonic modes with operators a , frequencies ω , and oscillation directions , 2mω a e −iω t + h.c. .
(26) If the frequencies are sufficiently well separated, so that the rotating wave approximation applies, and for amplitudes sufficiently small, so that field can be approximated locally, the Lindblad operators take on the form (App. B 3) e. Oscillating dipole-The corresponding calculation for an oriented dipole with time-independent dipole vector p whose center of mass is oscillating, yields the Lindblad operators f. Oscillating quadrupole-Likewise, for an oscillating oriented quadrupole with time-independent quadrupole tensor Q, one obtains g. Librating dipole-A linear rotor with moment of inertia I and dipole operator p Ω = pn 3 parallel to the symmetry axis n 3 is assumed to librate at R cm in the two directions α , β perpendicular to its equilibrium orientation eq . Denoting the frequencies of the libration modes by ω α , ω β and the bosonic operators by a α , a β , the rotor symmetry axis evolves as For sufficiently well-separated libration frequencies, the Lindblad operators take the form h. Librating quadrupole-Likewise, if the linear rotor is a quadrupole Q Ω = Q 33 [3n 3 ⊗ n 3 − 1] /2 with symmetry axis n 3 , librating as in (30), the Lindblad operators are given by i. Rotating planar dipole-A linear rotor located at R cm is assumed to rotate in a plane spanned by 1 and 2 . For sufficiently sharp angular momentum p α the symmetry axis n 3 rotates approximately at a frequency ω rot = p α /I determined by the moment of inertia I, (33) with Euler-angle operator α. For a dipolar charge distribution p Ω = pn 3 with symmetry axis n 3 the Lindblad operators read as j. Rotating planar quadrupole-Likewise, for a quadrupolar charge distribution with quadrupole tensor k. Freely rotating dipole-Finally, we derive the Lindblad operators of a freely rotating linear rotor at a position R cm with the dipole operator specified by the Euler angle operators α and β, The energy eigenstates entering (17) and (18) are given by the spherical harmonics α, β| , m = Y m (β, α), with associated energies E = 2 ( + 1)/2I. Expanding the dipole operator in the angular momentum basis m m =(−1) m (2 + 1)(2 + 1) The expansion coefficients (18) thus depend on four quantum numbers, which are nonzero only for | − | = 1 and m−m = ±1, 0 due to the Wigner 3-j symbols in (38). Since (17) is restricted to positive energy spacings, only the Lindblad operators with frequencies ω := (E − E )/ = ( + 1)/I contribute. We will use the operators (40) in Sec. V D to simulate the rotational decoherence of a dipole above a superconducting surface covered with a dielectric layer.
General charge distribution-The dissipator for a charge distribution characterized by several multipole moments does in general not decompose into a sum of pure multipole dissipators. However, this sum may be a good approximation if the characteristic frequencies associated with the motion of the distinct multipole components are sufficiently well separated, so that all crossterms effectively average to zero.

IV. SURFACE-FLUCTUATION KERNEL
Having specified the Lindblad operators, we now discuss the surface-fluctuation kernel, which depends on the geometry of the surface and its dielectric response, as characterized by the relative permittivity ε r (r, ω). The latter is subject to the Kramers-Kronig relations, implying that high frequency resonances described by peaks in the imaginary part of ε r (r, ω) contribute to the low frequency behaviour of the real part ( Fig. 2). As we will see below, even weak resonances due to surface-adsorbed molecules and imperfections can strongly affect the energy loss function at the typical trapping frequencies of ions and nanoparticles.
A simple model for the dielectric response is given by a sum of Drude-Lorentz oscillators [102] The imaginary part exhibits peaks at the resonance frequencies ω n , whose relative weights are determined by the (real-valued) coefficients f n , and whose widths are determined by the damping rates γ n . The special case of a single resonance at ω n = ω r with f n = ω 2 p /ω 2 r yields the Drude-Lorentz model for a single oscillator In the limit that the resonance frequency vanishes, this is referred to as the Drude metal, We will also discuss decoherence close to a superconducting surface, with critical temperature T c , at temperature T < T c . Its permittivity can be written in the two-fluid model as [103,104] where Λ 0 is a material-specific constant. For temperatures T approaching T c from below, this superconductor turns into a normal metal. At much smaller temperatures it behaves like an ideal electron gas without damping.

A. Homogeneous half-space
We consider a particle levitated in vacuum at r · e 3 > 0 above a semi-infinite dielectric half-space at r·e 3 < 0 with homogeneous permittivity ε b (ω). The Green function (5) for this configuration can be found in App. C 4. For the special case that r and r are both in vacuum, as appearing in the master equations (12) and (22), it is given in Tab. II. The imaginary part of the Green function decays with the inverse distance between the location r and the mirror image of r at the location Mr , with the mirror tensor M = 1 − 2e 3 ⊗ e 3 . In the special case of the Drude-Lorentz model (41) one obtains and for the Drude-metal half-space We will use this surface-fluctuation kernel to calculate below the decoherence rate of electron beams near a metallic surface. For a superconductor, the surface-fluctuation kernel follows from (44). In the static limit one obtains h(r, r ) = 0 since n(ω)Im [2/ [ε sc (ω) + 1]] vanishes. This means that there is no charge-induced decoherence of a slowly moving charge distribution in front of a pristine superconductor.

B. Half-space covered by thin surface layer
In typical experimental situations, the surface may be covered by a thin dielectric surface layer. The presence of this layer can strongly affect the charge-induced decoherence and heating, even though its impact on the overall dielectric response of the surface is rather weak. Denoting the layer permittivity by ε s (ω) and its thickness by d s , one finds the imaginary part of the Green function for r and r in vacuum where Here, J 0 (·) denotes the zeroth order Bessel function. The full Green function is calculated in App. C 5.
For large distances, r·e 3 d s and r ·e 3 d s , Eq. (47) can be expanded to first order in kd s 1, yielding the expression given in Tab. II. In case of a low-permittivity dielectric layer on a metallic (or superconducting) bulk, |ε b (ω)| |ε s (ω)| 1, the metal (or superconductor) effectively acts as a mirror. One then obtains For a dielectric layer characterized by the surface permittivity ε s (ω) = ε DL (ω) on top of a Drude metal, Here ε DL (0) = 1 + n≥0 f n is the static value of the surface response.
Replacing the metal by a superconductor, ε b (ω) = ε sc (ω), one obtains which is equivalent to (51) for negligible damping γ m = 0 in the bulk. Like in (50), the superconductor acts like a mirror for the fluctuations generated in the dielectric layer. This approximation is justified if the surface layer is thinner than the London penetration depth, which is the case for typical superconducting materials.
C. Electric-field power spectral density The surface noise kernel also characterizes the electric field fluctuations observed above thermal surfaces. As we shall see below, the resulting electric field power spectral density (PSD) is closely related to the decoherence and heating rates of electric monopoles and dipoles.
The PSD at position R 0 is defined as the Fourier transform of the autocorrelation function, (53) where E(r, t) denotes the electric field operator, described by the scalar potential (7), and {·, ·} is the dyadic anti-commutator. For a thermal distribution one finds in the quasistatic limit which reduces for ω = 0 to For instance, a thin dielectric layer, as described by Eq. (50), gives rise to In the large temperature limit k B T / ω 1, this is equivalent to the expression given in [89]. On the other hand, the Drude-metal half-space with (43) yields in the limit of small frequencies which is consistent with the calculations in [89,105][106]. As a final example, a Drude-Lorentz layer on a superconductor with (52) yields These examples illustrate that the presence of a thin dielectric layer induces the experimentally observed scaling with d = R 0 · e 3 the distance to the surface. Importantly, this contribution from the layer dominates for typical electrode and surface layer materials [89]. We will see in the following section that this distance behavior translates to nanoparticle decoherence and heating rates.

D. Relevance of low-frequency surface excitations
Inserting the surface-fluctuation kernel into the master equation (22), for instance using Eq. (50), demonstrates that the magnitude of surface-induced decoherence is determined by the energy loss function Im[ε r (ω)]/|ε r (ω)]| 2 and by its thermal version The latter dominates the former for small frequencies (see Fig. 2b), and thus determines the high-temperature frequency behaviour of the electric field PSD (56). For instance, for a Drude-Lorentz dielectric (41) (Fig. 2a), the absolute value in the denominator of (60) can be approximated by its static value ε DL (0), so that the function plateaus at for frequencies much less than ω n , ω 2 n /γ n . Typically, high-frequency electronic transitions exhibit significantly larger weights f n and thus dominate the static value of the dielectric permittivity. However, Eq. (60) can be enhanced by broad low-frequency contributions to the spectrum, each contributing with relative strength f n γ n /ω 2 n . Such weak and broad absorption bands could for instance be due to weakly bound surface adsorbates [66,99].
Moreover, let us assume that these low lying transitions give rise to a broad plateau in the absorption spectrum with mean imaginary dielectric response ε im (see Fig. 2a). It then follows that the function l(ω) exhibits the experimentally observed 1/ω frequency scaling (Fig. 2b). Specifically, one obtains This illustrates the decisive role played by the lowfrequency dielectric response of the electrode surface in ion traps. As we shall see below, this frequency scaling carries over to the surface-induced decoherence and heating rates. We note that non-flat combinations of lowlying resonances can give rise to other power-law scalings [107].

V. HOW TO USE THE TOOLBOX
Having discussed how the particle charge distribution and the surface properties enter the master equations (12) and (22), we now explain how one can identify the suitable master equation for a given experimental setting. This will be followed by two examples and a short discussion of relevant limiting cases, in which the master equations simplify to known results.

A. Finding the right master equation
To identify the appropriate master equation for a given physical situation requires comparing the involved timescales. These are determined by the setup, as described by the form and location of the particle chargedistribution, its free dynamics, as well as the surface geometry, dielectric response, and temperature. Given this information, one proceeds as follows: 1. Specify the body-fixed charge distribution of the particle, e.g. in terms of its multipole moments, and estimate the characteristic timescales of the particle dynamics in absence of the surface.
2. Specify the surface response function, e.g. by experimental measurements or by a suitable model such as Eq. (41).
3. Check that electric-field retardation is negligible and calculate the quasistatic Green function (5).
If the surface is clean/covered by a thin dielectric layer, one can use the first/second row in Tab. II.
5.1 If the particle dynamics are much slower than the surface-fluctuation correlation time, use the slowparticle master equation (12).
If the particle charge distribution is characterized exclusively by its monopole, dipole, or quadrupole moment, one can use the master equations a., b., or c. in Tab. I. Otherwise, insert the charge distribution into (12).
5.2 If both the free particle dynamics and the decay of surface-fluctuation correlations are much faster than the surface-induced system dynamics, use the resonant master equation (22).
If the particle centre-of-mass is oscillating and its charge distribution is characterized solely by its monopole, dipole, or quadrupole moment, one can use the master equations d., e., or f. in Tab. I. If the particle orientation is librating and its charge distribution is characterized solely by its dipole or quadrupole moment, one can use the master equations g. or h. in Tab. I. If the particle is rotating in a plane with approximately constant angular velocity and its charge distribution is characterized only by its dipole or quadrupole moment, one can use i. or j. in Tab I. If the particle is simultaneously oscillating and librating/rotating or for a general charge distribution, calculate the operators (17) and insert them into (22).
It is worth noting that for a flat surface-fluctuation spectrum, l(ω ) l(0), the resonant master equation describes the same dynamics as the slow-particle master equation if a rotating wave approximation with respect to the system's surface-free dynamics is justified.
The next subsection shows that the Markov and quasistatic approximations are justified for typical state-ofthe-art setups, ranging from trapped ions, to molecules and nanoparticles. We will then illustrate how to use the master equations in two examples, (i) a slowly rotating two-ion Coulomb crystal and (ii) a freely rotating polar molecule, both close to a surface covered by a thin dielectric layer. Finally, this section closes by discussing limiting situations in which the derived master equations reduce to known results.

B. Experimental adequacy of the description
Quasistatic approximation-The quasistatic approximation applies if any retardation in the wave-propagation from the source in the surface to the particle can be neglected. The experimentally observed distance scaling ∝ d −4 of the electric field noise [65,66,[69][70][71] indicates that the predominant noise sources are longitudinal fields originating from a thin surface layer. Retardation effects can be usually neglected given that the time it takes a photon to travel through the surface layer of thickness d s to the particle at distance d is much smaller than the characteristic timescale 1/ω 0 of the particle motion in the absence of the surface, ( |ε s |d s + d)ω 0 /c 1. Since |ε s | 40 and d s 5 nm for typical surface contaminants [89], retardation is negligible for frequencies up to ω 0 = 30 GHz at d = 1 mm, as applies to typical setups with atoms [2,65,66], molecules [7][8][9], and nanoparticles [25,34,[108][109][110][111][112]. In the absence of a thin surface layer, the quasistatic approximation is still justified as long as the skin depth δ is greater than the trapping height. For a Drude metal, one finds δ 2 = γ m c 2 /ω 0 ω 2 pm , yielding δ ≈ 4.3 mm for gold at kHz frequencies and δ ≈ 140 µm at MHz frequencies.
Born-Markov approximation-The weak-coupling Born-Markov approximation is justified if the temporal correlations in the surface fluctuations are irrelevant for the particle dynamics. Rigorously, this requires the surface-induced particle relaxation to take place on a timescale much greater than that of temporal correlations. The surface-noise correlation timescale is dominated either by the thermal correlation time /k B T or by the width of the spectral density (B10). For temperatures above 100 mK the thermal correlation time falls below 100 ps, while the correlation timescale due to e.g. the THz-resonance in Fig. 2 reaches 1 ps at most. Similar timescales are obtained for typical metals, though low-lying overdamped surface resonances may dominate (γ n /ω 2 n ≈ 10 µs for Fig. 2). These timescales can be compared to the particle-relaxation time, which lies in state-of-art setups in the range of up to milli-seconds. In addition, the rigorous Born-Markov approximation requires the surface-free particle dynamics to be either much slower than the decay of surface correlations (slow-particle limit) or much faster than particle-relaxation timescale (resonant limit). For instance, a rotating nanoparticle at 100 kHz is described by the slow-particle limit given the above surface properties, while the resonant limit applies for an oscillating ion at MHz frequencies. We emphasize that even if the Born-Markov approximation is not justified formally, it often yields an adequate description given that non-Markovian modifications of the dynamics are often marginal [113].

C. Example 1: Two-ion Coulomb crystal
We consider two equally charged ions with q = e, rotating with constant angular velocity ω rot at a constant interparticle distance d = 5 µm around their center of mass R · e 3 = 100 µm above a flat surface with surface normal e 3 , motivated by the experiment of Ref. [25]. Choosing the reference charge distribution such that the interparticle axis is aligned with the surface normal, one has It is characterized not only by the monopole and quadrupole moment, but also exhibits higher charge moments. The direction of the interparticle charge axis at a given crystal orientation Ω is given by n 3 = cos α sin βe 1 + sin α sin βe 2 + cos βe 3 . For a flat surface layer of thickness d s = 5 nm at T = 300 K with the exemplary response of Fig. 2, exhibiting spectrally flat and weak features in the kHz to GHz range on a bulk metal, the quasistatic approximation is justified. This is because (i) retardation between the surface and the particle is negligible on the timescale of the particle motion and (ii) because even the weak features with ε im ≈ 1 × 10 −5 dominate the total noise, i.e. noise from the gold bulk as well as blackbody radiation from infinitely distant surfaces, at the particle frequencies and its small distance above the surface. In addition, the surface-fluctuation correlations (B8) and (B9) decay on a timescale that is determined by the spectral widths of l(ω) and Im[ε(ω)]/|ε(ω)| 2 . In Fig. 2 the width of l(ω) is approximately given by ω 2 2 /γ 2 = 10 5 s −1 . The width of the energy loss function is approximately given by the width of the underdamped high frequency peak of γ 5 = 10 12 s −1 . The Markov approximation and Eq. (12) are justified since for rotation frequencies ω rot 100 kHz the crystal moves (and relaxes) much slower than those surface correlation timescales.
The resulting orientational decoherence rate (15) of the crystal with the surface-fluctuation kernel (52) is depicted in Fig. 3. We display the decoherence of a superposition between a fixed orientation (black dumbell and dots) and all other possible orientations. The rate as a function of the varying orientation is depicted in a contour plot in Mollweide projection. In superposition experiments, such as in [25], the varying orientation typically rotates on great circles (red).
Here we defined the operators n 0,i = n The coefficients n (i) m m are given in Eq. (38). Inserting (64) into (22) yields the dissipator with the surface-fluctuation kernel The particle is located at a distance of R cm · e 3 = 100 nm [7] from the superconducting surface, with normal vector e 3 , which is covered by a thin dielectric layer. We assume a layer thickness of d s = 4 nm and ε s (ω ) = 3 × (1 + 0.001i), as is realistic for e.g. an oxide layer [89,114]. This permittivity is of the same order of magnetitude as the response shown in Fig. 2. At these frequencies, its value is dominated by the high frequency peak rather than the low-frequency features. The quasistatic approximation of the electric fields is justified since (i) retardation between the surface and the particle is negligible and since (ii) the noise from the layer dominates over blackbody radiation and noise from a realistic metal or superconductor bulk. As discussed in the previous example, the surface-fluctuation correlations of the spectrum shown in Fig. 2 decay at a rate of 10 5 s −1 . Consequently, the Markov approximation is justified because the particle decoheres much slower than surface correlations decay and than the particle rotates, as we will show next.
In Fig. 4, we solve the master equation (8) with the dissipator (66), neglecting the coherent surface interaction. Using Eq. (50) we can approximate expressing that the noise is predominantly due to the surface layer. The initial state ρ(0) = |Ψ Ψ| generates a mean dipole moment that points along the 3 axis. In contrast, pure angular momentum states, or an incoherent mixture thereof, lead to a vanishing mean dipole. As displayed in Fig. 4, the decoherence rate crucially depends on the alignment between the dipole and the surface normal. The order of magnitude of the decoherence rate can be roughly estimated as p 2 i h ii (R cm , ω )/ ≈ 4 × 10 2 /s. This rate will increase for higher temperatures and a dielectric with a larger loss tangent, becoming relevant for hybrid quantum systems with polar molecules and superconducting devices [7,10]. . Rotational decoherence of a point dipole near a superconducting surface. (a) The particle is initially prepared in an angular momentum superposition state with initial dipole p 3 pointing perpendicular (blue) or parallel (red) to the surface normal vector e3. The noise originating from the surface is dominated by a thin layer covering the superconducting device (see Sec. V). The dynamics governed by (66) are displayed via the density operator in the basis of the lowest angular momentum states | m (quantization axis 3) for the times t = 0, t1 = 10 ms and t2 = 30 ms. Surface-induced dissipation relaxes the particle to its rotational ground state. (b)(top) While the groundstate population approaches unity for large times, the excited state population decreases at a constant rate. The latter differs in the two scenarios, marked in dashed blue and solid red, respectively. (b)(bottom) The coherence decay in (a) occurs also at constant rate. Decoherence in the red case is stronger, since the field fluctuations are stronger in the direction of the surface normal.

E. Relation to previous results
Electron beams near a metal surface-In the slowparticle limit the decoherence rate of a monopole reduces to that calculated in Ref. [88]. Specifically, the decoherence rate follows from Tab. Ia. and (14) together with the surface-fluctuation kernel of a Drude metal (46) as For superposed electron beams traveling parallel to a metal surface with a fixed distance between the superposition branches, this is identical to the path decoherence predicted in Ref. [88] and measured in [74].
Motional decoherence of ions in quadrupole traps-For a point monopole, the resonant master equation (22) reduces to the one used for describing heating of ions levitated in quadrupole traps. Specifically, using row d. in Tab. I at a trapping frequency ω k = ω 0 and oscillation direction k = 0 , yields the groundstate heating rate Γ h (ω 0 ) = q 2 n(ω 0 )h 0 (R eq )/mω 0 [71,105]. For k B T / ω 0 1, as is the case for typical trapping conditions [71], the heating rate can be approximated by means of the electric-field PSD (54) which exhibits the characteristic distance and frequency scaling discussed in Sec. IV. Rotational heating of polar molecules-Groundstate heating of dipolar rotors near surfaces has been studied in Ref. [87]. In order to verify that the two results agree, we evaluate the master equation (22) for a point dipole, neglecting all other multipole moments, at a fixed center-of-mass position with non-degenerate rotational states. Considering the quasistatic approximation (C5) in Ref. [87] shows that both theories coincide.

VI. DISCUSSION
Given the dielectric properties and geometry of the surface, the presented toolbox serves to predict quantitatively the induced heating and decoherence for a multitude of particle shapes and charge distributions. By disentangling the role of the particle from that of the surface, it provides novel tools for systematically studying the problem of surface-induced noise and decoherence in electric traps. In the following, we discuss a possible roadmap for future investigations, then comment on potential future extensions of the framework, and finally provide our conclusions.

A. Taming the surface noise
Characterizing the surface-As a main feature of the presented toolbox, the material properties of the surface enter only via the dielectric response functions of bulk and surface layer. Ex-situ measurements of the lowfrequency behavior of the dielectric material can therefore be used to predict its impact on the particle dynamics. Performing such measurements for different materials and surface contaminants will serve to identify the origin and strength of broad, low-frequency features in the dielectric permittivity.
Sensing field fluctuations-The levitated particle acts as a highly sensitive probe of electric surface noise since the heating force felt by a monopole and the torque experienced by a dipole are fully characterized by the electric field PSD (54). Specifically, the decoherence rate of a slow monopole at position R 0 in a superposition of size Such a quadratic dependence on ∆R is accompanied by momentum diffusion, and thus heating, with diffusion tensor q 2 S EE (R 0 , 0)/2. Likewise, for orientational superpositions of a slow point dipole in a superposition of orientations Ω and Ω , the decoherence rate reads with ∆p = p Ω − p Ω . Again, this decoherence is associated with angular momentum diffusion and thus heating [115]. Sensing gradient fluctuations-In more general situations, the electric field PSD does not suffice to describe the induced heating. For instance, the orientational decoherence and heating of a point quadrupole probes the PSD of the field gradient, a fourth-rank tensor, at its center-of-mass position. Likewise, the centerof-mass decoherence of a point dipole is determined by field gradient fluctuations characterized by the same tensor. For instance, the orientational decoherence of a point quadrupole at position R cm reads Γ Ω Ω (R cm , R cm ) = 1 36 Note that this equation exhibits a surface distance dependence distinct from the scaling in Eq. (56), showing that large-quadrupole particles might yield further insight into the origins of anomalous heating [66].
Probing non-local field correlations-It follows from Eq. (15) that the center-of-mass decoherence rate of a spatially delocalized wavepacket depends on non-local noise correlations, rather than on the local electricfield PSD. In particular, the decoherence rate remains bounded for large distances in the superposition, as observed in experiments with electron beams [74]. Creating and monitoring spatially extended superposition states thus enables deducing information on the surface noise that is not encoded in the heating rates, e.g. by utilizing modern parameter estimation techniques [116].
Modelling the dielectric response-Identifying the microscopic origins of surface noise will require input from atomistic theoretical models. This encompasses accounting for the role of microscopic processes as sources and as mediators of electromagnetic noise. Ab-initio calculations [97,98] and molecular-dynamics simulations [99] link microscopic models to the macroscopic dielectric response, which in turn serves as one input for the presented toolbox and can thus be tested experimentally.
Tailoring the surface response-Once the mechanism of surface noise is understood, surface-induced decoherence may be suppressed by choosing materials with beneficial dielectric properties. Specific surface coatings may help suppressing such noise at the relevant particle frequencies or prevent adsorption of surface contaminants. Plasma cleaning [117] and ion milling [118] of electrode surfaces in ion traps have shown that even slight surface modifications can have great impact on the noise level.
Geometry design-As a further handle to reduce the impact of unavoidable surface fluctuations, one may adjust the surface geometry to modify the electric fields at the particle position. For instance, when working with dipolar particles, one can minimize field noise, as compared to field gradient noise, by a mirror-symmetric arrangement of short-circuited electrodes.

B. Extending the toolbox
The presented toolbox is based on a number of assumptions, namely the quasistatic approximation of the surface-induced field, the Born-Markov approximation for the particle dynamics, the rigid description of the particle charge distribution, and the neglect of magnetic surface effects, all of which is justified for typical ion and nanoparticle trapping experiments. That said, the toolbox could be extended in following ways: Transverse electric field-Taking the transverse field into account requires two steps: First, the Green tensor, which describes both the longitudinal and the transverse field propagation, must be calculated instead of the scalar Green function. Second, the coupling between the particle velocity and the transverse fields must be accounted for. Note that the impact of Thomson scattering, which describes scattering from the purely transverse field at vanishing particle velocity, is discussed in App. D, where it is shown that its impact is negligibly small compared to decoherence due to the longitudinal field and other experimentally relevant decoherence processes.
Anisotropic and non-local dielectrics-To incorporate the possible anisotropy and non-locality of surface materials requires using the tensorial and non-local dielectric response in the definition of the surface polarization operators [100] and the scalar Green function. While the subsequent calculation is presumably tedious, we expect the derivation presented in this paper to carry over to the generalized case.
Non-rigid charge distributions-The fact that the bound charges inside the particle can fluctuate and react to fields emanating from the surface can be included in terms of the particle polarization. This adds to the dynamics (i) the Casimir-Polder interaction between the particle and the surface [91] as well as (ii) the associated decoherence channel [96]. This effect is expected to become relevant for weakly-charged micron-sized dielectrics.
Magnetic effects-Magnetic contributions to the particle-surface interaction can be taken into account by adding to the Hamiltonian the coupling of fluctuating surface magnetization currents to the charge distribution (see App. A 2). Likewise, if the particle exhibits permanent or fluctuating magnetic moments, their coupling to the surface polarization and magnetization fields can be included [91,95]. These extensions may serve to model the ro-translational heating of levitated magnets or superconductors close to surfaces [119][120][121][122], which are an alternative to electrically levitated particles [27,38].
Ro-translational thermalization-The dynamics of rotranslational damping are fully accounted for in the resonant limit. In the slow-particle limit, they can be described by including to first order the particle motion during the decay of surface-fluctuation correlations [101]. This adds to the master equation terms that are linear in the center-of-mass momentum and angular momentum operators. The resulting master equation describes ro-translational thermalization if it is of a specific form [123].
Non-Markovian effects-If the surface-fluctuation correlations decay on a timescale comparable to that of the particle motion, the weak-coupling Born-Markov approximation loses justification, even though it may still provide a good description for most situations. The impact of colored noise and related non-Markovian effects can be taken into account via the method of influence functionals [101] or by using quantum Langevin equations [124]. The latter feature an operator-valued noise to describe the dynamics of selected observables. If necessitated by empirical evidence, such a treatment is possible in principle, though significantly more involved and less versatile.

C. Conclusions
The framework described in this article predicts quantitatively surface-induced rotational and translational decoherence and heating of extended objects with arbitrary shape and charge distribution. The quantum master equations (12) and (22) incorporate the surface properties in terms of (i) the dielectric permittivity, which is empirically accessible and microscopically computable, and (ii) the scalar Green function, which describes the field propagation in arbitrary geometries. Our work applies to a variety of state-of-the-art experimental setups, ranging from trapped ion setups [66], to hybrid quantum devices with charged atoms and molecules [10], to nanoparticles in electric and optical traps [38]. We expect the presented framework and the worked-out toolbox Tab. I to be instrumental for a wide class of future quantum experiments with charged objects of finite extension. This appendix derives the coupling Hamiltonian between the charged particle and the surface. We first derive the minimal coupling Hamiltonian (1) from the Lorentz force and torque laws for a moving and revolving rigid charge distribution. Inserting into this Hamiltonian the electromagnetic fields in terms of surface excitations then yields the surface-particle interaction Hamiltonian, which serves as a starting point for the derivation of the master equations.

Rotor-field coupling
To describe the classical dynamics of a particle in the electromagnetic field, we introduce the position R and the orientation via Euler angles Ω = {α, β, γ} in the z-yz convention. The particle charge distribution at a given point r reads (r) = 0 R T Ω (r − R) . Here, 0 describes the charge distribution in the reference position and the rotation tensor R Ω relates the space-fixed axes e i to the principal axes n i = R Ω e i . The angular velocity is defined via the kinematic relationṘ Ω = ω × R Ω . The coupling Lagrangian, involves the charge distribution and the current density, due to the rigid motion of the body. We first show that this Lagrangian is consistent with the Lorentz force and torque laws. To demonstrate this equivalence, we use that ω =αe 3 +βe ξ +γn 3 (A3) with e ξ = − sin αe 1 +cos αe 2 . We employ (A3) to express the kinetic energy E kin and the current density in terms of the Euler angles. We then transform the integrals in (A1) to the body-fixed frame yielding The Lagrange equation for the center of mass involves the total time derivative of the vector potential which describes the center-of-mass dynamics of an extended charge distribution in presence of the electric and magnetic fields The Lagrange equations for the Euler angles are more involved. First, we deduce some useful relations from the definition of the rotation matrix and (A3), with µ ∈ {α, β, γ} and e µ ∈ {e 3 , e ξ , n 3 }. Note that de 3 /dt = 0 since the space-fixed axis does not move. Let us start with the angle α, and use (A10a) to write d dt ∂ ∂α Since I 0 = R T Ω I Ω R Ω is independent of the Euler angles, we find ∂ ∂α where we used (A10c). Combining (A12) and (A13), one thus obtains the e 3 -component of the angular-momentum time-derivative The term involving the scalar potential can be rewritten using (A10b), again with R 0 = R + R Ω r. Likewise, we evaluate the term involving the vector potential. In a first step, we compute where we used (A10a). Noting Eq. (A16) yields Subtracting the α-derivative gives With the help of (A10b) one finds that where (A6) and the Graßmann identity have been used. We combine (A14), (A15) and (A21) with the Lagrange equation for α to finally obtain This is the e 3 -component of the Lorentz torque acting on an extended charge distribution. Since (A10a)-(A10c) equally hold for β and γ, the corresponding Lagrange equations yield the analogue of (A22), but with e 3 replaced by e ξ and n 3 , respectively. Thus, we demonstrated that the Lagrangian (A1) is indeed consistent with the Lorentz torque acting on an extended charge distribution, In order to obtain the Hamiltonian we first define the canonical linear and angular momenta and then perform a Legendre transformation. Here we used (A10a) and related the angular momentum vector to the canonical Euler momenta via p α = e 3 ·J, p β = e ξ ·J and p γ = n 3 · J. The Legendre transformation yields the Hamiltonian (1).

Macroscopic quantum electrodynamics
One central ingredient of macroscopic quantum electrodynamics is the Green tensor G(r, r , ω), giving the electromagnetic field at the position r due to a source current located at r . In the presence of dispersive media with permittivity ε r (r, ω) and permeability µ r (r, ω), it fulfills the defining equation along with the boundary condition G(r, r , ω) = 0 for |r − r | → ∞. The electric field due to the current density j(r, ω) is then given by If the current density obeys the fluctuation-dissipation theorem, one obtains [91], The theory is quantized by expressing the current density in terms of quantum harmonic oscillators [91,100], the so-called polarization and magnetization oscillators f λ (r, ω), λ ∈ {e, m}. They obey the canonical commutation relations [f λ (r, ω), f λ (r , ω )] = 0 and [f λ (r, ω), f † λ (r , ω )] = 1δ λ,λ δ(r − r )δ(ω − ω ). These oscillators are related to the (noise) polarization and magnetization of the medium through so that the associated noise current density acts as a source for the electromagnetic field. The free field dynamics are governed by the Hamiltonian In Coulomb gauge the electromagnetic potentials originating from the dielectric can be expressed through the (Schrödinger-picture) polarization and magnetization oscillators as [91,100] Here is the left-transverse part of the tensors G λ , which are related to the Green tensor, The transverse delta function is defined as Expressing the Green tensor in terms of the scalar Green function g(r, r , ω) (see App. C) shows that the vector potential vanishes, while the scalar potential is of the form Eq. (7). In the following, we will formulate the theory in terms of the Green tensor to keep the discussion general. Simplifications associated with using the scalar Green function g are discussed along the way.
The quantized electric field (A29) thus reads (A41) We note that in presence of a particle, the field (A41) is not yet the physical electric field, since the Coulomb field of the particle charge in absence of the dielectric medium has still to be added [91]. The dynamics generated by the minimal coupling Hamiltonian (1) and the free-field dynamics (A34) are consistent with the Maxwell equations, the Lorentz force and torque, and thermal field fluctuations in absorbing media.
In the following, we neglect the vector potential to get the total Hamiltonian where the free Hamiltonian involves the external potential V ext (R, Ω) and the particle-field interaction reads as Here R, P, Ω = {α, β, γ}, p α = e 3 · J, p β = e ξ · J and p γ = n 3 · J are operators satisfying the canonical commutation relations of translation and rotation.

Appendix B: Born-Markov approximation
This section derives the two main master equations from the above discussed interaction Hamiltonian. In a first step, we will perform the weak-coupling approximation, and then evaluate the resulting Markovian master equation in the slow-particle limit and in the resonant limit.

Weak-coupling approximation
In the weak-coupling approximation, the time evolution of the interaction-picture density operatorρ(t) = can be written as [101] ∂ρ ∂t where the bath density operator ρ B determines the state of the field degrees-of-freedom, which are assumed to be uncorrelated with the particle state, while the bath correlation (surface-fluctuation correlation) is assumed to decay much faster than the particle-relaxation time. The interaction Hamiltonian (A44) in the interaction picture decomposes as where we defined the Coulomb field operator We assume the bath state ρ B to be thermal with respect to the free field Hamiltonian, implying that f λ (r , ω) B = 0, where · B = tr B {·ρ B }. The second moments of the bath oscillators are finite, where we introduced the Bose-Einstein distribution n(ω) = 1/[exp( ω/k B T ) − 1]. The second moments determine the bath-operator autocorrelations, These autocorrelations contain the tensorial dissipation kernel [101] D(s, s , τ ) = 2 (4π) 2 ∞ 0 dω J(s, s , ω)sin(ωτ ) (B8) and the tensorial noise kernel Both are determined by the spectral density as follows by using the relation [91], In addition, the Onsager reciprocity [91] of the Green tensor, also applies to D, N and J, while the Schwarz reflection principle of G leads to In summary, Eq. (B1) can be rewritten as This non-Markovian evolution equation serves as a starting point to derive the Markovian master equations in two different limiting cases: the slow-particle limit and the resonant limit.

Slow-particle limit
If the free particle motion is slow, so that the autocorrelations D and N drop on a much shorter timescale than the timescale over which the Coulomb field operator C(s , t) changes, one can approximate In this case, the time-integration in (B14) can be carried out, yielding Transforming equation (B14) back to the Schrödinger picture one obtains where we introduced the coherent surface interaction (Lamb shift) The frequency integral included in (B19) can be rewritten as an integral over the whole frequency axis, (B20) where P denotes the Cauchy principal value. Using contour integration, exploiting that the Green tensor is holomorphic in the upper complex frequency plane, that it approaches G(s, s , ω) ∼ −c 2 1δ(s−s )/ω 2 for large ω [91], and that it is longitudinal with (C5) for ω towards zero, the integral can be evaluated Consequently, the surface-interaction Hamiltonian contains the negative electrostatic self energy of the rigid particle charge distribution and an additional contribution due to the electrostatic interaction of particle and the dielectric, The self energy, given by the first term, is independent of the particle operators R and R Ω and thus does not have an impact on the particle dynamics in the slow-particle limit.
The second term gives the image-charge interaction. It is the potential energy of a classical charge distribution due to the presence of a dielectric surface. For instance, a monopole with 0 (r) = qδ(r − R) in front of a Drudemetal half space (see Tab. III) with permittivity (43), experiences the image potential The master equation (B18) can be brought into Lindblad form by using the relation (B11), In the quasistatic approximation, the Green tensor is longitudinal and takes the form (C5). As a result, L m (r) = 0, so that the master equation simplifies to with where we performed a partial integration. The charge distribution enters through as a direct consequence of the definition (B3). Having shown that the master equation has Lindblad form, we can again use (6) to simplify the master equation in the quasistatic limit, yielding Eq. (12). This equation may also be obtained directly by inserting the quasistatic Green tensor (C5) into (B18) and subsequently carrying out partial integrations.

Resonant limit
A Markovian master equation can also be derived in the limit that the free particle dynamics are much faster than the bath-induced particle-relaxation time. We start from the non-Markovian evolution equation (B14) and expand the Coulomb field operator (B3) in eigenoperators C (s) of the free evolution, The index covers all positive energy differences ω . We express through the eigenstates |Ψ n of H 0 with energies E n and coefficients (B33) The operators c (s) in (17) expand the charge distribution, and are thus related to (B30) via where the coefficients read We define the Fourier integrals of the damping and noise kernel, and in order to write Here we assumed that the time-independent term in (B30), C 0 (s), is a c-number, as it is the case in all situations discussed in the main text. The case of an operatorvalued C 0 (s) could be included along the lines of the derivation in the slow-particle limit in Sec. B 2.
With (B39) we find Γ N (s, s , −ω )C † (s )e iω t ,ρ(t) . (B40) If the free particle dynamics are fast, so that all frequencies and frequency differences are much larger than the rate at which the density operatorρ(t) relaxes, we can carry out the rotating wave approximation and cancel all rapidly oscillating terms. This yields after a rearrangement that uses the symmetries in ω, s and s . Here describes the coherent surface interaction in the resonant limit. The two principal-value integrals can be evaluated using contour integration. We express the imaginary part of the Green tensor through the total Green tensor and use the Schwarz reflection principle to then exploit that the Green tensor is holomorphic in the positive imaginary plane. Accounting for the pole at ω = ω , all the poles of n(ω) in the positive imaginary plane and at ω = 0, and for the behaviour of the integrand at infinity, the two integrals yield and (ω th m) 2 ω 2 + (ω th m) 2 G(s, s , iω th m), (B45b) Here ω th = 2πk B T / denotes the thermal Matsubara frequency. In the limit T → 0 the Matsubara frequency summation transforms into an ordinary integral that is independent of temperature. In the quasistatic limit (C5) one then obtains (21).
Similar to the the slow-particle limit we now briefly give the result for a particle in front of a perfectly conducting half space. The quasistatic Green tensor (C5) in front of a perfectly conducting half space (see Tab. III) reads yielding for (B45a) and (B45b) temperature-independent expressions, The surface-interaction Hamiltonian (B44) can now be evaluated. The transverse delta-function does not contribute and one finds after a partial integration by use of (B35). A smallamplitude oscillating monopole with (27) and a single frequency ω 0 , oscillation direction 0 and annihilation operator a 0 is consequently exposed to the potential Here e 3 is the surface normal. This quantum potential could also be derived by taking the classical image-charge potential of an extended charge distribution, promoting the charge distribution to an operator, accounting for the small-amplitude oscillations of the monopole and carrying out a rotating wave approximation.
To simplify the Lindblad master equation (B53) to a scalar form we can insert the quasistatic Green tensor (C5), Here the Lindblad operators do not depend on the C (s), but instead involve the operators c (s) used to expand the charge distribution. We use relation (6) to finally arrive at (22). Similar to the slow-particle limit, this equation may also be obtained directly by inserting the quasistatic Green tensor (C5) into (B43) and subsequently carrying out partial integrations.
Appendix C: Quasistatic approximation of the Green tensor Depending on the material response ε r (r, ω), the frequency, and the distance |r − r |, retardation effects in the propagation of the electromagnetic field to the particle can be neglected. In this quasistatic limit, the Fourier transformed electric field E(r, ω), generated by an arbitrary charge distribution ρ(r, ω) follows as the solution of the quasistatic Maxwell equations The charge distribution is related to the current density by the continuity equation Expressing the electric field in terms of the Green tensor and using (A29), the quasistatic Maxwell equations (C2) can be rewritten as They are solved by where the scalar Green function g(r, r , ω) is the solution of ∂ ∂r · ε 0 ε r (r, ω) ∂ ∂r g(r, r , ω) = −δ(r − r ).
Comparing Eqs. (C6) and (A27) shows that the complexity of the problem was reduced considerably. From the defining equation (C6) and the Schwarz reflection principle for the dielectric response follows the Schwarz reflection principle for the Green function g * (r, r , ω) = g(r, r , −ω * ). The boundary conditions of the problem can be deduced from inserting (C5) into the electric field (A29) and carrying out a partial integration, To ensure that the electric field vanishes infinitely far from the source, one requires that ∂g(r, r , ω)/∂r = 0 for |r − r | → ∞. Equation (C6) determines g(r, r , ω) up to a constant, so that we can set g(r, r , ω) = 0 for |r − r | → ∞, as is beneficial for partial integration. From Eq. (C7) we can see that g(r, r , ω) yields the quasielectrostatic potential at r of a test charge located at r oscillating in magnitude with ω. Inserting (C5) into the scalar potential (A36) and a subsequent partial integration yields, after transforming to the interaction picture, the simplified form (7).

Continuity at an interface
At an interface ε r (r, ω) changes discontinuously, so that the defining equation (C6) implies continuity conditions. These conditions will be used in Sec. C 4 and Sec. C 5 to derive the scalar Green function of the dielectric half-space and the dielectric layer configuration.
For example, a flat interface between two homogeneous dielectrics with ε r (r, ω) = ε < for r · e 3 < 0 and ε r (r, ω) = ε > for r · e 3 > 0, volume integration at the boundary gives the condition where r > and r < approach the surface from above and below. This expresses the continuity of the displacement field's normal component. Likewise, by integrating over a loop through the interface we find the continuity for the electric field's parallel component e · ∂ ∂r g(r > , r , ω) = e · ∂ ∂r g(r < , r , ω), where e is perpendicular to e 3 .

Symmetry of the Green function
To prove the relation (6), we first show that g(r, r , ω) = g(r , r, ω). (C9 )  Table III. The Green function for a dielectric half-space with εr(ω) for r · e3 < 0 and vacuum for r · e3 > 0.
For the half-space geometry depicted in Tab. II, where we assume vacuum for r · e 3 > 0, while for r · e 3 < 0 we set ε r (r, ω) = ε r (ω) to be the response of an arbitrary homogeneous dielectric. With r · e 3 > 0 in vacuum the solution can be obtained from the ansatz with the mirror tensor M = 1−2e 3 ⊗e 3 . The second term in Eq. (C20a) is a homogeneous solution, since e 3 · Mr < 0. Continuity at the interface (C8) implies .
The homogeneous solution can be interpreted as the potential due to the presence of an image charge. The same calculation for r ·e 3 < 0 in the dielectric yields the entire Green function given in Tab. III.

Half-space covered by a surface layer
We next derive the Green function g(r, r , ω) in a geometry consisting of a vacuum half-space with ε r (r, ω) = 1 for r · e 3 > 0, a dielectric layer with ε r (r, ω) = ε s (ω) for −d s < r · e 3 < 0, and a dielectric half-space with ε r (r, ω) = ε b (ω) for r · e 3 < −d s , see Tab. II.
In a geometry with two interfaces, an infinite number of image charges appears. We will now use a compact formulation by expanding the particular solution (C19) in exponential functions [125], with the in-plane distance ∆r 2 = [(r − r ) · e 1 ] 2 + [(r − r ) · e 2 ] 2 . One can show that ∂ ∂r · ∂ ∂r e ±kr·e3 J 0 (k∆r) = 0, so that we can make the ansatz g(r, r , ω) = 1 4πε 0 ∞ 0 dk g k (r, r , ω)J 0 (k∆r) , (C24) with the g k (r, r , ω) given in Tab. IV. Then (C24) solves the defining differential equation (C6) in all regions of constant ε r (r, ω) and accounts for the boundary condition at infinity and the symmetry (C9). The coefficients c 1 to c 10 appearing in Tab. IV can be determined by using the continuity at the interface (C8a) and (C8b). Defining the coefficients read some of which can also be found in [89,125]. Expanding the denominator in (C26) using the geometric series yields an infinite sum of image potentials for the layer configuration. Table IV. The Green function for a geometry consisting of a vacuum half-space with εr(r, ω) = 1 for 0 < r · e3, a dielectric layer with εr(r, ω) = εs(ω) for −ds < r · e3 < 0 and a dielectric half-space with εr(r, ω) = ε b (ω) for r · e3 < −ds. For the table we defined z = r · e3 and z = r · e3. For the dependence of the Green function on g k (r, r , ω), see Eq. (C24).
In order to calculate the decoherence rate, one has to determine the dissipator L Th from the double commutator in (D7). In the slow particle limit and after transforming back to the Schrödinger picture we find where C i1i 1 (s, s ) := C i1i 1 (s, s , 0). The Born-Markov approximation assumed here is justified if the particle operator C(s, s , t) evolves slowly in comparison to the timescale on which the dissipation kernel decays.
The integral of the noise kernel without the time-independent contribution yields Here Tr [·] is the tensor trace. The dissipator L Th describes decoherence of a particle with arbitrary charge distribution near a dielectric medium due to the A 2 -term in the coupling Hamiltonian. The properties of the dielectric enter through the spectral density J(s, s , ω), while the particle properties enter through the operator C(s, s ). with 1 , 2 forming an orthonormal basis with n (like the two polarization vectors of a plane wave traveling in direction n ). With these definitions and (D3) one can write Here, one can identify the orientation-dependent scatting amplitude of Thomson scattering of an electromagnetic wave with incoming wave vector k and polarization j , elastically scattered at a rigid charge distribution to an outgoing direction with wave vector k, as can be derived from classical electrodynamics [127]. With this (D27) can be simplified to L j (kn, kn ) · ρ L † j (kn, kn ) − 1 2 L † j (kn, kn ) · L j (kn, kn ), ρ . (D32) The dissipator (D32) is of the typical form to describe scattering decoherence. The Lindblad operators are given by plane waves multiplied by the scattering amplitude, while the pre-factor accounts for the statistics of the scattering incidents. The dissipator is consistent with the dissipator derived in the Monitoring Approach [126,128].

Thomson-scattering decoherence of a point monopole
To estimate decoherence of a point monopole 0 (r) = qδ(r), we use that the transverse electromagnetic field, and thus the vector potential, dominate the particlesurface interaction for distances larger than the wavelength, e.g. the thermal wavelength 2π c/k B T . For simplicity, we thus assume the particle to be located in free space, where (D32) applies. For the point particle, the scattering amplitude (D28) reduces to which fulfills × d 2 nd 2 n 1 + (n · n ) 2 1 − e ik(n −n)·(R−R ) .

(D35)
We express the spherical integrals through the spherical Bessel functions of the first kind j n (x), where x = k|R − R |. The rate vanishes for |R − R | = 0, while for separations |R − R | much larger than the thermal wavelength it saturates at the value (D37) The impact of decoherence due to the A 2 -term can now be evaluated explicitly. For example a single electron in vacuum at temperature T = 300 K is subject to a decoherence rate of Γ ∞ ≈ 1.5 × 10 −5 s −1 , which is negligible in comparison to the decoherence due to the longitudinal field. This is a consequence of the tiny scattering cross section of Thomson scattering, which is even smaller for nanoparticles with their reduced ratio of q 2 /m.