Non-contact friction in ultracoherent nanomechanical resonators near dielectric materials

NCF for an infinitesimally thin string

We calculate the NCF damping coefficient for a one-dimensional tensioned string interacting with a substrate. The schematic is shown in Extended Data Fig. 1a. The string has length L and a linear charge density of ρe,l. The motion of the beam is given by a mode shape U(x) (motion along the z axis). Generally, for situations where the electric energy dissipated within the volume of the resonator itself can be neglected, g is equivalent to the electric field generated by a polarization density \({\bf{P}}({\bf{r}})={\rho }_{e}({\bf{r}}){{\bf{U}}}_{m}({\bf{r}})\) (Supplementary Information). Due to the zero volume of the string, to calculate the field variation g we can use this approximation and replace the linear charge density, with a linear polarization density \({{\bf{P}}}_{{\rm{l}}}(x)={\rho }_{e,{\rm{l}}}U(x){{\widehat{\mathbf{z}}}}\). Moreover, for modes where \(d\ll | {U}^{{\prime} }(x){| }^{-1}\) we can approximately reduce the problem to two dimensions in the cross-section and compute the field for an infinite line (Supplementary Information). We consider two different types of substrates and compute the NCF-limited quality factors for each case.

NCF due to a conductive substrate with a thin oxide layer

Even for high resistivity silicon, the characteristic time constant for the charges to decay is given by ϵ/σ ≈ 10−8 s, rendering silicon as a conductor for our frequency range. Hence, we have to treat the silicon substrate as a conductor in our models. For modelling the effect of the native silicon oxide layer, the schematic is shown in Extended Data Fig. 1b. For a thin oxide layer (Δ ≪ d), we can neglect the effect of the dielectric layer on the electric field and compute it as if the string was in front of a conductor. Hence, the electric field variation is given by the field generated by the original dipole density \({{\bf{P}}}_{{\rm{l}}}(x)\) in addition to its image \(-{{\bf{P}}}_{{\rm{l}}}(x)\) at z = Δ + d. Moreover, we approximate the field within the layer to be uniform and equal to the field generated by the two dipoles at the boundary, reduced by \({\epsilon }_{r}^{{\prime} }\). The resulting field within the dielectric only has a z component and is given by

$${\bf{g}}({\bf{r}},\omega ;U(x))\approx \frac{{\rho }_{e,{\rm{l}}}}{{\uppi }{\epsilon }_{0}{\epsilon }_{r}}\frac{\left({y}^{2}-{d}^{\,2}\right){{\widehat{\mathbf{z}}}}}{{\left({y}^{2}+{d}^{2}\right)}^{2}}U(x).$$

(6)

To calculate γNCF, the NCF damping for a harmonic mode of the string with frequency Ωm, we use equation (2) and integrate \(| {\bf{g}}{| }^{2}\) within the oxide layer

$${\gamma }_{\mathrm{NCF}}=\frac{{\rho }_{e,{\rm{l}}}^{2}{\epsilon }_{r}^{{\prime\prime} }}{| {\uppi }{\epsilon }_{0}{\epsilon }_{r}{| }^{2}{\Omega }_{m}}\times {\int }_{0}^{\varDelta }{\int }_{-\infty }^{\infty }{\int }_{0}^{L}{\rm{d}}x{\rm{d}}y{\rm{d}}z\frac{{\left({y}^{2}-{d}^{\,2}\right)}^{2}}{{\left({y}^{2}+{d}^{2}\right)}^{4}}| U(x){| }^{2}.$$

(7)

The integration along the x axis is evaluated using the definition of the effective mass. Given that the mode Un(x) is normalized to its maximum the effective mass is defined as \({m}_{\mathrm{eff}}={\rho }_{m,{\rm{l}}}{\int }_{0}^{L}{\rm{d}}x| U(x){| }^{2}\), where ρm,l is the linear mass density. Carrying out the integration in the yz plain yields

$${\gamma }_{\mathrm{NCF}}=\frac{{m}_{\mathrm{eff}}{\rho }_{e,{\rm{l}}}^{2}\varDelta }{4{\uppi }{\epsilon }_{0}{d}^{3}{\rho }_{m,{\rm{l}}}}\frac{\tan \delta ({\varOmega }_{m})}{{\epsilon }_{r}^{{\prime} }({\varOmega }_{m}){\varOmega }_{m}}.$$

(8)

Calculating the damping rate, ΓNCF = γNCF/meff, and the quality factor yields

$${Q}_{\mathrm{NCF}}=\frac{4{\uppi }{\epsilon }_{0}{d}^{\,3}{\rho }_{m,{\rm{l}}}}{{\rho }_{e,{\rm{l}}}^{2}\varDelta }\frac{{\epsilon }_{r}^{{\prime} }({\varOmega }_{m}){\varOmega }_{m}^{2}}{\tan \delta ({\varOmega }_{m})}.$$

(9)

To obtain the numerical estimates for QNCF shown in Fig. 1e, we have considered the dielectric permittivity and loss tangent of the SiO2 substrate to be 3.9 and 10−3, respectively40. We have also retrieved the designs and the chips from refs. 4,11,12, measured the distance d using an optical microscope and computed the corresponding QNCF for each geometry.

NCF due to a semi-infinite dielectric substrate

The schematic for the semi-infinite lossy dielectric substrate is shown in Extended Data Fig. 1c. The z > 0 space is filled by the dielectric with relative permittivity ϵr(ω). Using an elementary method of image charges, we calculate the field within the dielectric (z > 0) by replacing \({{\bf{P}}}_{l}(x)\) with \((1+\xi ){{\bf{P}}}_{l}(x)\). Following these steps, we find the electric field variation for z > 0 as

$${\bf{g}}({\bf{r}},\omega ;U(x))=\frac{(1+\xi \,){\rho }_{e,{\rm{l}}}}{2{\uppi }{\epsilon }_{0}{\epsilon }_{r}}\times \frac{2y(z+d){{\widehat{\mathbf{y}}}}+\left[{(z+d)}^{2}-{y}^{2}\right]{{\widehat{\mathbf{z}}}}}{{\left[{y}^{2}+{(z+d)}^{2}\right]}^{2}}U(x).$$

(10)

To calculate γNCF,n, the NCF damping for a harmonic mode of the string Un, we use equation (2).

$${\gamma }_{\mathrm{NCF},n}=\frac{| 1+\xi {| }^{2}{\rho }_{e,{\rm{l}}}^{2}}{4{{\uppi }}^{2}{\epsilon }_{0}| {\epsilon }_{r}{| }^{2}}\frac{{\epsilon }_{r}^{{\prime\prime} }}{{\varOmega }_{n}}\times {\int }_{0}^{\infty }{\int }_{-\infty }^{\infty }{\int }_{0}^{L}{\rm{d}}x{\rm{d}}y{\rm{d}}z\frac{| {U}_{n}(x){| }^{2}}{{\left[{y}^{2}+{(z+d)}^{2}\right]}^{2}}.$$

(11)

Carrying out the integrations similar to the previous case, one obtains

$${\gamma }_{\mathrm{NCF},n}=\frac{{m}_{\mathrm{eff},n}{\rho }_{e,{\rm{l}}}^{2}}{8{\uppi }{\epsilon }_{0}{d}^{2}{\rho }_{m,{\rm{l}}}}\frac{\mathrm{Im}[\xi ({\varOmega }_{n})]}{{\varOmega }_{n}},$$

(12)

where we have used the relation \(\mathrm{Im}[\xi ({\varOmega }_{n})]\approx \frac{2{\epsilon }_{r}^{{\prime\prime} }}{{(1+{\epsilon }_{r}^{{\prime} })}^{2}}\) with the assumption \({\epsilon }_{r}^{{\prime\prime} }\ll {\epsilon }_{r}^{{\prime} }\). Calculating the damping rate, ΓNCF,n = γNCF,n/meff,n, and the quality factor yields equation (3).

Measurement and calibration techniques

The measurement apparatus used to perform the ringdown and thermal noise measurements is the same as the one used in refs. 4,12. For completeness, we show a simplified schematic of the setup in Extended Data Fig. 2a. The sample is placed in a vacuum chamber with P < 10−8 mbar and room temperature. The mechanical motion is read out using a balanced homodyne interferometer, where laser light from a 780-nm diode laser is focused on the sample using a microscope objective. The reflected light is then interfered with a local oscillator and detected using a balanced photodetector. The optical path length difference and conequently the local oscillator–signal phase difference is stabilized by feeding back the low frequency (that is, <300 Hz) component of the homodyne signal to cascaded double acoustic–optic modulators (AOMs) in the local oscillator path, using a technique described in ref. 12. The high frequency component contains the mechanical motion signal and is used for the ringdown and thermal noise measurements. Using a lock-in amplifier unit (Zurich Instruments UHFLI), the signal is demodulated at the mechanical frequency and the amplitude of motion (that is, R in Fig. 3d) is obtained which is used for the ringdown measurements. The mechanical mode is excited using a piezoelectric actuator attached to the sample holder. To avoid optical back-action effects, the ringdowns are done stroboscopically using a shutter placed in the optical path. In our measurements the typical open to close duty cycle of the shutter is about 1:10. The calibrated thermal force measurement, shown in Fig. 3e, is performed using the technique described in ref. 4. The power spectral density of the homodyne signal is obtained using a spectrum analyser and it is calibrated in units of displacement spectral density, that is, m2 Hz−1 using the method described in ref. 4. To avoid the distortion of the spectrum due to slow drifts in the mechanical frequency, the mechanical mode is feedback cooled using the piezoelectric actuator by broadenning the mechanical linewidth from its intrinsic value of Γm to a value of Γeff. The thermomechanical spectrum shown in Extended Data Fig. 2b is acquired under this condition. The measured spectrum is given by \({S}_{x}^{\mathrm{meas}}(\omega )={S}_{{\rm{F}}}^{\mathrm{th}}| {\chi }_{\mathrm{eff}}(\omega ){| }^{2}+{S}_{x}^{\mathrm{imp}}\), where

$${\chi }_{\mathrm{eff}}(\omega )=\frac{1}{{m}_{\mathrm{eff}}\left({\varOmega }_{m}^{2}-{\omega }^{2}-i\omega {\varGamma }_{\mathrm{eff}}\right)},$$

(13)

is the effective mechanical susceptibility under feedback cooling, which is fully known since the effective mass is also known from the FEM simulations, and Γeff is obtained from a Lorentzian fit to the measured spectrum. \({S}_{x}^{{\rm{imp}}}\) is the imprecision noise floor of the measurement. We obtain the measured thermal force spectrum by deviding the measured displacement spectrum by the effective susceptibility, which yields \({S}_{{\rm{F}}}^{\mathrm{meas}}={S}_{{\rm{F}}}^{\mathrm{th}}+{S}_{{\rm{F}}}^{\mathrm{imp}}\). The resulting spectra are shown in Extended Data Figs. 2c and 3e, in all of which the flat component corresponds to \({S}_{{\rm{F}}}^{\mathrm{th}}\), and the rising branch corresponds to \({S}_{{\rm{F}}}^{\mathrm{imp}}\).

Ruling out other damping mechanisms

We discuss five mechanisms that could explain a distance-dependent damping in both the uniform strings and the nano-optomechanical system. We show that none of them can fully explain our observations.

Gas and squeeze-film damping

Although the samples are tested under high vacuum (P < 10−8 mbar), the small gaps can result in enhancement of gas damping due to squeeze-film damping. However, squeeze-film gas damping results in a frequency-independent damping rate41,42, different from our observations.

Local surface contamination

Local surface contamination means reduction of the intrinsic mechanical Q in a local area on the nano-beam. For example, this could be a result of the etching of the PhCs in the vicinity of the beam. Such a mechanism could in principle explain the different behaviour of modes with different parities. To study this effect, we use a modified version of the FEM simulations used for estimating the dilution factor for tensioned and high-aspect-ratio resonators43. These simulations are based on evaluating the kinetic and bending energies, Wkin and Wbend, for a given mode and then computing the quality factor as Q = Wkin/Wbend. For a two-dimensional resonator (in the xy plane) and an OP mode with frequency Ω and mode shape U(x, y), these energies are given by

$${W}_{\mathrm{kin}}=\frac{1}{2}{\rho }_{m}h{\varOmega }^{2}{\int }_{S}{\rm{d}}x{\rm{d}}yU{(x,y)}^{2},$$

(14)

$$\begin{array}{l}{W}_{\mathrm{bend}}={\frac{E{h}^{3}}{24(1-{\nu }^{2})}}\displaystyle \int _{S}{\rm{d}}x{\rm{d}}y\phi (x,y)\\ \times \left[{({U}_{xx}+{U}_{yy})}^{2}+2(1-\nu )\left({U}_{xy}^{2}-{U}_{xx}{U}_{yy}\right)\right],\end{array}$$

(15)

where ρm is the mass density, h is the thickness, E is the Young’s modulus, ν is the Poisson’s ratio of the resonator and Uij corresponds to partial derivatives of U(x,y). The integrations are carried over S, the area of the resonator. ϕ(x, y) defines the local mechanical loss angle of the resonator where in the usual cases, it is constant and given by the inverse of the material’s intrinsic quality factor \(\phi (x,y)={Q}_{{\rm{int}}}^{-1}\).

For the resonator design used to take the data in Fig. 3g, we use this FEM simulation to estimate the effect of local contamination. In a region on the resonator nearby the PhC’s location (Extended Data Fig. 3a), we set ϕ higher than the usual value so that we reproduce the observed reduced quality factor. We compute the Q for the first few OP modes and compute the added damping coefficient compared to the fully clean case (similar to equation (5)). We show these values for the modes with an anti-node in the centre in Extended Data Fig. 3d. We observe an increasing frequency dependence for the added damping coefficient which is unable to explain the experimental observations.

Mechanical coupling to the PhC

Another mechanism that can introduce a distance-dependent loss in the nano- optomechanical system is coupling of the high-Q mode to the low-Q mechanical modes of the PhC cavity22. Extended Data Fig. 3b shows an optical microscope image of the PhC cavities and the string suspended between them. One can clearly see that the cavities are also suspended and can support mechanical modes. We show an FEM simulation of the fundamental mode of the PhC cavity in Extended Data Fig. 3c. This mode has the frequency of 2.1 MHz and we also observe it in the experiments. It can in principle couple to the modes of the binary-tree resonator through multiple mechanisms such as electrostatic or van der Waals forces. Generally, for two mechanical modes with frequencies Ω1(2), damping rates Γ1(2) and effective masses m1(2), that are coupled by a spring constant kc, the additional damping for mode 1 is given by

$$\delta {\varGamma }_{1}=\frac{{k}_{{\rm{c}}}^{2}}{{m}_{1}{m}_{2}}\frac{{\varGamma }_{2}}{{\left({\varOmega }_{2}^{2}-{\varOmega }_{1}^{2}\right)}^{2}+{\varGamma }_{2}^{2}{\varOmega }_{1}^{2}},$$

(16)

having a resonance-like frequency dependence. For the measured parameters for the two mechanical modes, we set kc to a value to reproduce the additional damping for the fundamental mode of the binary-tree resonator at 160 kHz. The result is shown in Extended Data Fig. 3d, clearly inconsistent with our observed frequency scaling.

Finite electrical conductivity

Throughout our work, we have assumed that the dissipation within the dielectrics (SiO2 and Si3N4) is purely due to the dielectric loss (that is the imaginary part of the permittivity). However, having a finite charge mobility and consequently electrical conductivity of σ ≠ 0, can also introduce NCF. For a conductivity, small enough that does not disturb the electric field distribution at the mechanical frequency (that is, σ ≪ ϵ0ϵr/Ωm) it can be shown that the NCF coefficient arising from conductive losses is given by (Supplementary Information)

$${\gamma }_{\mathrm{NCF}}^{\mathrm{cond}}=\int \,{d}^{\,3}r\,| {\bf{g}}({\bf{r}},{{\varOmega }_{m}};{{\bf{U}}}_{m}({\bf{r}})){| }^{2}\frac{\sigma ({\bf{r}},{\varOmega }_{m})}{{\varOmega }_{m}^{2}}.$$

(17)

For a frequency-independent conductivity, evidently this form of NCF cannot explain our observed frequency scaling of 1/Ωm. The only possible way in which the conductive losses give rise to 1/Ωm is that the conductivity has a strictly linear frequency dependence in the 100 kHz to 1 MHz range. For SiO2, previous measurements of the conductivity of thermally grown oxide on silicon44, demonstrate \({\sigma }_{{\mathrm{SiO}}_{2}}\approx 1{0}^{-14}\,{\Omega }^{-1}\,{{\rm{m}}}^{-1}\) at d.c. and a weak frequency dependence. For stoichiometric silicon nitride a similar value has been reported45. With these values for the conductivity, the contribution to the loss tangent is on the order of 10−9, much smaller than our estimated values.

Fundamental thermal-electrodynamic damping

In the absence of static charges, the finite thermal fluctuations of the electromagnetic field gives rise to the ultimate form of NCF46. This mechanism that often goes by different names (for example, vacuum friction, van der Waals friction and Casimir friction) is a direct consequence of the fluctuation–dissipation theorem for the electromagnetic field in the presence of lossy dielectrics. When two dielectrics are in relative motion, the retardation of the fluctuating electromagnetic field in conjunction with the force exerted on the dielectrics through Maxwell’s stress tensor, gives rise to friction forces. Estimating the damping force for arbitrary geometries is complicated, but there is an analytical solution for the case of two semi-infinite bodies47. We estimate the damping coefficients, \({\gamma }_{{\rm{NCF}}}^{{\rm{th}}}\), for silicon nitride bodies for both parallel and perpendicular relative motions. As shown in Extended Data Fig. 3e, even for distances of tens of nanometres, \({\gamma }_{{\rm{NCF}}}^{{\rm{th}}}\) is orders of magnitude smaller than our observed values. Moreover, it has been shown that for a mechanical oscillator, this form of damping does not depend on the harmonic potential of the oscillator (that is, \({\gamma }_{{\rm{NCF}}}^{{\rm{th}}}\) is frequency independent)48.

Dissipation–dilution limit to force sensitivity

For a string with thickness h, stress σ and mass density ρm, the ultimate limit for the Q at a given frequency Ω from dissipation–dilution (also known as the clampless limit) is given by10

$$\frac{{Q}_{m}^{\mathrm{cl}}}{{Q}_{\mathrm{int}}}\le \frac{12{\sigma }^{2}}{{\rho }_{m}E{h}^{2}{\varOmega }^{2}},$$

(18)

where E is the Young’s modulus. To obtain a limit for the force sensitivity at a given string width, w, we set the effective mass to the value for the fundamental mode of a uniform string with frequency Ω, \({m}_{\mathrm{eff}}^{\mathrm{cl}}=\frac{{\rho }_{m}wh}{2\varOmega }\sqrt{\frac{\sigma }{{\rho }_{m}}}\). This is the minimal value for a given thickness and width of the string, and for soft-clamped modes, this value is up to a factor 2 larger. Together, the limit of force sensitivity is given by

$${S}_{{\rm{F}}}^{\mathrm{th},\mathrm{cl}}=4{k}_{{\rm{B}}}T\frac{E{h}^{3}w}{24{Q}_{\mathrm{int}}}{\left(\frac{{\rho }_{m}}{\sigma }\right)}^{\frac{3}{2}}{\varOmega }^{2}.$$

(19)