Nonlinear squeezing of stochastic motion

Stochastic mechanical systems, characterized by fixed mass m, are typically assumed to be immersed in a Markovian environment at constant temperature T > 0 with fixed linear damping Γ. Since, in the classical description, any noise can be eliminated trivially by cooling to zero temperature, useful theoretical thresholds will always require prior assumptions such as these.

The Markovian stochastic dynamics for a one-dimensional mechanical system in an arbitrary potential are described by the Langevin equation

$$\begin{equation} \ddot{x} = -\frac1m\partial_xV\left(x\right)-\gamma\dot{x}+\sqrt{\frac{2\gamma}{m\beta}}\xi\left(t\right)\,, \end{equation} \tag{ 1 }$$

where V(x) is a potential depending only on the mechanical position, $\beta = 1/k_{\mathrm {B}} T$ is the inverse temperature, with T the fixed temperature of the environment and $k_{\mathrm {B}}$ the Boltzmann constant, $\gamma = \frac{\Gamma}{m}$ is the medium damping with Γ the drag coefficient of the medium (e.g. air at low pressure) and $\xi(t)$ is the Langevin force of a thermal bath satisfying the Markov conditions $\langle{\xi(t)} = 0$ and $\rangle{\xi(t)\xi(t^{\prime})} = \delta(t-t^{\prime})$. Observe that both the damping and the Langevin force depend on γ, so that γ = 0 retrieves conservative motion.

For any potential locally described by an up-to-quadratic potential $V_2(x) = Fx+\frac12m\omega^2 x^2$, with ω the frequency and F arising from a constant force, the system is defined to be linear. Many nonlinear oscillators operate in a regime close to the equilibrium position which can be well approximated by such a linear potential. Such oscillators are referred to as linearized and are included in the following derivation of the threshold. In deriving the threshold we will assume that ω and F are constant in time during a single experimental run and may only vary randomly between independent experimental runs, as opposed to parametrically driven systems, in which frequency may be a periodic function of time. For many experiments, the position statistics are accumulated over time and the momentum statistics are well-estimated from numerical differentiation of the measured position [16] and so in reality, even when m, β and γ remain constant during measurement the potential might vary. This introduces a statistical mixture over the potential parameters ω, F of equation (1). However the thresholds for nonlinear squeezing will turn out to be independent of such variation.

The simplest nonlinearity taking the system beyond the linear dynamics of $V_2(x)$ is the pure cubic potential $V_3(x) = k_3x^3$. This potential cannot be approximated by a harmonic potential in any limiting cases. Motivated by the short-term deterministic dynamics of a massive particle in a cubic potential [23, 24, 26] we examine the nonlinear variable $p+\lambda x^2$. This variable is chosen as an exemplar, however the principles applied to derive the threshold may be applied to other variables such as $p+\lambda x^n$, providing a collection of thresholds for nonlinear motion. Here $\lambda\in\mathbb{R}$ is a parameter characterizing the choice of nonlinear variable, analogous to varying the angle in phase space for linear variables [7].

Nonlinear variables such as $p+\lambda x^2$ involve higher order moments. For Gaussian distributions, the higher order moments are fixed by the first and second, and thus cannot vary arbitrarily. We will leverage this property into a noise threshold that captures when the measured statistics of a dynamical system cannot be explained by any model of the form of equation (1) under up-to-quadratic potentials.

In order to describe how to construct this threshold, we begin by observing some general properties of equation (1) under linear potentials. For an initially Gaussian phase space distribution the Gaussian form is preserved [27]. Thus, we may assume that at any time t, even in the limit $t\rightarrow\infty$, the distribution ρ and its associated covariance matrix Σ take the following forms

$$\begin{align} \rho& = \frac{1}{2\pi}\frac{1}{\sqrt{\text{Det}\left(\Sigma\right)}}e^{-\frac12\left(X-\mu\right)^\top\Sigma^{-1}\left(X-\mu\right)} \end{align} \tag{ 2 }$$

$$\begin{align} \Sigma& = \begin{pmatrix}\sigma_x && \sigma_{xp}\\ \sigma_{xp} && \sigma_p\end{pmatrix}\,, \end{align} \tag{ 3 }$$

where $X = \begin{pmatrix}x & p\end{pmatrix}^\top$ and $\mu = \begin{pmatrix}\bar{x} & \bar{p}\end{pmatrix}^\top$ is the corresponding vector of mean values.

In this general form the variance $\text{Var}(p+\lambda x^2)$ can be evaluated via a Gaussian integral to be

$$\begin{equation} \text{Var}\left(p+\lambda x^2\right) = \sigma_p+2\lambda\left(2\bar{x}\sigma_{xp}+\lambda\sigma_x\left(2\bar{x}^2+\sigma_x\right)\right)\,. \end{equation} \tag{ 4 }$$

Observe the independence of this expression from $\bar{p}$, which is important in applying the threshold to mixtures of Gaussian distributions. Also note that, since $\text{Det}(\Sigma) = \sigma_x\sigma_p-\sigma_{xp}^2$, we may replace the correlations σ_xp with the determinant. Then it is straightforward to minimize the variance for Gaussian distributions. The result is the threshold

$$\begin{equation} \mathcal{T}_0 = \frac{3\left(\text{Det}\left(\Sigma\right)\lambda\right)^\frac23}{2^\frac13}\,. \end{equation} \tag{ 5 }$$

It is worthwhile to reflect on the meaning of this threshold. Formally, if the statistical quantity $\text{Var}(p+\lambda x^2)$ of a pair of random variables x and p is below the threshold then a certain class of multivariate normal distributions (those with bounded $\text{Det}(\Sigma)$) are excluded as possible models for the multivariate probabilistic system.

To connect this to the noisy damped linear dynamics of equation (1) we must associate a physical meaning with Σ. For conservative linear dynamics, γ = 0, the determinant of the covariance matrix is preserved and, via the decomposition of Williamson's theorem [28], always associated with a fixed thermal distribution of a fiducial harmonic oscillator, for which we have $\text{Det}(\Sigma) = \frac{1}{(\beta\omega)^2}$. Nontrivial noise thresholds can only be derived when bounds are placed on temperature, and so we have assumed constant β. What remains is to set a fiducial, or reference, frequency ω₀ associated with the fiducial harmonic oscillator. We then express the conservative threshold as

$$\begin{equation} \mathcal{T}_0 = \frac{3|\lambda|^\frac23}{2^\frac13\left(\beta\omega_0\right)^\frac43}\,. \end{equation} \tag{ 6 }$$

Any nonlinear squeezing below the threshold, resulting from conservative dynamics and under the conditions specified, i.e. fixed mass and temperature and an associated reference frequency ω₀ indicates conservative non-Gaussian motion.

A quantum mechanical treatment of this problem, replacing the fiducial states at frequency ω₀ with the minimum uncertainty states, recovers the threshold for quantum non-Gaussianity reported in [29]. The threshold of equation (6) for the transient dynamics does not depend on mass m and holds even for stochastic mixtures of the distributions arising from conservative linear dynamics with random linear potentials (see SM for details).

More generally, for γ > 0, $\text{Det}(\Sigma)$ is not preserved. This can be easily seen by the fact that for $t\rightarrow\infty$ a new equilibrium distribution is achieved which may have arbitrary $\text{Det}(\Sigma)$ when the dynamical frequency ω is different from the fiducial frequency ω₀. At first glance it appears that loosening the restriction on $\text{Det}(\Sigma)$ trivializes the threshold, but remarkably this is not the case. Consider that in the fully general case of noisy damped dynamics (equation (1)) the system relaxes from the state associated with the fiducial frequency ω₀ to a new equilibrium state defined by the frequency ω and a linear force F. During the relaxation process the threshold equation (6) is indeed violated, indicating that the conservative motion threshold can be passed merely by the application of damped linear dynamics. However thorough numerical analysis reveals that such damped dynamics are also subject to limitations due to their linearity. Since the distribution remains Gaussian it is not possible to achieve arbitrarily low values of $\text{Var}(p+\lambda x^2)$ simply by changing the potential parameters, even if no simple geometric picture is possible.

However the threshold must now be modified to accommodate any transient damping effects up to and including equilibration to the thermal bath. For the up-to-quadratic potential $V_2(x)$, the equilibrium distributions are described by the Gaussian distributions

$$\begin{equation} \rho_\text{eq} = \frac{e^{-\frac{F^2\beta}{2m\omega^2}}\beta\omega}{2\pi}e^{-\beta\left(\frac{p^2}{2m}+\frac12m\omega^2x^2+Fx\right)}\,. \end{equation} \tag{ 7 }$$

The noise in the nonlinear variable for these equilibrium distributions is given by

$$\begin{equation} \text{Var}\left(p+\lambda x^2\right) = \frac{m}{\beta}+\lambda^2\frac{4F^2\beta+2m\omega^2}{m^3\beta^2\omega^6}\,, \end{equation} \tag{ 8 }$$

irrespective of γ. For such equilibria, if the dynamical frequency ω is allowed to vary such that $\omega\rightarrow\infty$, then the nonlinear variance reduces to a constant noise floor given by $\frac{m}{\beta} = \text{Var}(p)$, representing an equilibrium threshold, independent of λ. No adjustment of the linear potential parameters F or ω can produce equilibrium distributions with $\text{Var}(p+\lambda x^2)$ below this threshold.

However, an appropriate selection of ω, F and transient time t will result in the equilibrium threshold being surpassed by transient conservative linear motion in the region $|\lambda|\lt\sqrt{2\beta}\omega_0^2\left(\frac{m}{3}\right)^\frac32$. Outside this region the noise floor can only be approached by the relaxation process induced by damping. Thus a candidate for the full threshold for noise in the nonlinear variable $p+\lambda x^2$, for noisy damped linear dynamics, is a piecewise function of the noise floor $\frac{m}{\beta}$ and the conservative motion threshold $\mathcal{T}_0$. In order to also accommodate the transient damping effects we must ensure that for no value of γ can the damped transient linear dynamics violate this candidate for the full threshold. Although $\text{Var}(p+\lambda x^2)$ can be expressed analytically (methods detailed in SM) optimization over the parameters γ, F and ω must be performed numerically. Our numerical testing provides evidence that damped noisy linear dynamics cannot surpass these combined thresholds (see SM). The full threshold is therefore finally expressed as

$$\begin{equation} \mathcal{T} = \begin{cases}\frac{3|\lambda|^\frac23}{2^\frac13\left(\beta\omega_0\right)^\frac43} & |\lambda| \lt \sqrt{2\beta}\omega_0^2\left(\frac{m}{3}\right)^\frac32 \\ \frac{m}{\beta} & \text{otherwise}\,. \end{cases} \end{equation} \tag{ 9 }$$

If, while mass m and temperature β are fixed, the statistics in the nonlinear variable $p+\lambda x^2$ go below the threshold in equation (9), then they cannot be achieved by noisy damped up-to-quadratic dynamics following the Markovian model in equation (1), starting from a Gaussian distribution with fiducial frequency ω₀. Additionally, given these fixed parameters, these statistics cannot be represented as a mixture of any Gaussian distributions dynamically available starting from the fiducial distributions and are non-Gaussian in this sense. This approaches refines definitions of classical non-Gaussian states for Markovian dynamics at a fixed temperature T. We observe that the threshold depends only on the mass m for the equilibrium noise threshold and this explains why it is required to be fixed (or at least lower bounded) whenever the transient frequency ω may be freely chosen. Alternatively, one may fix the transient frequency and allow mass to become a varying parameter (see SM). In this case the threshold is identical in form to the conservative threshold $\mathcal{T}_0$, where the fiducial frequency is replaced by $\Omega = \text{max}(\omega,\omega_0)$.

In fact, since the threshold $\mathcal{T}$ is monotonic in these variables, it is sufficient to bound the mass and temperature from below and the fiducial frequency from above. In order to surpass the threshold under these conditions nonlinear stochastic dynamics are required.

To nontrivially overcome the threshold in equation (9) we examine the transient effects of the cubic potential $V_3(x) = k_3x^3$ on an initially Gaussian distribution. The cubic potential is softer then the inverted quadratic potential on one side and harder than quadratic on the other side and so the particle motion quickly diverges in finite time [23, 24] irrespective of the linear damping strength. Therefore we examine the motion in the short transient, starting from a zero-mean Gaussian distribution with $\text{Var}(x_0) = \sigma_{x_0}$, $\text{Var}(p_0) = \sigma_{p_0}$ and $\text{Cov}(x_0,p_0) = \sigma_{x_0p_0} = \langle{x_0p_0}\rangle-\langle{x_0}\rangle\langle{p_0}\rangle$.

Such an initial distribution can be obtained by the linear dynamics of equation (1) from the values of $\sigma_{x_0} = \frac{1}{\beta m\omega_0^2}$, $\sigma_{p_0} = \frac{m}{\beta}$ and $\sigma_{x_0p_0} = 0$ corresponding to the fiducial equilibrium state. Here we distinguish an initial distribution, as a possibly engineered Gaussian distribution prepared for cubic dynamics, and the fiducial equilibrium state. The fiducial distribution is used to define the noise scale and contains our prior assumptions mentioned in the above introduction of linear stochastic motion. Any initial distribution, prepared by linear dynamics from the fiducial distribution, can be used to improve the generation of the nonlinear squeezing without exiting these conditions as it is guaranteed to begin from above the threshold and cannot surpass it through any further linear dynamics.

One of the primary candidates for implementing cubic-like potentials for a free mechanical particle is modern systems of levitating nanoparticles [30] and in what follows we use the parameters $\omega_0 = 120$ kHz, $m = 8\times10^{-17}$ kg and T = 300 K derived from state of the art experiments [15, 16, 31]. In the high pressure regime the Langevin dynamics of equation (1) can be approximated by the overdamped equation

$$\begin{equation} \dot{x} = -\bar{\kappa} x^2 + \sqrt{\frac{2k_{\mathrm {B}}T}{\Gamma}}\xi\left(t\right), \end{equation} \tag{ 10 }$$

where $\bar{\kappa} = k_3/\Gamma$ is the stiffness of the cubic potential normalized to the damping Γ. In this approximation the noise in the nonlinear variable $p+\lambda x^2$ for the short transient, with vanishing mean position $\langle x_0 \rangle = 0$ (see SM for the derivation), reads

$$\begin{align} \text{Var}\left(p+\lambda x^2\right)& \approx \frac{2m^2 k_{\mathrm {B}}T}{\Gamma}+\frac{2m^2\bar{\kappa}^2 \sigma_{x_0}^2}{\Gamma^2} -4m\bar{\kappa}\sigma_{x_0}^2\lambda \left[1 +\frac{6\bar{\kappa}^2 \sigma_{x_0} t^2}{\Gamma^3}\right] \nonumber\\ &\quad +2\lambda^2 \left[\sigma_{x_0}^2 + \frac{4k_B T t \sigma_{x_0}}{\Gamma} \left(1+\frac{k_{\mathrm {B}}Tt}{\Gamma \sigma_{x_0}} \right) + \frac{42 \bar{\kappa}^2 t^2 \sigma_{x_0}^3}{\Gamma^2} \left(1+\frac{4\Gamma k_{\mathrm {B}} T t}{7\sigma_{x_0}} + \frac{8\bar{\kappa}^2 t^2\sigma_{x_0}}{7\Gamma^2} \right) \right].\nonumber \end{align} \tag{ 11 }$$

Note that in the limit of large damping Γ, the surviving term from equation (11) is $2\sigma_{x_0}^2\lambda^2$ which indicates nonlinear squeezing is always below the threshold of linear dynamics (equation (9)). At realistic high pressures of the order of γ = 200 kHz ≈ 200 mbar [32, 33], and evolution in the cubic potential for $t = 9~\mu$s (figure 1 panel (h), red dots), the short time approximation of equation (11) (red solid) accurately predicts a noise suppression below the threshold of linear dynamics (dashed blue). Note how the strong damping decreases the fluctuations of both x and p (panel (d)) compared to the initial thermal distribution (panel (a)). This noise reduction due to damping occurs before the noise begins to diverge due to the cubic dynamics, and aids the observation of nonlinear squeezing. In fact, if the Gaussian approximation to this non-Gaussian distribution is constructed from the covariance matrix it also shows nonlinear squeezing below the threshold, due to this initial noise reduction. In contrast, the overdamped regime with up-to-quadratic potentials does not pass the threshold and therefore, the nonlinear squeezing visible from this non-Gaussian distribution is more correctly attributed to the cubic dynamics rather than the large damping.

Zoom In Zoom Out Reset image size

As time increases (green dots, panel (h)), the nonlinear variance moves closer to the threshold, as the increasing noise in position (panel (e)) hinders the noise suppression. However, it is still below the threshold. This highlights the benefit of large damping as a tool to observe nonlinear squeezing for the first time below the threshold of linear dynamics for slower experiments with timeframes of up to t = 40 ms. Of course, it is more practically important to observe nonlinear squeezing with decreasing damping approaching conservative dynamics, therefore we focus on this case in what follows.

In the levitodynamics operational regime, with $~10^{-11}$ Hz $\unicode{x2A7D} \gamma \unicode{x2A7D} 20$ kHz [15, 31, 34–37], the fiducial thermal distribution does not display nonlinear squeezing below the threshold when evolved for transient times in a cubic potential, as visible in figure 1 (panel (i), red dots). The approximate expression for the noise in the nonlinear variable $p+\lambda x^2$ for an arbitrary initial Gaussian distribution is too large to report here, however for vanishing means ($\langle{x_0}\rangle,\langle{p_0}\rangle = 0$) it corresponds well in the short-time approximation to

$$\begin{align} \text{Var}(p+\lambda x^2)& \approx \sigma_{p_0} (1- m \gamma t )^2 + 2 m \gamma k_{\mathrm {B}} T t + 2m^2\kappa^2\sigma_{x_0}^2t^2 - 2\lambda \left[2m\kappa \sigma_{x_0}^2 t+\frac{\kappa \sigma_{p_0}\sigma_{x_0}t^3}{m} \right.\nonumber\\ &\quad \left. \ + 7 \kappa \sigma_{x_0} t^2 \left(1- \frac{5 m\gamma t}{7} \right) \sigma_{x_0p_0}+ \frac{4 \kappa t^3\sigma_{x_0p_0}^2}{m}\right] + 2\lambda^2 \left[\sigma_{x_0}^2 + \frac{2\sigma_{p_0}\sigma_{x_0}t^2}{m^2}(1- m\gamma t ) \right.\nonumber\\ &\quad \left.\,+ \frac{4 \sigma_{x_0} t}{m} \left( 1-\frac{\gamma t}{2} + \frac{\sigma_{p_0} t^2}{m^2 \sigma_{x_0}} \right)\sigma_{x_0p_0} + \frac{4 t^2}{m^2}(1-m\gamma t )\sigma_{x_0p_0}^2 \right],\nonumber \end{align} \tag{ 12 }$$

where $\kappa = k_3/m$ is the stiffness of the cubic potential k₃ normalized to the mass of the particle. For initial thermal equilibrium distributions we have $\sigma_{x_0p_0} = 0$, and in the low pressure regime with γ < 200 kHz such distribution (red line) does not exhibit nonlinear squeezing in the short time approximation given by equation (12). Outside the short time evolution the system begins to exhibit divergent behavior and the noise in position and momentum increases rapidly, and nonlinear squeezing is not achievable.

In order to observe nonlinear squeezing, a more suitable initial distribution can be prepared from the fiducial distribution, by allowing it to evolve it in a harmonic potential with fixed temperature and damping prior to the cubic dynamics. In general, out-of-equilibrium dynamics involving both a new frequency and a linear force F is required to create such a distribution. Merely displacing the fiducial state or the generation of nonzero correlations σ_xp appears to be insufficient (see SM).

One example of such a sensitive initial distribution can be created by evolving the fiducial equilibrium distribution in a tilted quadratic potential of frequency $\omega = 1.5\omega_0$ and linear force $F = 2k_{\mathrm {B}}T/\mu m$ (figure 1 panel (c)). This state, after transient cubic dynamics, shows nonlinear squeezing (figure 1, yellow dots, panel (i)) below the threshold. Here the full form of equation (12), including arbitrary initial Gaussian distributions, still accurately describes the noise suppression (solid yellow). Such behavior displays a high tolerance to damping (up to γ < 200 kHz) where the short-time approximation holds. The preparation of suitable initial distributions can also be accomplished via postselection of recorded trajectories [16].

Levitated platforms, as well as many NEMS experiments [38–46], more commonly experience nonlinearities in the form of a Duffing oscillator. A Duffing oscillator $V_4(x) = k_2 x^2/2 + k_4 x^4/4$ ($k_4 \ll k_2$) adds a softening character ($k_4 \lt 0$) to the harmonic oscillator. To best access the nonlinearity at transient times, a linear tilt $k_1 x$ is added to the Duffing potential. Extensive numerical simulations in the underdamped regime $\gamma = 10^{-2}$ Hz, performed with $k_1 = 10^{-13}$N [31], $\omega = \sqrt{k_2/m} = 120$ kHz, and $k_4 = 2.7 \mu {\mathrm {m}} ^{-2}$ [16], reveal a nonlinear squeezing below the threshold of equation (9) only for engineered distributions initially squeezed and displaced to negative momentum and positive positions (see SM).

In the opposite regime of hardening character ($k_4 \gt0$), and strong ($k_4 \gg k_2$) quartic bound, the dynamics is dominantly characterized by the linear drift at short transients, where the nonlinear traits are encoded in the higher moments of position. To achieve noise suppression below the threshold equation (9) for such transients $t \approx 6\,\mu$s, an initially squeezed distribution displaced to negative momentum and positive positions is required (see SM).