Subfemtonewton force fields measured with ergodic Brownian ensembles

TL;DR

This paper demonstrates a method to measure and reconstruct radiation pressure force fields with sub-femtonewton resolution using an ergodic ensemble of colloidal particles, enabling ultra-sensitive force measurements in fluids.

Contribution

The authors introduce a novel approach leveraging ergodic Brownian ensembles to achieve unprecedented force measurement resolution in colloidal systems.

Findings

Force resolution of 0.3 fN achieved

Ergodic ensemble allows for high-precision force reconstruction

Method applicable to various force fields in fluids

Abstract

We demonstrate that radiation pressure force fields can be measured and reconstructed with a resolution of $0.3$ fN (at a $99.7\%$ confidence level) using an ergodic ensemble of overdamped colloidal particles. The outstanding force resolution level is provided by the large size of the statistical ensemble built by recording all displacements from all diffusing particles, regardless of trajectory and time. This is only possible because the noise driving the particles is thermal, white and stationary, so that the colloidal system is ergodic, as we carefully verify. Using an ergodic colloidal dispersion for performing ultra-sensitive measurements of external forces is not limited to non-conservative optical force fields. Our experiments therefore give way to interesting opportunities in the context of weak force measurements in fluids.

Equations

$$\gamma {\dot{y}}_i(t)=F(y_i,z_i,t)+\sqrt{2k_{\mathrm{B}}T\gamma }{\xi }_i(t)$$

$$\gamma \frac{\Delta y_i(t_k)}{\Delta t}=F(y_i,z_i,t_k)+\sqrt{2k_{\mathrm{B}}T\gamma }\cdot \frac{w_i(t_k)}{\sqrt{\Delta t}}$$

$$⟨F(y_i,z_i,t_k)⟩=\gamma \frac{⟨\Delta y_i(t_k)⟩}{\Delta t}$$

$$\overline{\delta y_i^2(\Delta )}=\frac{1}{{\mathcal{T}}_i-\Delta }{\int }_0^{{\mathcal{T}}_i-\Delta }(y_i(t^{\prime }+\Delta )-y_i(t^{\prime }){)}^{2}d{t}^{\prime },$$

$$⟨\delta y^2(\Delta )⟩=\frac{1}{N}\sum _i(y_i(t+\Delta )-y_i(t){)}^2$$

$$\underset{\mathcal{T}/\Delta \to \infty }{\lim }\overline{\delta y_i^2(\Delta )}=⟨\overline{\delta y_i^2(\Delta )}⟩$$

$$ϵ(\Delta )=\frac{{\sigma }^2\left(\overline{\delta y_i^2(\Delta )}\right)}{⟨\overline{\delta y_i^2(\Delta )}{⟩}^2}=\frac{⟨{\left(\overline{\delta y_i^2(\Delta )}\right)}^2⟩}{⟨\overline{\delta y_i^2(\Delta )}{⟩}^2}-1$$

$$\underset{\mathcal{T}/\Delta \to \infty }{\lim }ϵ(\Delta )=\frac{4\Delta }{3\mathcal{T}}.$$

$$F_y(z)=F_0\exp \left(-\frac{2(z-z_0{)}^2}{w_0^2}\right),$$

$$F_0⟨\exp \left(-\frac{2(z_n-z_0{)}^2}{w_0^2}\right)⟩\equiv \frac{\gamma }{\Delta t}⟨\Delta y_n⟩.$$

$$⟨F{⟩}_{min}=m\cdot \frac{\gamma }{\Delta t}\sigma \left(⟨\Delta y_n⟩\right)$$

$$\sigma (\Delta y_n)=\sqrt{\frac{2k_{B}T\Delta t}{\gamma }},$$

$$⟨F{⟩}_{min}=\sqrt{2k_{B}T\gamma }\frac{m}{\sqrt{N\Delta t}}.$$

$$\Delta y_i(t_k)=\sqrt{\frac{2k_{B}T\Delta t}{\gamma }}\cdot w_i(t_k),$$

$${\sigma }^2(\tau )=\frac{1}{2}⟨\Delta {\overline{\Omega }}^2(\tau )⟩$$

$${\overline{\Omega }}_k(\tau )=\frac{1}{\tau }{\int }_{t_k}^{t_k+\tau }\Omega (t^{\prime })d{t}^{\prime }$$

$${\sigma }^2(\tau )=\frac{1}{2(N-2n)}\sum _{k=1}^{N-2n}{\left({\overline{\Omega }}_{k+n}(\tau )-{\overline{\Omega }}_k(\tau )\right)}^2.$$

$$\theta (\tau )={\int }_0^{\tau }\Omega (t^{\prime })d{t}^{\prime }$$

$${\sigma }^2(\tau )=\frac{1}{2{\tau }^2}⟨{\left(\theta ((i+2)\tau )-2\theta ((i+1)\tau )+\theta (i\tau )\right)}^2⟩,$$

$$C_{\theta }(\Delta )=\frac{1}{2\pi }{\int }_{-\infty }^{+\infty }{S}_{\theta }(\omega )e^{i\omega \Delta }d\omega .$$

$${\sigma }^2(\tau )$$

$$\displaystyle~{}~{}~{}+(e^{2i\omega\tau}+e^{-2i\omega\tau})-2(e^{i\omega\tau}+e^{-i\omega\tau})\big{]}$$

$$\displaystyle=\frac{1}{4\pi\tau^{2}}\int_{-\infty}^{+\infty}d\omega S_{\theta}(\omega)\big{[}6-8\cos(\omega\tau)+2\cos(2\omega\tau)\big{]}$$

$$=\frac{4}{\pi {\tau }^2}{\int }_{-\infty }^{+\infty }d\omega S_{\theta }(\omega ){\sin }^4\left(\frac{\omega \tau }{2}\right).$$

$${\sigma }^2(\tau )=\frac{4}{\pi {\tau }^2}{\int }_{-\infty }^{+\infty }\frac{d\omega }{{\omega }^2}{S}_{\Omega }(\omega ){\sin }^4\left(\frac{\omega \tau }{2}\right),$$

$$\sigma (\tau )=\sqrt{\frac{2k_{\mathrm{B}}T}{\gamma \tau }}.$$

$${\theta }_k=\Delta t\sum _{j=1}^k{\Omega }_j$$

$${\overline{\Omega }}_k(\tau )=\frac{{\theta }_{k+n}-{\theta }_k}{\tau }$$

$${\overline{\Omega }}_{k+n}(\tau )=\frac{{\theta }_{k+2n}-{\theta }_{k+n}}{\tau }.$$

$${\sigma }^2(\tau )=\frac{1}{2(N-2n){\tau }^2}\sum _{k=1}^{N-2n}({\theta }_{k+2n}-2{\theta }_{k+n}+{\theta }_k{)}^2.$$

$$⟨\overline{\delta y^2(\Delta )}⟩=2D\cdot \left(\Delta -\frac{t_E}{3}\right)+2{\sigma }^2$$

$$\sigma ={\sigma }_0{\left(1+\frac{D{t}_E}{s_0^2}\right)}^{1/2}$$

Full text

I Introduction

That light can exert a pressure force on an illuminated object is one central prediction of Maxwell’s theory Maxwell (1954) that immediately challenged a few experimentalists (see Worrall (1982) and references therein). But their initial attempts in measuring radiation pressure were mostly hindered by thermal effects induced by the illuminating light, such as convective and radiometric forces. The first positive experimental demonstration of light radiation pressure was achieved by P.N. Lebedev in 1899, with results published in a 1901 article Lebedev (1901). In this tour-de-force experiment performed with reflecting winglets suspended on torsional balances, Lebedev was able to measure radiation pressure forces down to $3\times 10^{-10}$ N, with an accuracy better than $6\%$.

Today, the interest in measuring radiation pressure remains. The non-conservative character of radiation pression indeed plays a central role in the generation and control of complex optical force fields He et al. (1995); Cuche et al. (2012); Lehmuskero et al. (2014); Ma et al. (2015); Canaguier-Durand et al. (2013a); Liu et al. (2016); Schnoering et al. (2019), with non-trivial optomechanical effects recently discussed Roichman et al. (2008); Wu et al. (2009); Brzobohatỳ et al. (2013); Gloppe et al. (2014); Canaguier-Durand and Genet (2015); Sukhov et al. (2015); Sukhov and Dogariu (2017); Magallanes and Brasselet (2018); Mangeat et al. (2019). Usually, experiments measuring external force fields involve Brownian probes and monitor, through different means and methods, the shifts of the statistical distributions of displacement of the probes from their stable equilibrium positions induced by the force fields.

In this paper, we show that an ensemble of colloidal particles diffusing in the overdamped regime within an external non-conservative force field can be advantageously used for ultra-sensitive radiation pressure measurements. We demonstrate indeed an unprecedented level of force resolution of $0.3$ fN (within a thermally limited $99.7\%$ confidence interval) by exploiting the large size of the statistical ensemble of Brownian displacements supplied by the colloidal dispersion. We emphasize that such large statistical ensembles gather displacements recorded on different trajectories and at different times, so that our force measurement method crucially relies on the ergodic hypothesis. In order to qualify our colloidal dispersion as an ergodic system, and thus validate our strategy, we demonstrate the thermal and stationary nature of the noise driving our Brownian colloidal ensemble.

One cannot fail to highlight the importance of the ergodic hypothesis which is needed for averaging out the contribution of the thermal when determining the force, as we explain below. The hypothesis is pivotal in many recent experiments, in particular those involving optical traps in fluids, but it is scarcely verified. In this work, we close this loophole by resorting to an appropriate measure of ergodicity. This measure yields a precise confirmation of the ergodic hypothesis for our experiments, all systematic and tracking errors being accounted for.

II Experimental setup

Our experiments consist in illuminating with a horizontal laser beam a colloidal dispersion of micron-sized melamine spheres, diffusing and sedimenting inside a cuvette filled with water. Coming from one side of the cuvette, the laser beam induces a radiation pressure on the particles that modifies their diffusion dynamics. As discussed in details below, this mechanical action can be analyzed by looking at colloidal trajectories in real-time, recorded by tracking the successive positions of the particles, which are fluorescent (weakly dye-doped). The large number of particles leads to collect a large number of trajectories, and thereby to provide a large statistical ensemble of single-step displacements on which a high-precision motional analysis can be performed. To do so, the optical setup, described in details in Fig. 1, has important features.

First, the cuvette has large dimensions compared to the size of the imaging region-of-interest that only extends over a small central region far away from all walls. This, together with the small-volume fraction of the colloidal dispersion allows us to neglect the influence of possible boundary-wall and particle-interaction effects on the colloidal diffusion dynamics. Importantly, such dimensions also insure that sedimentation and possible laser-induced convection effects (see Appendix A) are laminar, only acting on the vertical $z$ axis. These effects thereby are perfectly projected out by looking at Brownian motion along the $y$ axis only. This is a central point in our methodology, as discussed below.

Second, the illumination conditions are set so that, within the imaging region-of-interest, the Gaussian profile of the laser beam is uniform along the horizontal optical $y$ axis (i.e. a Rayleigh length much larger than the width of the cuvette) with a waist much larger than the diameter of a particle. Such conditions, close to plane-wave illumination, minimize any gradient contribution to the optical force field that turns out to be only determined by radiation pressure.

Third, our setup can be operated in two illumination modes. In the single-beam mode, the laser beam illuminates the dispersion from one direction and pushes the particles along the optical axis $y$ –as described in Fig. 1(c). This mode is used for measuring radiation pressure forces, whose effects, along $y$, are decoupled from sedimentation and convection along the $z$ axis. The possibility to access a free diffusion regime is particularly important for carefully assessing the properties of the noise driving the colloidal system and for verifying its ergodic nature. But this regime is unreachable in the single-beam mode where the laser, even at its lowest intensity at the threshold of detection level, always pushes the particles. For this reason, we implement a dual-beam mode that consists in illuminating the dispersion from both sides along the $y$ axis. Because the laser beam is splitted by a polarizing beam splitter, the two counter-propagating beams are crossed polarized. They do not interfere and yield therefore a uniform intensity profile inside the cuvette. The intensities in each beam are carefully balanced so that radiation pressures coming from both sides of the cuvette can be perfectly compensated. In such conditions, the colloidal particles freely diffuse along the $y$ axis.

III Brownian dynamics under radiation pressure force field

With gravity and convection acting along the $z$ axis only, the Brownian motion projected on the $y$ axis for each micron-sized colloidal particle can be described independently from these effects. Therefore, the motion of each particle, subjected inside the cell to the external optical force field $F$ directed along the $y$ axis and depending on both $y,z$ coordinates of the particle 111The dependence on $x$ of the force field can be safely neglected considering that our imaging depth of field $\sim 8~{}\mu$m is much smaller than the beam waist $w_{0}\sim 65~{}\mu$m inside the cuvette. The particles imaged on our camera therefore diffuse along the $x$ axis within a constant laser intensity throughout the experiments, is described by the simple overdamped Langevin equation:

[TABLE]

where $y_{i}(t)$ corresponds to the position of the $i^{\rm th}$ particle measured at a time $t$ along the $y$ axis, $k_{\rm B}$ the Boltzmann constant, $T$ the temperature of water, and $\gamma$ the Stokes friction drag. The stochastic Langevin force is modeled as a Wiener process that satisfies $\langle\xi_{i}(t)\rangle=0$ and $\langle\xi_{i}(t)\xi_{j}(t^{\prime})\rangle=\delta_{ij}\delta(t-t^{\prime})$, with $\delta_{ij}=1$ if $i=j$ and [math] if $i\neq j$, and where $\langle\cdots\rangle$ stands for an ensemble average performed over all the realizations of the stochastic variable $\xi_{i}(t)$.

Experimentally, we acquire images at a given frame rate $f$ and thus implement a discrete version of the equation involving successive displacements $\Delta y_{i}(t_{k})=y_{i}(t_{k+1})-y_{i}(t_{k})$. The discrete Langevin equation reads as:

[TABLE]

where $t_{k}=k\cdot\Delta t$ with $\Delta t=1/f$ and $k$ an integer number. Thermal fluctuations are set by $w_{i}(t_{k})$ as a discrete random number, $\mathcal{N}(0,1)$ normally distributed, which satisfies $\langle w_{i}(t_{k})\rangle=0$ and $\langle w_{i}(t_{k})w_{j}(t_{l})\rangle=\delta_{ij}\delta_{kl}$.

With a statistical distribution of displacements recorded over the ensemble $\{i\}$ of colloidal particles and for all times $\{k\}$ at a given laser intensity, we perform an ensemble averaging of Eq. (2)

[TABLE]

that allows to measure the strength of the force field, as soon as the Stokes friction drag is determined, as discussed in details below.

We emphasize that this simple ensemble average relation relies on averaging out of the thermal fluctuation term in Eq. (2). That $\langle w_{i}(t_{k})\rangle=0$ in the ensemble of all displacements measured at different times and over different trajectories can only be true if one first confirms that all displacements are driven by the same source of thermal (white), stationary noise regardless of the chosen trajectory $i$ and the selected time $k$.

IV Thermal noise and stationarity

In order to verify the thermal and stationary properties of the noise at play in our colloidal system, we illuminate the colloidal dispersion in the dual-beam mode corresponding to two external force fields that compensate each other. Under such conditions, we perform an Allan variance-based analysis on a reconstructed single long trajectory composed of all concatenated displacements $\{\Delta y_{i}(t_{k})\}$ over the ensemble of trajectories $\{i\}$. This calculation, detailed in Appendix B, is performed over 15 different experiments and the result is displayed in Fig. 2(a). The Allan variance clearly shows, in log-log scale, the $-1/2$ exponent expected for a system driven by a white thermal noise Allan (1966). Remarkably, this behavior spans ca. 3 decades of time lag $\Delta$, showing no drift in the experimental thermal noise and thereby indicating that one can exploit the full ensemble of displacements indeed.

As reminded in Appendix B, the Allan variance is related to the displacement covariance $C(\Delta,t)=\langle\Delta y_{\Delta+t}\Delta y_{t}\rangle$ which can also be directly evaluated from all the displacement available in the concatenated trajectory. The result is displayed in Fig. 2(b) showing that $C(\Delta,t)$ fluctuates around a zero mean for all time lags $\Delta>0$, as expected from the white noise Allan deviation. The degree of stochastic independence of each successive steps can be quantified using the correlation of covariance $\eta={C(\Delta)}/{\sigma_{\Delta y}^{2}}$ which, as expected again for a white noise, fluctuates around a zero mean with an amplitude much smaller than $1$.

It also interesting to look at the time ensemble average mean square displacement (MSD) $\langle\overline{\delta y_{i}^{2}(\Delta)}\rangle$ evaluated by, first, averaging over their duration $\mathcal{T}_{i}$ every single trajectory $i$ MSD

[TABLE]

and, then, taking the mean over the ensemble $\{i\}$ of trajectories of such time averages.

The results, evaluated under the same dual-beam mode of illumination for different laser intensities in each beam, are displayed in Fig. 2(c). The clear linearity of $\langle\overline{\delta y_{i}^{2}(\Delta)}\rangle$, for all time lags $\Delta$ corresponds to two important characteristics of the stochastic process driving the particles: first, a zero mean (linearity) and second, a slope $2D$ that measures the diffusion coefficient $D$ of the free Brownian motion along the $y$ axis. As shown in Fig. 2(d) for different laser intensities injected in each path of the dual-beam mode, this measurement falls in good agreement with the value $D=k_{\rm B}T/\gamma$ expected from the Langevin equation, including systematic tracking errors, temperature value and particle size dispersion errors, as discussed in details in Appendices C and D.

The zero mean, fixed variance and covariance independent of the time lag, all together manifest the white stationary character of the noise driving our experiments, in full agreement with the stochastic description of the Langevin force in Eq. (2).

V Ergodicity

This analysis implies that all the successive displacements drawn from different single trajectories at different times within the concatenated trajectory are driven by the same white noise. Under such conditions, the colloidal dispersion must behave like an ergodic system, at the level of which it is possible to collect displacement values acquired from different $i$ trajectories at different $k$ times and to perform large time ensemble averages necessary to evaluate Eq. (3) with great precision. In the context of high resolution force measurements therefore, ergodicity is an important property to verify.

Stationarity implies that the ensemble average mean-square displacement (MSD)

[TABLE]

is independent from the choice of the initial time $t$. As a consequence, $\langle\delta y^{2}(\Delta)\rangle$ equals its time average that simply corresponds to the ensemble mean $\langle\overline{\delta y_{i}^{2}(\Delta)}\rangle$ of single $\{i\}$ trajectory time average MSD given by Eq. (4). This property is clearly seen in Fig. 3 with a ratio $\rho(\Delta)=\langle\overline{\delta y_{i}^{2}(\Delta)}\rangle/\langle\delta y^{2}(\Delta)\rangle$ equal to one for all time lags.

Ergodicity per se then requires any single-trajectory time average MSD $\overline{\delta y_{i}^{2}(\Delta)}$ drawn from the trajectory ensemble $\{i\}$ to become equal to the ensemble mean of such time average MSD

[TABLE]

for a sufficiently long averaging time $\mathcal{T}$. This time corresponds to the full duration of the shortest trajectory recorded and hence defines the duration on which all trajectories are limited or possibly, if long enough, subdivided.

It is clear that this sufficient condition for ergodicity can be quantified by looking at the evolution of the statistical distribution of the single-trajectory MSD $\overline{\delta y_{i}^{2}(\Delta)}$ which, according to Eq. (6), is expected to reach a mean value equal to $\langle\overline{\delta y_{i}^{2}(\Delta)}\rangle$ with a variance $\sigma^{2}\left(\overline{\delta y_{i}^{2}(\Delta)}\right)$ decreasing as $\mathcal{T}/\Delta$ increases. Accordingly, and following Metzler et al. (2014), we introduce the ergodic parameter

[TABLE]

that can be explicitly evaluated for the case of the one-dimensional free Brownian motion expected to take place in our experiments along the optical axis $y$ in the dual-beam mode. The calculation is detailed in Appendix E and leads to the simple ergodic law

[TABLE]

As shown in Fig. 3, this law is clearly followed by our data. We stress in Appendix F that apparent deviations of $\epsilon$ from the ergodic law at short time lags are in fact induced by tracking errors on displacement measurements. Once these errors are accounted for, our measurements eventually agree well with Eq. (8), clearly demonstrating the ergodic nature of our colloidal dispersion.

VI Radiation pressure force measurement and profile reconstruction

White noise stationarity and ergodicity ensure that Eq. (3) is experimentally valid for our system. It can therefore be evaluated by collecting all displacements in order to perform the ensemble averages. Because radiation pressure is determined by fixed laser intensities, and therefore constant in time, $\langle F(y_{i},z_{i},t_{k})\rangle=\langle F(y_{i},z_{i})\rangle$ in Eq. (3). In addition, the Rayleigh range of our Gaussian illumination mode is large compared with the imaging region-of-interest. This gives a radiation pressure force field invariant along the $y$ axis and only dependent on $z$ with an expected Gaussian profile centered on a maximal force $F_{0}$

[TABLE]

where the central position $z_{0}$ of the laser beam and its waist $w_{0}$ are determined optically with a high level of precision (see below). Under such conditions, $\langle F_{i}(y_{i},z_{i})\rangle$ can be evaluated as an average performed on the successive $z$ positions that all diffusing colloids cross within the laser beam while being pushed along the $y$ axis. Experimentally therefore, we record the $(y,z)$ positions of each colloid simultaneously and build two statistical ensembles gathering, one, the successive $n^{\rm th}$ $\Delta y_{n}=y_{n+1}-y_{n}$ displacements measured along the $y$ axis, and the other, the corresponding $z_{n}$ positions of the colloids on the $z$ axis recorded for the $\Delta y_{n}$ displacement. We will note $N$ the total number of such displacements available from all the trajectories and at all times in one experiment.

Such two statistical ensembles give the possibility to estimate the maximal force $F_{0}$ from the agreement between the ensemble averages performed respectively over the ${\{z_{n}\}}$ and ${\{\Delta y_{n}\}}$ ensembles

[TABLE]

The estimators of $F_{0}$ evaluated for different laser intensities in the single-beam configuration are presented in Fig. 4(a), including $\pm\delta F_{0}$ error bars associated with the different sources of errors described in detail in the Appendix G. The results, linearly dependent on the laser intensity, are in good agreement with radiation pressure force values expected from a Mie calculation Canaguier-Durand and Genet (2014).

Once the estimator $F_{0}$ is determined, the full force profile (9) can be reconstructed, as shown in Fig. 4(b) for a given laser intensity. Associated uncertainties here include the residual errors in determining optically the $(z_{0},w_{0})$ values. The waist $w_{0}$ is measured by imaging, with a camera facing the laser, the transverse intensity of the laser beam. Analyzing the profile with a Gaussian model, the waist value ($65~{}\mu$m) is estimated with a relative error of $4.9\%$. The beam center $z_{0}$ is determined by superposing all the frames (gathering $\sim 3\times 10^{6}$ fluorescent spots) from all the experiments, by analyzing the upper and lower limits and the distribution of all points in the imaging field of view, leading to estimate the center position $z_{0}$ with a relative error of $0.2\%$.

A force resolution level of our thermally limited system can be derived, for each experiment providing $N$ displacements, as the minimal measurable force

[TABLE]

where $m$ fixes the confidence interval chosen, and $\sigma\left(\langle\Delta y_{n}\rangle\right)=\sigma\left(\Delta y_{n}\right)/\sqrt{N}$ is the standard error evaluated from the variance $\sigma\left(\Delta y_{n}\right)$ of the displacement ensemble $\{\Delta y_{n}\}$. This variance can be evaluated directly from Eq. (2) in the absence of any force field as

[TABLE]

giving a resolution level of our force measurement as:

[TABLE]

As clear from this relation, the resolution can be improved by increasing the size of the statistical ensemble (i.e. the number $N$ of displacements) within the thermal, stationary and ergodic limits set in Sec. B. For our experiments, the thermally limited force sensitivity is $\sqrt{2k_{B}T\gamma}\simeq 8.56$ N$/\sqrt{\rm Hz}$. Within such limits, a typical experiment yields $N\sim 2\times 10^{5}$ successive displacements, which corresponds to a remarkable $\sim 0.3$ fN resolution level (within a $m=3$ $99.7\%$ confidence interval). Such levels of resolution are displayed in Fig. 4(a) and they vary just as $N$ depends on the experiment performed at a given laser intensity. These data clearly show that in the low laser power regime, radiation pressure can be determined down to the thermal limit, with a $\sim 0.5$ fN force actually measured at a $99.7\%$ confidence level.

Finally, we exploit the fact that video-tracking microscopy also gives the possibility to measure the external optical force field directly within limited horizontal layers chosen from both sides of the optical axis of the illumination beam. We define $18.4~{}\mu$m thick layers at given $z$ positions and within which we collect all displacements along the $y$ direction. We use again Eq. (3) inside such selected layers to give a measurement of the average radiation pressure force zone-by-zone. The results are reported in Fig. 5 and falls in very good agreement with the reconstruction method.

VII Conclusions

Our results demonstrate that colloidal dispersions can be used as highly sensitive Brownian probes for measuring external force fields. Because such systems remain thermally limited, with stationary and ergodic dynamics, large statistical ensembles mixing particle displacements recorded from different trajectories and times become available and thus make sub-fN force resolution levels accessible. The main experimental uncertainty stems in our work from the particle size dispersion in the colloidal dispersion commercially available. But this source of systematic error could be reduced working with better monodisperse colloidal dispersions. Our method of statistical reconstruction of the force profile, here demonstrated in two-dimensions for a non-conservative optical force field, can obviously be extended in three-dimensions, involving efficient and available techniques in three-dimensional tracking Martinez-Marrades et al. (2014); Verrier et al. (2015). Validating colloidal systems for ultrasensitive force detection strategies, our work opens interesting perspectives. For instance, colloids can be involved in reconstructing complex topological force fields, in particular in the near field where momentum exchanges are enhanced by near-field inhomogeneities Cuche et al. (2012); Canaguier-Durand et al. (2013a). As recently proposed, colloidal systems can also become pertinent tools for weak force measurements in the context of Casimir physics Ether et al. (2015). The advent of nano- and micro-structured and functionalized colloids McPeak et al. (2014); Kumar et al. (2016); Ma et al. (2017) can lead to new types of dynamical responses to external fields, as exemplified with chiral optical forces Canaguier-Durand et al. (2013b); Cameron et al. (2014); Ding et al. (2014); Zhao et al. (2017); Kravets et al. (2019a, b). The outstanding stability of the statistical properties of our system offers new possibilities for deciphering non-trivial force fields at a genuine sub-fN resolution level. We finally stress that our approach is also relevant for weak force experiments that do not necessarily involve optical fields, experiments found for instance in the burgeoning field of mechanical chiral resolution Marichez et al. (2019).

VIII Acknowledgments

This work was supported in part by Agence Nationale de la Recherche (ANR), France, ANR Equipex Union (Grant No. ANR-10-EQPX-52-01), the Labex NIE projects (Grant No. ANR-11-LABX-0058-NIE), and USIAS within the Investissements d’Avenir program (Grant No. ANR-10-IDEX-0002-02).

IX Appendix A: Sample preparation

The samples are prepared from an initial dispersion ($2.5\%$ mass-volume ratio) of melamine micro-particles of diameter $d=0.940\pm 0.05~{}\mu$m purchased from Micro-Particles GmbH, weakly doped with a fluorescent dye for a most efficient detection in water. We dilute the dispersion $\sim 10^{4}\times$ with ultra-pure water and fill a cuvette with the colloidal dispersion to a $1.3\times 10^{-6}$ low-volume fraction. The dimensions of the cuvette –see Fig. 1 in the main text– are large enough so that within imaging region-of-interests, boundary wall effects can be safely neglected. The filled cuvette is covered and sealed with vacuum grease in order to prevent water evaporation and to isolate the fluid from other environmental influence. The cuvette and its cover are exposed at least one hour to UV light in order to ensure the absence of any bacterial contaminant. The sample is grounded in order to remove any electrostatic charge on the surface of the cuvette. Before performing our experiments, we leave for about $1$ hour the sample relaxing in its holder until well thermalized with the environment.

One sees in our system some convective drag effects at play, induced by the illuminating laser depending on its intensity. These effects can be simply understood. When the laser is turned on, the water is heated. Although this effect is minute, it is sufficient for reducing the density of water within the laser beam. As a consequence, buoyancy driven flows are induced that can eventually drag the particles upwards at the highest laser intensities. This is what is seen, for instance, in Fig. 1 (c) where the convection induced at high laser intensity drags the melamine spheres against sedimentation. The large dimensions of the cuvette insure that such a convective drag is strictly laminar and performed along the $z>0$ axis. In such conditions, both sedimentation and convective effects are perfectly decoupled from the diffusive dynamics performed on the optical $y$ axis, along which external force fields are measured.

It is interesting to note that for smaller laser intensities, the drag resulting from sedimentation and convection is directed downwards. This means that it is actually possible to find a condition of illumination where sedimentation and convective drags can compensate each other. In the dual-beam illumination mode, this condition corresponds to a laser intensity of $0.148$ kW/cm2 in each beam. Remarkably then, the Brownian motion of the ensemble becomes totally free in three dimensions within the imaging field of view.

X Appendix B: Allan variance analysis of the noise

The Allan variance is a statistical tool particularly relevant for quantifying noise sources because of its direct relation with noise power spectral densities Allan (1966); Barnes et al. (1971) . It has recently been used in the context of optical trapping in order to assess limits of stability of traps and to provide optimal measurement bandwidths Gibson et al. (2008); Czerwinski et al. (2009); Schnoering et al. (2019).

We first explain how the Allan variance is measured in our experiments, performed in the dual-beam mode with two laser beams of equal intensity illuminating the dispersion. This leads to a precise compensation of external radiation pressure forces –i.e. to free Brownian motion along the $y$ axis, described by the simple Langevin equation

[TABLE]

only driven by the noise $w_{i}(t_{k})$. In such conditions, the successive displacements $\Delta y_{i}(t_{k})$ recorded from all the trajectories (over the shortest observation time lag $\Delta t$) are concatenated one-by-one in order to form a single time series $\Omega(\mathcal{\tau})$ of displacements with a total time $\mathcal{T}=N\Delta t>\tau$.

The Allan variance

[TABLE]

is calculated on this long time series using successive and non-overlapping series of shorter duration $\tau$ with $\Delta\overline{\Omega}(\tau)=\overline{\Omega}_{k+n}(\tau)-\overline{\Omega}_{k}(\tau)$ and for which time average values $\overline{\Omega}_{k}(\tau)$ of $\Omega(\tau)$ are defined as:

[TABLE]

with $t_{k}=k\Delta t$ and $\tau=n\Delta t$. Here, $\langle\cdots\rangle$ defines an average performed over all accessible series, implying that $n<N/2$ and $k\leq N-2n$. The Allan variance is then explicitly calculated as:

[TABLE]

The connection to noise power spectral density can be seen by defining an intermediate estimator $\theta$

[TABLE]

that corresponds in Allan (1966) to instantaneous phase fluctuations for $\Omega(\tau)$ corresponding to instantaneous frequency fluctuations. With this estimator, the Allan variance can be written as

[TABLE]

which gives a sum of different covariance function $C_{\theta}(\Delta)=\langle\theta(\tau+\Delta)\theta(\tau)\rangle$ related to a power spectrum densitiy (PSD) as

[TABLE]

Eq. (19) can thus be simply expressed as

[TABLE]

Considering the relation $\omega^{2}S_{\theta}(\omega)=S_{\Omega}(\omega)$ between the PSD of the noise and the PSD of its estimator, we have

[TABLE]

as the important relation that shows why the Allan variance is an appropriate tool for quantifying all kind of noise sources characterized by their respective PSD.

Let us evaluate the Allan variance for the case of a white noise $\Omega(t)=w(t)$ characterized by the well-known covariance $\langle w(t)w(0)\rangle=h_{0}\delta(t)$ and the associated PSD $S_{w}(\omega)=h_{0}$. For such a PSD, we can directly evaluate Eq. (22) to find that $\sigma^{2}(\tau)=\frac{h_{0}}{\tau}$. For the thermal white noise driving our experiments, $h_{0}=2k_{\rm B}T/\gamma$, leading to

[TABLE]

This type of Allan variance is clearly verified in our experiments (Fig. 2 (a), main text) with a linear evolution observed in the log-log scale with a $-1/2$ slope and an amplitude fixed at $2k_{\rm B}T/\gamma$.

Experimentally, the Allan variance is calculated from a discretized version of the estimator

[TABLE]

where $\theta_{k}=\theta(k\Delta t)$ and $\Omega_{j}=\Omega(j\Delta t)$. The time averages of $\Omega$ –Eq.(16)– performed over the temporal length of the trajectory $\tau=n\Delta t$ for two successive and non-overlapping displacement samples, is then simply expressed as:

[TABLE]

[TABLE]

Following Eq. (17), we therefore calculate in this work the Allan variance as:

[TABLE]

XI Appendix C: Tracking error data correction

As discussed in Michalet and Berglund (2012), recording the trajectory of a diffusing particle under the microscope is accompanied by tracking errors that can affect the Brownian statistical analysis. Such errors combine localization errors for single-shot positional measurements and blurring effects caused by the finite exposure time $t_{E}$ of such a measurement over which the particle keeps diffusing with a diffusion coefficient $D$. For free Brownian motion, the one-dimensional time average mean-square displacement (MSD) acquired along a given axis (the $y$ axis, following our experimental setup) has to be corrected. In order to include these effects, the MSD writes as:

[TABLE]

where $t_{E}$ is the exposure time of the camera and $\sigma$ is the dynamic localization error. This error

[TABLE]

is related to the static localization error $\sigma_{0}$ and the full-width at half maximum $s_{0}$ of the Gaussian profile used to approximate the microscope point-spread function.

For our experimental setup, we have: $\sigma_{0}\sim 15$ nm and $s_{0}\sim 1.5~{}\mu$m, extending over $3$ pixels of the CCD camera. Our exposure time is fixed at $t_{E}=7$ ms. Using the expected value of $D$ in water at room temperature for a particle of diameter $940$ nm, these values implies that $s_{0}^{2}\gg Dt_{E}$, so that $\sigma\sim\sigma_{0}$.

Therefore, the single-particle time average MSD $\overline{\delta y_{i}^{2}(\Delta)}_{exp}$ experimentally acquired

[TABLE]

departs from the real MSD $\overline{\delta y_{i}^{2}(\Delta)}_{0}$ by a random localization error $\varepsilon_{i}$. At the level of a colloidal dispersion, we perform an additional ensemble averaging giving

[TABLE]

with $\langle\overline{\delta y_{i}^{2}(\Delta)}\rangle_{0}=2D\Delta$ and $b=\langle\varepsilon_{i}\rangle=2\sigma_{0}^{2}-2Dt_{E}/3$. Within our experimental conditions, this error parameter is estimated to be $b\sim-2.6\times 10^{-3}~{}\mu{\rm m}^{2}$, negative as expected looking at the experimental MSD data point at small $\Delta$.

Such corrections impact the determination of the experimental MSD under an external force field, hence the diffusion coefficient, and the determination of the ergodicity parameter $\epsilon$ that can gives the impression, if not corrected, to display signatures of ergodicity breaking at short time scales, as discussed in Sec. F in details.

We can also get an estimate of the tracking error by linearly fitting the 15 experimental time ensemble average MSD. We obtain an average error parameter $b_{\rm fit}\sim-2.8\times 10^{-3}~{}\mu{\rm m}^{2}$ including a standard deviation of $\sim 0.4\times 10^{-3}~{}\mu{\rm m}^{2}$, showing that $b_{\rm fit}$ is in good agreement with the estimated value $b$.

XII Appendix D: Determination of diffusion coefficients

Determining diffusion coefficients is central to our experiments since it allows to extract a value for the Stokes friction drag $\gamma$ once the temperature of the system is known. Under an external force field, the diffusion coefficient is related to the experimental variance

[TABLE]

given by the statistical ensemble of displacements measured along one $r$ axis with a time difference $\Delta$, including the tracking error $b$ discussed above in Eq. (31).

This variance can be calculated from the difference between the experimental time ensemble average MSD and the square of the mean displacement measured over $\Delta$ as:

[TABLE]

where one recovers on $\langle\overline{\delta r_{i}^{2}(\Delta)}\rangle$ the ballistic contribution of the external force field.

Fitting the variance respectively in the $y$ and $z$ directions with a chosen time lag $\Delta$ measures the diffusion coefficient $D_{y}$ and $D_{z}$. The results are shown in Fig. 6 for different illumination powers, hence different strengths of radiation pressure exerted along the $y$ axis. The results clearly show, within error bars, that $D_{y}\sim D_{z}$ and that they fall in good agreement, within error bars (see caption), with the theoretical value of the diffusion coefficient $D=k_{B}T/\gamma$ expected with a measured temperature $T=302.15$ K and a Stokes friction drag $\gamma$ evaluated from the corresponding viscosity of water and the known diameter of the colloidal sphere. This result is important since it provides a valid value (including uncertainties) for the Stokes friction drag, necessary in Eq. (3) for converting measured spatial displacements into force strength signals within well determined confidence intervals.

XIII Appendix E: Ergodicity for free Brownian motion

We here consider a one-dimensional free Brownian motion $x(t)$ described, in the overdamped limit, by the simple Langevin equation:

[TABLE]

where $U(t)$ corresponds to a thermal noise with $\langle U(t)\rangle=0$ and $\langle U(t_{1})U(t_{2})\rangle=\frac{2k_{B}T}{\gamma}\delta(t_{1}-t_{2})$ where $k_{\rm B}$ is the Boltzmann constant, $T$ the temperature of the noise, and $\gamma$ the Stokes friction drag acting on the Brownian object considered. The ensemble average $\langle\cdots\rangle$ is done over all the realizations of the stochastic variable. Having:

[TABLE]

and setting $x(0)=0$, the mean-square displacement (MSD) is given by the Einstein equation:

[TABLE]

The position correlation function $C_{x}(t_{1},t_{2})$ can then be calculated:

[TABLE]

If we suppose that $t_{2}>t_{1}$:

[TABLE]

since the two integrals do not overlap.

Similarly, if $t_{1}>t_{2}$, then $C_{x}(t_{1},t_{2})=2Dt_{2}$. Both cases can be described in a unified manner with:

[TABLE]

We now calculate explicitly the ergodic parameter $\epsilon$ defined in Eq.(7) in the main text. First, we remind that $\langle x(\Delta)^{2}\rangle=\langle\overline{\delta x^{2}(\Delta)}\rangle$ where $\overline{\delta x^{2}(\Delta)}$ is the time average MSD taken over the duration $\mathcal{T}$

[TABLE]

and $\langle\overline{\delta x^{2}(\Delta)}\rangle=2D\Delta$ the corresponding ensemble average performed over all available $\overline{\delta x^{2}(\Delta)}$.

The ergodic parameter $\epsilon$ is built on the variance of time average MSD $\sigma^{2}(\overline{\delta x^{2}(\Delta)})=\langle\overline{\delta x^{2}(\Delta)}^{2}\rangle-\langle\overline{\delta x^{2}(\Delta)}\rangle^{2}$. Looking at

[TABLE]

we will use Wick’s relation for normally distributed ensembles of displacements:

[TABLE]

Using Eq. (36) and Eq. (39), one gets:

[TABLE]

with:

[TABLE]

The term $4D^{2}\Delta^{2}$ is cancelled in the variance by $\langle\overline{\delta x^{2}(\Delta)}\rangle^{2}$ and one is therefore left to calculate

[TABLE]

To do so, we use the change of variables described in Fig. 7 that splits the integration over the two distinct sectors $t_{2}<t_{1}$ and $t_{2}>t_{1}$. For the $t_{2}<t_{1}$ sector, we have

[TABLE]

and for the $t_{2}>t_{1}$ sector

[TABLE]

The absolute value integrands are finally evaluated using a last change of variable $v=t^{\prime}+\Delta$ on Eq.(46) and $w=t^{\prime}-\Delta$ on Eq.(47) that restricts the integration to $v>0$ and $w<0$ values, and eventually leads to:

[TABLE]

In the long-time limit $\mathcal{T}/\Delta\rightarrow\infty$:

[TABLE]

so that Eq.(45) simplifies to:

[TABLE]

The ergodic parameter $\epsilon=\sigma^{2}\left(\overline{\delta x^{2}(\Delta)}\right)/\langle\overline{\delta x^{2}(\Delta)}\rangle^{2}$ then simply writes in the long-time limit as:

[TABLE]

This scaling is discussed already in Metzler et al. (2014) as a sufficient condition for ergodicity. As we demonstrate in the main text, this scaling is clearly verified in our experiments, proving the ergodic character of our colloidal system.

XIV Appendix F: Correcting the ergodic parameter from tracking errors

With such tracking-error corrections, we can give a relation between the experimentally estimated ergodicity breaking parameter $\epsilon(\Delta)_{exp}$ and the real one $\epsilon(\Delta)_{0}$, starting from the definition of the $EB$ parameter that involves the MSD according to:

[TABLE]

Defining the ratio

[TABLE]

which, from Eq. (29), can be written as a function of $\Delta$:

[TABLE]

allowing us to rewrite Eq. (52) as:

[TABLE]

Assuming that $\varepsilon_{i}$ and $\overline{\delta y_{i}^{2}(\Delta)}_{0}$ are not correlated leads to $\langle\varepsilon_{i}\overline{\delta y_{i}^{2}(\Delta)}_{0}\rangle=b\langle\overline{\delta y_{i}^{2}(\Delta)}_{0}\rangle$, so that:

[TABLE]

considering that:

[TABLE]

This relation has an important consequence when discussing the time evolution of the ergodic parameter $\epsilon$. Indeed, as seen from Eq. (56), for sufficiently large time lags $\Delta\gg\Delta t$ with respect to the sampling time, the error term vanishes as $\phi\rightarrow 1$. In such conditions, $\epsilon(\Delta)_{exp}\sim\epsilon(\Delta)_{0}$. But when evaluating $\epsilon$ in the regime of small time lags with $\Delta\sim\Delta t$, the experimental result $\epsilon(\Delta)_{exp}$ turns out to be larger than $\epsilon(\Delta)_{0}$ precisely because of the influence of the tracking errors. What can be taken for the signature of some breaking of ergodicity at small time lags is simply, in this case, an artefact of the experiment which is simply accounted for using Eq. (56).

In order to evaluate Eq. (56) as presented in Fig. 2 (b) in the main text, we estimate $\phi$ using Eq. (54) and evaluate $\sigma^{2}(\varepsilon_{i})$ in two steps. First, we assume that $\overline{\delta y_{i}^{2}(\Delta)}_{exp}$ and $\overline{\delta y_{i}^{2}(\Delta)}_{0}$ are uncorrelated, so that according to Eq. (30):

[TABLE]

We then take Eq. (7) in the smallest $\Delta\sim\Delta t$ limit, for which $\sigma^{2}\left(\overline{\delta y_{i}^{2}(\Delta t)_{0}}\right)\sim\langle\overline{\delta y_{i}^{2}(\Delta t)_{0}}\rangle^{2}\times 4\Delta t/3\mathcal{T}\sim(\Delta t)^{2}/\mathcal{T}\sim 0$ considering that $\langle\overline{\delta y_{i}^{2}(\Delta)_{0}}\rangle^{2}=2D\Delta t$. This leads eventually to:

[TABLE]

XV Appendix G: Radiation pressure force uncertainty

The global uncertainty related to a radiation pressure force measurement $F_{0}$ involves different sources of errors: $(i)$ of systematic nature regarding the determination of the colloidal particle size, the temperature (hence viscosity) of water and the laser beam waist $w_{0}$ and center $z_{0}$, and $(ii)$ of statistical nature related to the standard error associated with the ensemble averaging of the displacement distribution. These different errors combine in a global uncertainty as:

[TABLE]

From Eq. (10), we know that the $F_{0}$ can be written as:

[TABLE]

where $\alpha=\langle\exp\left(-2(z_{n}-z_{0})^{2}/w_{0}^{2}\right)\rangle$. The error on $\alpha$ comes from the determination error of $z_{0}$ and $w_{0}$, which yields a relative error of $\varphi_{\alpha}=\delta\alpha/\alpha\sim 0.4\%$. Within the Stokes drag term $\gamma=6\pi\eta a$, both systematic errors are at play. The $\pm 1$ K error on temperature around a mean value of $302.15$ K leads to a relative error in the viscosity $\varphi_{\eta}\sim 2.4\%$, and the $50$ nm standard deviation in particle size given by the manufacturer corresponds to a relative error in particle diameter of $\varphi_{a}\sim 5.3\%$. The standard error in the displacement distribution is calculated from the ensemble of $\{\Delta y_{n}\}$, so that, using the maximum error estimation method, the uncertainty writes as follows:

[TABLE]

where $\sigma(\Delta y_{n})$ is the sample standard deviation defined as:

[TABLE]

XVI Appendix H: Mie calculation uncertainty

Following Canaguier-Durand and Genet (2014), radiation pressure forces are calculated in the Mie regime as $F_{pr}=(n_{w}I/c)\cdot Q_{pr}$, where $n_{w}$ is the refractive index of water, $I$ is the intensity of the illuminating laser, $c$ the velocity of light and $Q_{pr}$ is the calculated pressure cross section. Here, the uncertainty comes from the laser intensity because of inevitable errors made in measuring $(i)$ the laser intensity, and $(ii)$ the optical waist of the laser beam inside the cuvette. With $I=2P/\pi w_{0}^{2}$, the laser intensity uncertainty is therefore given by:

[TABLE]

In our experiments, we have a relative error in the power determination of ca. $5\%$ and a relative error in the waist measurement of ca. $4.9\%$ through the Gaussian fitting of the intensity profile with $95\%$ confidence level (see main text).