Optically levitated nanoparticle as a model system for stochastic bistable dynamics

Nano-mechanical resonators have gained an increasing importance in nanotechnology owing to their contributions to both fundamental and applied science. Yet, their small dimensions and mass raises some challenges as their dynamics gets dominated by nonlinearities that degrade their performance, for instance in sensing applications. Here, we report on the precise control of the nonlinear and stochastic bistable dynamics of a levitated nanoparticle in high vacuum. We demonstrate how it can lead to efficient signal amplification schemes, including stochastic resonance. This work contributes to showing the use of levitated nanoparticles as a model system for stochastic bistable dynamics, with applications to a wide variety of fields.

The full experimental set-up is schematically shown in Fig. 1. The light source is a an ultra-stable low noise Nd:YAG laser, and is used to both trap the silica nanoparticle (d ∼ 177nm in diameter) and track its motion. After passing a Faraday isolator (FI), the beam is split into a high power branch (trapping beam) and a low power branch (detection beam) by means of a polarizing beam splitter (PBS). To decouple the two branches and to avoid interferences in the laser focus we use cross polarized beams and we shift the frequency (∆f ≈ 110MHz) of the detection beam with an acousto-optic modulator (AOM). The trapping beam passes through the electro-optic modulator (EOM) that allows the parametric actuation of the particle. After convenient beam shaping with telescope systems (TS) the two paths are recombined with a second PBS before entering in the vacuum chamber where a high NA microscope objective (OBJ) focuses the beams -creating the optical trap -and an aspheric lens (AL) collects and collimates the light, further sent out of the chamber carrying the information about the particle's motion. After the vacuum chamber the trapping beam is dumped, while the detection beam is analysed with a split detection scheme. This scheme relies on the fact that the interference pattern between the transmitted light and the light scattered from the particle is position dependent. Therefore, by half-splitting the beam with a sharp edge D-shaped mirror (DM) and sending the two halves to a balanced photodetector, it is possible to measure the oscillating component of the input signal, which is proportional to the particle's position. Figure 4(a) shows the ability of such a scheme to track particle position with a noise floor on the order of 1 pm/ √ Hz [1]. Below ∼ 6 mbar pressures, the particle enters in the ballistic regime, and a spectral analysis of its motion shows the presence of three peaks at the eigenfrequencies of the oscillation modes. Further decreasing the pressure we reach an highly underdamped regime (P 10 −6 mbar), where the particle oscillates with very small amplitude and corresponding sub-kelvin effective temperatures. At these temperatures the three mechanical modes of the oscillation are fully decoupled and the restoring gradient force depends linearly on particle's displacement by means of the linear trap stiffness κ i = αE 2 0 /w 2 i , where α 0 is the polarizability of SiO 2 , E 0 is the electric field intensity at the focus and w i is the beam waist radius (i = x, y) or Rayleigh range (i = z). Due to the asymmetry of the trapping potential, the three mechanical modes along the three main axes have different oscillation frequencies: (Ω x , Ω y , Ω z )/2π (125, 140, 40) kHz, as we show in Fig. 2. Under this non-degentrate conditions, the feedback scheme can therefore act without cross-coupling between the three signals.
In the absence of feedback cooling, or for coherent driving of one of the modes, the amplitude of the oscillation becomes large and the particle explores the anharmonic part o the optical potential, which features Duffing nonlinearities. In this regime, the restoring force is well approximated by F grad i = −κ i 1 + j=x,y,z ξ j x 2 j x i , where x i is the displacement of the particle from the trap centre and the nonlinear coefficients are given by ξ j = −2/w 2 j [2]. The coefficients ξ j can be calculated with high accuracy measuring the particle motion in the nonlinear oscillation regime, as explained in the main text and shown in Fig. 1d.

SUPPLEMENTARY NOTE 2: Feedback electronics optimization
The high stability of our system takes partially merit from a fine optimization of the feedback electronics, with particular regards to the phase shifter and frequency doubler module, responsible for the implementation of an efficient cooling. In fact the relative phase between particle oscillation and feedback signal is a crucial parameter: a phase shift of π/2, for example, turns cooling (nonlinear Supplementary Figure 2: Power spectral densities of the decoupled modes Experimental power spectral densities S i (i = x, y, z) of the particle's 3D motion at P ∼ 6 mbar (dots), and corresponding Lorentzian fitted function (solid lines). The peaks appear at non-degenerate frequencies, which allows the feedback to act separately on each oscillation mode. Also, decoupling is ensured in our system, as almost no trace of orthogonal oscillation is present in each considered mode.
damping term η > 0) into amplification (η < 0). Therefore, a fine tuning of the feedback phase is of uttermost importance. As explained in the main text, one way to retrieve the value of the nonlinear damping coefficient η is to measure and fit the nonlinear amplitude response of the particle during parametric driving frequency down-sweeps, which reads: Fig. 3 sketches the procedure that allows to optimize the feedback efficiency for the three channels (oscillation modes) separately. For each value of phase φ we measure the particle response inside the Arnold tongue ( 3(a)). Then, by fitting the oscillation amplitude at different values ( 3(b)), and discarding the fits with an R-square lower than a certain threshold (=0.995 for the case in Fig. 3(c)) we retrieve η and the relative standard deviation σ η from a weighted averaging over the fits outcomes. The procedure is repeated for each mode x, y, z. We clearly from Fig. 3(d) that the non-linear damping coefficients (η x , η y , η z ) display a maximum in their trend, which consequently provides the optimal phase delays (φ x , φ y , φ z ) to be set.

SUPPLEMENTARY NOTE 3: Frequency stability
Frequency stability of the resonator is monitored by tracing the resonance peak position in time.
Examples of relative frequency time traces are shown in Fig. 4(b). At high pressure, and with no feedback (purple trace), relative fluctuations are of the order of σ Ω /Ω 0 ≈ 10 −2 . However when high vacuum is reached and the feedback cooling is turned on (orange trace), we find σ Ω /Ω 0 ≈ 10 −4 , and a net fluctuation reduction of two orders of magnitude is observed. The remarkably small values obtained correspond, in absolute terms, to fluctuations of ∼ 30 Hz, smaller than the effective linewidth of the resonator (∼ 100 Hz, considering feedback broadening effects). The considerable improvement in frequency stability, as compared to previous experiment [2,3], has been achieved by reducing the two main sources of disturbances, namely the optical noise introduced by the feedback cooling, and air turbulences affecting pointing stability of both the trapping beam and the detection scheme. A comprehensive description of feedback optimization is given in Supplementary SUPPLEMENTARY NOTE 2:, while air disturbances problems have been circumvented by hardware screening the optical set-up.

SUPPLEMENTARY NOTE 4: Nonlinear response
Parametrically driving the particle into the nonlinear regime leads to an expected amplitude response given by [3,4] Performing frequency down-sweeps and up-sweeps also provides interesting spectral features of the two oscillation modes. These are shown in Fig. 5a. The low amplitude mode A low consists of a frequency pulled main mode at Ω P , together with the drive sideband at Ω S . The position of these peaks depends of the driving parameters (δ m , ) according to the relation Ω P, , with the ± sign corresponding to P and S subscripts respectively. Instead, the high amplitude mode A high is a locked steady state oscillation, with fixed frequency Ω 0 (1 + δ m /2) = Ω m /2. These frequency relations along up-and down-sweeps are depicted in Fig. 5b, where we show experimental data (dots), together with analytical prediction (solid lines). These comparisons, namely both amplitude fits and iso-amplitude lines of Fig. 1d (main text), and the frequency curves in Fig. 5c demonstrate the high level of accuracy of our measurements and the perfect agreement with the model given by Eq. (1) of the main text. The analytical prediction of the iso-amplitude lines shown in Fig. 1d (main text) directly follows from eq. 2 after representing the amplitude response as a 2-dimensional surface in the (δ m , , A) space (cross-cuts in Fig. 1d actually give a hint of such a surface). By setting A 2 (δ m , ) = A 2 ∆ = const, one selects an horizontal cut of the surface, having the property of maintaining a constant amplitude , that projected in the (δ m , )−plane identifies the iso-amplitude line. Mathematical derivation is as follows:

SUPPLEMENTARY NOTE 5: Noise Calibration
The system behaviour after noise injection is shortly described in the main text and here supported by the graphic data in Fig. 6. The 1 MHz bandwidth Gaussian noise process ζ(t), used to drive the particle in the stochastic switching regime, satisfies ζ(t) = 0 and ζ(t)ζ(t ) = δ(t − t ). When noise amplitude N V is varied from 0 to 8 V, we obtain the discussed increase in oscillation amplitude that provides the noise temperature calibration T N ≈ N 1.44 V (see Fig. 6a). The histogram in the inset shows the energy distribution of the particle derived from a sample of 10 5 measurements. The

ρ(E)
Supplementary Figure 6: Noise temperature calibration without parametric driving. a, When injecting optical noise in the system the particle's oscillation amplitude increases (blue dots), together with its standard deviation (shaded area). Converting amplitudes into effective temperatures it is possible to map noise levels into noise temperatures and extract a power law calibration curve. The inset shows the measured energy distribution ρ(E) at the thermomechanical noise level.
dashed yellow line is a fitted analytical curve of the form: where, given the high accuracy in the calculated value of η (see main text Fig. 1d), γ 0 remains the only fitting parameter.
Tracking the time evolution of the amplitude provides information about its spectral properties. By applying the slow varying envelope approximation to the equation of motion of the system in the absence of coherent driving (Eq. (1) main text, having set = 0), one finds that the amplitude behaves as an overdamped variable, in spite of being the result of the average over a fast (underdamped) oscillating variable [3]. This behaviour is depicted in Fig. 6b where the amplitude power spectral density S A is shown for different noise temperatures. We fit the measured power spectral densities (PSDs) with a typical overdamped Lorentzian S A = S 0 is the height of the plateau and f c is the corner frequency. The inset shows the power law dependency of these fitting parameters (in terms of relative change with respect to the no-injected-noise case) as a function of noise temperature. As expected, the low-frequency plateau raises with the noise intensity, in agrement with the observation in Fig. 6a [5,6]. The corner frequency also increases with the injected noise, leading to a larger bandwidth of amplitude fluctuations, this likely being due to the fact that the optical noise becomes multiplicative after the EOM [7] SUPPLEMENTARY NOTE 6: Symmetric effective potential When measuring Kramers' transition rates a static symmetric configuration of the potential is needed. The same applies when performing the Stochastic Resonance (SR) experiment, where on top of a symmetric potential the harmonic perturbation contributes with a modulation of the shape and consequently of the switching probability. If such a starting symmetric condition is not ensured, unwanted effects of probability quenching are encountered, and proper SR phenomena are not observed. Therefore, a preliminary fine characterization to know where this symmetry configuration is found inside the bistable regime is of utmost importance. To do so we sweep the detuning δ m horizontally in a small region around this symmetry condition. During each (discretized) sweep we measure stochastic switching time traces and, identically to what explained in the main text, we retrieve the effective potential from from the amplitude histograms. Such a scan provides the data shown in Fig. 7. Fig. 7a,b show the amplitude of stable (A low and A high ) and unstable A un states, together with the amplitude gap A ∆ , as a function of the detuning δ m . As the scan inspects only a narrow region of the bistable tongue the amplitude of the stable states, and consequently A ∆ are kept constant, while the unstable state shift towards lower amplitudes due to the dynamical changes in the potential shape (see main text   Fig. 2d of the main text, state populations Σ low and Σ high invert at the symmetric configuration δ * m . As a result, fitting with sigmoid functions σ low (δm) and σ high (δm) (see eq. 5)provides information on the position δ * m of the symmetric configuration and on the width τ of the bistable regime.
also evident in Fig. 7c, where the hight of the potential barriers separating respectively the lower from the higher state and the vice versa, invert their values passing through the quasi-symmetric configuration. Correspondingly, the state populations in Fig. 7d also undergo a similar dynamics. The fitting function in the latter two plots are typical sigmoid function of the form where the parameter δ * m and τ respectively represent the horizontal position of the symmetric configuration inside the bistable tongue, and the effective width of the bistable regime. Repeating these measurements for higher levels of noise in the system yields to a clear dependency of these parameters from the noise temperature T N . This trends are displayed in Fig. 8 (circles) and width τ of the bistable regime (shaded region) as a function of noise temperature T N . The linear fit (solid line) provides an intrinsic characterization of the system that allows to follow the symmetric potential configuration when higher levels of noise are injected. Red crosses correspond to data points excluded from the fit.
these measurements allows to conclude that when noise temperature is changed in the system, an adjustment to the driving parameters, and more precisely to detuning δ m is necessary to lock on to the symmetric configuration. We believe the reason for this resides in two shift-contribution. On the one hand there is the noise-induced particle frequency drift (also observed in [8,9]), which automatically requires the detuning to follow. On the other hand, instead, the multiplicative nature of the noise [7] could also play a role when entering in the fully noise-driven dynamics of the particle.

SUPPLEMENTARY NOTE 7: Linear perturbation theory of stochastic resonance
According to linear perturbation theory, the expected amplitude power spectral density (PSD) has the form: This equation represents the analytical prediction of the powr spectral density(PSD) shown in the inset of Fig. 3a (main text). Here the first term corresponds to the typical PSD of an overdamped variable, while the second one containing the δ-function is related to the harmonic perturbation. It is