Coherent oscillations of a levitated birefringent microsphere in vacuum driven by nonconservative rotation-translation coupling

Experiment

Figure 1A shows the trapping geometry and coordinate axes used throughout this article. Vaterite microspheres have a spherulitic structure in which the optic axis of this highly birefringent material traces out nested hyperbolae arranged around a unique symmetry axis (section S1). A single vaterite microsphere with a radius of 2.2 μm is held in a Gaussian trap, linearly polarized in the x direction and propagating in the positive z direction with a wavelength of 1070 nm and beam power of 10 mW (measured at the back aperture of the microscope objective). Its motion is tracked with a fast complementary metal-oxide semiconductor (CMOS) camera (see Materials and Methods).

Figure 1B shows time-lapse images of the particle trapped at a residual gas pressure of 1.2 mbar [see movie S1, which is rendered at 15 frames per second (fps) from 5000 fps]. First, the particle aligns its symmetry axis (

in Fig. 1A) with the electric polarization (i.e., the x direction). The center of mass (CoM; yellow crosses) preferentially oscillates along the polarization direction with an amplitude that can be as large as ±1 μm compared to that of ±0.1 μm in the orthogonal y direction (Fig. 1C).

Figure 2 describes the strong dependence of the oscillatory x motion on the gas pressure. Column (A) shows the CoM positions of the vaterite particle in the x–y plane (red dots) measured at every 0.2 ms for the duration of 1 s (5000 data points in total) at different gas pressures. The oscillation amplitude in the x direction increases as the residual gas pressure is reduced by a few millibars, from 4 to 1.2 mbar. The motion is compared with that of an isotropic silica microsphere with a radius of 2.5 μm under the same trapping conditions (Fig. 2A, blue dots), confirming that the oscillatory motion is unique to the birefringent particle. Distributions of the x coordinate are given for vaterite and silica particles in Fig. 2B. In accordance with the equipartition of energy (Eq. 1), the distributions for silica are Gaussian and approximately independent of pressure, while those for vaterite broaden with decreasing pressure, finally changing form altogether. The power spectrum of the x coordinate clearly indicates the trap or oscillation frequency at around f_x (=ω_x/2π) ∼ 530 Hz, which remains the same at lower pressures (see red curves in Fig. 2C). Notably, the trapped vaterite particle exhibits much larger Q factors than those of the silica particle (blue curves), and kinetic energy is increasingly concentrated at a single resonant frequency.

Autocorrelation 〈x(t)x(t + τ)〉 (Fig. 2D) measures the rate at which the translational x(t) motion loses coherence because of thermal fluctuations. Vaterite exhibits an order of magnitude longer decay time than that of the silica (Fig. 2D, 4). These features correspond to greatly increased Q values for vaterite that, as discussed below, can be further enhanced through parametric driving. We note that as the gas pressure is decreased, inertial forces tend to exceed gradient forces, and the particle may leave the trap typically at a residual gas pressure of <1 mbar, depending on the optical power and the mass of the particle.

Last, Fig. 3A shows the position variance of the trapped particle in the x direction (red crosses) compared with that in the y direction (blue plus signs), where 〈y²〉 is directly associated with the thermal energy through

, where K_yy is the trap stiffness in the y direction and 〈y²〉 = 4.49 × 10⁻³ μm² at the pressure P = 4 mbar, assuming that the particle is at room temperature. The ratio 〈x²〉/〈y²〉 therefore approximates the relative energy in the x motion (Fig. 3A on the right axis). Figure 3B shows the reciprocal of the position variance against P, where experimental data are fitted with linear regression. Here, we clearly show the linear relationship of 1/〈x²〉 with P, which is predicted by the theory (see Eq. 13).

Simulation

To better understand the experimental observations, we performed direct numerical simulations of the thermal motion of a homogeneous birefringent sphere with default parameters (radius of 2.2 μm and density ρ = 2.54 g/cm³ with ordinary and extraordinary refractive indices corresponding to those of bulk vaterite n_o = 1.55, n_e = 1.65) in a linearly polarized Gaussian beam in vacuum (see section S3 for further details). In summary, we numerically integrate the equation of motion

(2)

f^optis the generalized optical force (i.e., forces and torques) acting on the particle, and q are the generalized coordinates specifying the CoM and orientation [i.e., q = (x, y, z, θ, ϕ), where θ and ϕ are the polar and azimuthal angles specifying the symmetry axis of the particle,

]. Ξ is the pressure-dependent friction matrix with diagonal elements Ξ_ii = 6πμa and 8πμa³ for translational and rotational motion, respectively. f^L is the stochastic, Langevin force (and torque), uncorrelated, with zero mean and amplitude fixed by the fluctuation-dissipation theorem, i.e.,

,

. M is a diagonal matrix whose elements are given by the mass (m) of the sphere and its moment of inertia (I), and

is the weight. The beam propagates in the positive z direction with its axis coincident with the z axis. Although the model is idealized (ignoring the inhomogeneity and roughness of the particle, aberrations in the optical beam, etc.), the behavior represented in Fig. 2 is convincingly reproduced (see fig. S3). In particular, fig. S3 shows simulated distributions of the CoM, power spectra, and autocorrelation functions, which should be compared with Fig. 2 (A, C, and D, respectively).

In addition, the simulations reveal strong coupling between the rotational and translational degrees of freedom of the sphere. As the viscosity is reduced, the growing oscillation of the x coordinate of the CoM is increasingly accompanied by an oscillation of the angular coordinate, θ. Figure 4 shows this correlation as a scatterplot (left-hand side) and as the time-dependent cross-correlation of a simulated Brownian trail. The underlying mechanism is clear. Angular oscillations produce oscillating force components,

, through rotation-translation coupling, which drive linear oscillations in x. Similarly, oscillations in x generate torque components,

, which feed into the angular oscillations. As the pressure (viscosity) is reduced, the coupled oscillations increasingly reinforce one another and the motion becomes more coherent and self-sustaining. This effect is very clearly seen in movies S2 and S3, which show the simulated stochastic motion in optical traps whose parameters correspond to the three cases represented in Fig. 4 (A and B).

Analysis

These observations can be quantified with reference to a simple analytical model. Linearizing the force field relative to the trapping configuration (in which external forces and torques vanish) and restricting attention to the coupled coordinates, θ and x, give

(3)

The first expression is the linearized equation of motion. The external (i.e., optical and gravitational) forces on the particle are given by

, where

is a stiffness matrix that includes coupling coefficients e.g.

(4)and q = (x, θ) are the coupled coordinates. Ξ is the friction matrix for a sphere of radius a in a fluid with viscosity μ, with diagonal components ξ_xx = 6πμa for translation and ξ_θθ = 8πμa³ for rotation. The second expression in Eq. 3 is the frequency domain equivalent, with Q = (X, Θ), the Fourier transform of q. The eigenfrequencies of A, denoted by ω_i, determine qualitative features of the trap including its linear stability, power spectral density (PSD), and correlation functions. In particular, if K is symmetric, then the force is locally conservative with (∇ × f) = 0 and the trap is at thermal equilibrium. However, this need not be the case for optical fields or other momentum flows (8). Complete calculations of these quantities are provided in the Supplementary Materials, and the main results are quoted and discussed below.

The real parts of the eigenfrequencies,

, describe oscillatory behavior. The imaginary parts,

, relate to relaxation [when

] or linear instability

]. For isotropic spheres, coupling is absent and K_xθ = K_θx = 0. The eigenfrequencies are then

(5)so that the sphere oscillates with the natural frequency of a harmonic trap, while its relaxation is determined by the Stokes’ drag and the mass; the trap is linearly stable with

. These traps are necessarily conservative within the linear regime and remain at equilibrium. When motional degrees of freedom are coupled, this constraint is lost. It is helpful to consider two extremes: the overdamped or low–Reynolds number regime (negligible inertia) and the underdamped limit (negligible viscosity).

In the low–Reynolds number regime, the eigenfrequencies are purely imaginary and positive unless K_xθK_θx > K_xxK_θθ (see section S4). This is a demanding condition, equivalent to the requirement that the coupling forces and torques exceed the restoring terms, and unlikely to be satisfied.

In the low-pressure limit, where viscosity is negligible, we have

(6)with

(7)

Setting the K_xθ = K_θx = 0 results in the familiar natural frequencies for harmonic oscillators, e.g.,

. Characteristic frequencies with negative imaginary parts emerge whenever Δ_LP < 0, destabilizing the trap. This requires that K_xθ and K_θx have opposite signs in Eq. 7. Physically, this condition ensures that the rotational and translational vibrations feed into one another constructively so that, for example, the rotational vibration generates an oscillating force that is in phase with the translational vibration and vice versa. Providing that K_xθK_θx < 0, the condition Δ_LP < 0 is relatively easily satisfied, especially when the natural rotational and translational frequencies are close as they are in this case. Calculated forces and torques are shown in Fig. 4 (C and D) for the default parameters given above, showing that K_xθK_θx < 0 as required. We note that calculated values of Δ_LP are negative for the default system parameters and for surrounding regions of parameter space (see section S4).

Whenever the trap is linearly stable at low Reynolds number (as usual), but unstable in the low-pressure limit (Δ_LP < 0 in Eq. 7), one of the values of

must pass through zero as the ambient effective viscosity is reduced. Requiring

yields the following expression for the critical effective viscosity, μ_X, at which one of the

(see section S4)

(8)where the friction coefficients ξ_xx = 6πμa ≡ s_xμ and ξ_θθ = 8πμa³ ≡ s_θμ are proportional to the viscosity and ω_X is the corresponding real part of the frequency

(9)

In the absence of rotational-translational coupling (K_xθ = K_θx = 0), Eq. 8 gives μ_X = 0, and when the condition

cannot be satisfied (because K_xθK_θx > 0, for example), a negative value is obtained for

, indicating a nonphysical solution. A survey describing typical variations in μ_X and ω_X with varying properties of the sphere and beam is included in section S4. The properties of the experimental system fall within realistic ranges given uncertainties in the beam power, for example, and the idealized nature of the model.

As μ → μ_X, the oscillations become more coherent and grow in amplitude. This is expressed in the PSD, which describes the distribution of motional power in frequency space and is given by the ensemble average of the squares of the Fourier-transformed coordinates, 〈Q(ω)Q*(ω)〉. As

approaches zero, the PSD for the vaterite sphere, with rotation-translation coupling, consists of a single peak with the maximum value, 〈Q(ω)Q*(ω)〉_max, and width, Δω (section S4)

(10)

(11)

Thus, the height of the peak in the PSD increases rapidly as μ_X is approached and its width, Δω, decreases toward zero. Time correlations take the following form (section S4)

(12)

As

, the amplitude

increases so that, for example, the instantaneous position variance, 〈x²〉, grows. Since the relationship between pressure and viscosity is approximately linear in this regime (19), we can write

, where P_X is the pressure corresponding to μ_X. Combining with Eqs. 10 and 12 gives the scaling behavior for the PSD maxima and linewidths as well as covariance and relaxation times. For example, 〈x²〉 is

(13)where P_X is the pressure corresponding to the critical viscosity, μ_X. In Fig. 3B, we confirm that this relation holds for the experimental data and use it to extract a value for P_X = 1.2 mbar. Previous measurements (see section S2) show that the corresponding viscosity is μ_X = 0.2 μPa s, which falls within the range predicted in the parametric survey (section S4). Similarly, the relaxation time,

becomes large as μ → μ_X, and the oscillations of the birefringent sphere become increasingly coherent.

The behavior described above, in Eqs. 10 to 12, is qualitatively different from that of isotropic spheres, for which the position variance, 〈x²〉, is independent of viscosity (Eq. 1), while the value of

remains finite for finite viscosity, ensuring that the linewidth and relaxation times also remain finite.

In summary, whenever Δ_LP < 0, the stochastic motion within the linear approximation is characterized by a critical viscosity, μ_X. For μ > μ_X, the motion is driven by thermal fluctuations, without which the particle would remain motionless in the trapping configuration. As the viscosity is reduced toward μ_X, spontaneous oscillations emerge with a linewidth, Δω, approaching zero and quality factor

tending to infinity. In reality, this limit is never reached since the position variance, 〈x²〉, increases rapidly, taking the particle beyond the linear range of the trap. For viscosities μ < μ_X, the linear approximation breaks down completely and the oscillations become self-sustaining. In essence, the particle takes energy from the optical field and dissipates it into the surrounding gas, forming a nonequilibrium steady state. The oscillations may be regarded as being driven, but the external force field is time invariant. Thermal fluctuations introduce phase errors that the system cannot correct since it has no knowledge of what the phase should be at any particular time. These errors accumulate and limit how narrow the linewidth can become.

Ultrahigh Q oscillators

Although the extremely high Q values predicted by the linear approximation (Eq. 11) cannot be obtained in practice, parametric modulation of the trap intensity is observed to drastically improve the quality factors of the free-running oscillator. This remarkable effect is due, in part, to the oscillations being confined to a single spatial dimension, analogous to the case of the tethered pendulum (20). In addition, the parametric drive introduces a periodic time variation against which phase errors can be measured, allowing the oscillator to synchronize with the drive (21, 22).

Figure 5A shows the PSD of the oscillating particle at around the resonant frequency f₀ = 832 Hz (centered at 0 Hz). When the particle is trapped with a stationary light field both spatially and temporally, the oscillation yields the resonant peak (orange dots), which is fitted with a Lorentzian (blue curve) with a linewidth of 0.90 ± 0.13 Hz, yielding Q = 924 ± 136 (Fig. 5A). Because of collisions with residual gas molecules and CoM excursion of the particle into regions of different light intensity, the oscillation frequencies exhibit a broad distribution.

However, when the particle is trapped by a light field with an amplitude periodically modulated at a frequency f = 2f₀ (i.e., parametrically driven), the mechanical mode becomes more coherent, leading to an ultranarrow linewidth. In this regime, there are two driving protocols: (i) The external modulation can be locked to the phase of the particle motion (either inward or outward, relative to the trap center), i.e., “phase locking.” (ii) The particle motion can be locked to the external clock, transducing the frequency of the time standard to the particle motion, i.e., “frequency locking” (see Materials and Methods). Figure 5B shows the PSD when phase locked, exhibiting a linewidth of 10.5 ± 3.0 μHz at f₀ ≈ 580.5 Hz and the corresponding Q = 5.53 ( ± 1.72) × 10⁷. These mechanical properties are further improved when frequency locked. In this case, the particle acts as a precision micromechanical transducer, allowing long-time measurements required for high-Q experiments. Figure 5C shows the PSD with a linewidth of 2.20 ± 0.62 μHz at f₀ ≈ 539.6 Hz and Q = 2.45 ( ± 0.75) × 10⁸ when frequency locked. These are, to the best of our knowledge, the narrowest linewidth and the highest Q factor reported to date for a mechanical oscillator at room temperature.

To demonstrate weak-force sensing capabilities of the driven oscillators, the phase-locked particle is used as a probe for real-time frequency measurements in a dilute gas (see Materials and Methods). A fractional change in gas pressure affects the damping rate of the particle, which translates gas pressure values onto the oscillation frequency. Figure 6 shows the oscillation frequencies (orange dots) depending on the residual gas pressure in the pressure range of 4 to 7 mbar. A linear fit to the data (blue solid line) yields a rate of frequency change of 0.21 Hz/mbar with a relative pressure sensitivity of 0.75% (2σ), which can be improved by increasing the length of the measurement (currently 100 s). We note that frequency locking can also be used, where the phase lag between the drive and the oscillator depends on nonconservative forces, such as light or gas scattering (15). This demonstrates the great potential of this self-sustained oscillator as a weak-force sensor.