Stability of long-sustained oscillations induced by electron tunneling

.


I. INTRODUCTION
Coupling of quantum devices to mechanical degrees of freedom can be exploited for high-precision measurements [1][2][3][4][5][6] and may serve as a platform for quantum and classical information processing [7,8].At low temperatures, carbon nanotube (CNT) devices can be operated as extremely sensitive mechanical oscillators which are strongly coupled to single electron tunneling [4,[9][10][11][12][13][14][15][16][17][18].The interplay between single electron tunneling and mechanical motion, in the absence of a mechanical drive, can give rise to self-sustained oscillations.Such self-oscillations were observed to be either present or absent depending on the electron transport regime, both in theory [19][20][21] and experiments [22][23][24][25][26][27][28].In this paper, we report that, at the boundary between different electron transport regimes, self-oscillations can appear or vanish spontaneously.These bistable states of motion of an undriven oscillator, until now unexplored, are of particular interest for applications of these devices in quantum and stochastic thermodynamics [29][30][31][32], in collective dynamics and synchronization [33], and in emulating neural behavior [34][35][36].We measure bistable self-oscillations both for single and double quantum dot configurations and present a theoretical analysis that provides a complete characterization of the stability of self-oscillations at different bias voltages.We find that, once started, undriven oscillations can be self-sustained and may decay over time scales of the order of 10 8 mechanical periods.In this way, our system explores timescales of electronic and mechanical origin that are separated by several orders of magnitude.
Our electromechanical device consists of a fully suspended CNT (see Fig. 1(a)) [14,26,37,38].Applying a bias voltage V SD between the source and drain electrodes, we measure a current I through the CNT.The gate electrodes, to which we apply voltages V G1-5 , are located beneath the CNT.These gate voltages define an electrostatic potential for the confined charges, and depending on the combination of gate voltage values, we can define single or double quantum dots within the CNT (Figs. 1(a,b)) [39,40].
The mechanical motion of the CNT can be excited by applying a radio-frequency (rf) signal to one of the gate electrodes.Sharp changes in the current through the CNT indicate its mechanical motion [9,41,42], and we identify the natural mechanical resonance frequency Ω/2π = 270 MHz with an approximate quality factor Q of 2000 (see Supplemental Material [43]).Experiments were performed at 60 mK.At charge degeneracy points, strong coupling between electron tunnelling and mechanical motion is evidenced by the softening of the mechanical resonance frequency [10-12, 14, 22, 44].
Self-oscillations are identified as sharp switches in current appearing in the absence of an rf excitation [22][23][24][25][26]28].We observe such switches both for single and double dot configurations (Fig. 1(c,d)).While in a double-dot configuration, these sharp switches are visible within the bias triangles (Fig. 1(d)).

A. Single-dot configuration
In the single-dot configuration, we observe the bistability through the hysteretic behavior of the current when sweeping -180 560 680 V SD in and out of the self-oscillation area, following the vertical white dashed line in Fig. 1(c), which corresponds to a constant gate voltage V G3 = 1.8 V.As V SD is swept upward, a sharp increase of current around V SD = 2.3 mV indicates the onset of self-oscillations (red curve in Fig. 2(a)).However, when V SD is swept in the opposite direction (blue curve), the sharp decrease of current is found at a significantly lower voltage V SD ≃ 1.3 mV.The sharp changes in current are reproducible over several sweeps of V SD , with a small variation in threshold voltages due to the stochastic nature of this process.This hysteresis cycle shows that there is a range of bias voltages, V SD ≃ 1.3 − 2.3 mV, labeled as (II) in Fig. 2, for which the system exhibits bistability.
We can explain the bistability regime observed using a model describing the motion of the CNT as a single vibrational mode of displacement x(t), frequency Ω, and effective mass m, whose evolution equation reads [19,20,30] Here γ = mΩ/Q is the friction coefficient affecting the CNT motion and n(t) = 0, 1 is the occupation number of the dot, which is a stochastic variable, and ξ(t) is thermal Gaussian white noise with zero average and ⟨ξ(t)ξ(t ′ )⟩ = 2γk B T δ(t − t ′ ), T being the CNT's temperature and k B the Boltzmann constant.Gate voltages exert a force on the CNT when the dot is occupied, i.e. n(t) = 1.This force depends on several parameters like the distance between the dot and the gate electrodes.For small oscillations, this force can be considered constant.The constant g m determines the strength of the coupling between the dot charge and the mechanical degree of freedom , arising from the electrostatic potential induced by the gates [14].
The stochastic occupation number n(t) undergoes Poissonian jumps when electron tunneling occurs [20,30].The rates of jumps from the left and right electrodes to the dot are given by Γ L,R (x)f L,R (ϵ(x)) and from the dot to the electrodes by where Γ L,R (x) are the tunneling rates, and f L,R (ϵ(x)) are the Fermi functions of each lead evaluated at the energy ϵ(x).The electrochemical potential of the dot, ϵ(x) = ϵ 0 + g m x, depends on the displacement x of the oscillation.The constant ϵ 0 is the electrochemical potential of the dot in the absence of any mechanical motion.Assuming that the energy of the dot reaches the chemical potential µ S at the onset of self-oscillations and transport, ϵ 0 = eV * S , where V * S is the value of V SD at the border between regions (II) and (III) in Fig. 2(a), i.e., ϵ 0 = 2.25 meV.The value of g m can be estimated as in Ref. [14].Notice that the tunnelling rates Γ L,R (x) depend on the energy of the dot and, consequently, on the displacement x of the oscillations.This inhomogeneity of the tunnel barriers have been found to be a necessary condition for the occurrence of self-oscillations [20,21,26,30].
In our experiments, the tunnelling rates are much larger than the mechanical frequency of the CNT, Γ L,R (x) ≫ Ω.Hence, the dynamics of the random variable n(t) is much faster than the motion of the CNT.Moreover, the thermal noise is small (see Supplemental Material [43] for a more detailed discussion of the effect of temperature in the system).This enables us to approximate the behaviour of the system through deterministic dynamics that are influenced by minor fluctuations arising from thermal noise and the random nature of the occupation number n(t).The deterministic dynamics determines the stability of self-oscillations, as explained below, whereas fluctuations are responsible for their generation and decay [27] (see also Appendices A and B).If we neglect the thermal noise in Eq. ( 1) and replace the random occupation number n(t) by its average n(t) we obtain: Equation ( 3) is the master equation for the average occupation number with transition rates 2) and ( 3) are deterministic and predict the appearance of self-oscillations [26].
Here we analyze the stability of self-oscillations calculating the exchange of energy between the quantum dot and the CNT.
To simplify the analysis, we further assume that the leads are at zero temperature, yielding sharp Fermi functions f R,L (ϵ) and that the tunneling rates are given by Γ L,R (x) = Γ L,R e αL,Rx [12], α L,R being parameters that quantify the inhomogeneity of the barriers.Self-oscillations were so far explained by proving that the last term in Eq. ( 2) acts as negative damping that can counterbalance the energy dissipation term [26].An insightful approach to understanding the stability of self-oscillations is to consider the energy ∆E that the CNT gains in an oscillation of a given amplitude A. The mechanical energy, see Appendix A 2 for details.Notice that the device functions as an engine: the second term in Eq.( 4) is the energy that the electrical charge transfers to the mechanical motion by electrostatic interaction, whereas the first term is the energy dissipated as heat through friction.
To calculate ∆E, it is enough to consider harmonic oscillations x(t) = A cos(Ωt).We solve Eq. ( 2) for the average occupation number n(t) using this harmonic approximation and then insert the solution into Eq.( 4) to calculate ∆E as a function of A and the other relevant parameters.This approximation is accurate for a high Q and a relatively small g m .The stability of an oscillation of amplitude A is given by the sign of ∆E.If the oscillation gains energy, that is, if ∆E > 0, then the amplitude increases.On the other hand, if ∆E < 0, the amplitude decreases.The magnitude of ∆E does not affect the stability of the oscillations, but it is relevant for the spontaneous generation and decay of the self-oscillations (see Appendix B).In Fig. 2(b) we show the areas where ∆E is positive (dark blue) or negative (light blue), depending on the value of the bias voltage V SD and the amplitude A of the oscillation.For the rest of the parameters, we use realistic values based on a fitting of the Coulomb peaks (see Appendix C).If the system is oscillating with an amplitude A located in the dark blue region, then the oscillator gains energy in each oscillation and the amplitude increases, as indicated by the vertical arrows, until it reaches the light blue region.On the other hand, the amplitude of the oscillation decreases in the light blue region.We can also distinguish three different regions depending on V SD .In region (I), ∆E is negative for any value of A; therefore, oscillations lose energy and fade out.In region (II), there are both positive and negative values of ∆E, in correspondence with the bistability observed in Fig. 2(a).In region (III), ∆E is mostly positive, and A reaches a saturation value A ≈ 0.25 nm.At the boundary between regions (II) and (III), for small values of A, the model also predicts unstable amplitudes.The shape of the dark blue region is given by dissipation.
Our model also provides an estimate of I, We estimate I for the x(t) resulting from the saturation value of A. The result, plotted in Fig. 2(a) (orange dashed curve), shows good agreement with the data (blue and red curves).We also explore the occurrence and duration of selfoscillations in region (II).The sample is subjected to a protocol where V SD is modified in sequence as shown in Fig. 3(a).We start the protocol with a low voltage, V KILL = 0 mV, for which self-oscillations are absent, and perform a rapid quench to the value V SD = V PROBE that we want to probe.This voltage is kept constant for about 10 seconds to observe the spontaneous onset of self-oscillations revealed by a sudden increase in I. V SD is then changed for about a second to a high pump bias voltage V PUMP = 3 mV, where self-oscillations are induced.We then move back V SD to V PROBE for the rest of the sequence (about 20 s) to measure the persistence of the selfoscillations.Fig. 3(a) shows I during the protocol for three representative regions.Sudden changes in the current indicate the presence of self-oscillations.
Fig. 3(b) summarizes our experimental results for a wide range of V PROBE .We identify four different regimes associated with the three regions in Fig. 2. For 0 < V PROBE < 1.3 mV, region (I), the absence of current indicates that self-oscillations are not stable.They do not appear spontaneously at V SD = V PROBE and vanish immediately after the pumping step.In re- gion (IIa), 1.3 mV < V PROBE < 2.0 mV, self-oscillations do not start spontaneously but endure during a random time after being triggered by the pumping step.In a small bias voltage range around V PROBE ≈ 2.0 mV, region labelled as (IIb), selfoscillations can start spontaneously after some time at V PROBE .Finally, for V PROBE > 2 mV, region (III), self-oscillations are always stable: they start spontaneously at V PROBE and are maintained during the whole protocol.
In the bistable regions, (IIa) and (IIb), self-oscillations appear and vanish due to fluctuations that allow the system to access areas in the stability diagram of Fig. 2(b) that have the opposite sign of ∆E than that dictated by our deterministic model.The probability of these excursions is very small but they can occur after a large number of oscillations.For instance, a self-oscillation of frequency Ω/2π ∼ 270 MHz can last for 20 seconds or 5.4 × 10 9 oscillations.This mechanism explains the huge separation of time scales in the device, as explained in detail in Appendix B.

B. Double-dot configuration
We perform a similar study when the device is in the double-dot configuration by measuring the hysteresis as a function of V SD (Fig. 4(a)) at the point designated by the white star in Fig. 1(d).The sharp changes in current are reproducible over several sweeps of V SD for a different thermal cycle of the device (see Supplemental Material [43]), with a small variation in threshold voltages due to the stochastic nature of this process.In the double dot configuration, the hysteretic switches define regions (i) to (iv).In region (iii) we observe current values evidencing self-oscillations, while in regions (i) and (v) these self-oscillations seem not to be present.
These current switches can be explained by our model under the assumption that the motion of the CNT mostly affects the electrochemical potential of one of the dots.This assumption is justified by the estimation of capacitive coupling of the dots to the gate electrodes which indicates that one of the dots is primarily controlled by V G3 , whilst the other is mainly controlled by V G5 .The electrochemical potential of the second dot can thus be considered aligned with µ D = 0.In this case, the system is similar to the one-dot case but replacing the transport between the dot and the rightmost electrode by tunneling between the two dots occurring when their electrochemical potentials are aligned, i.e., when x(t) is close to zero; in the model, this interdot exchange is represented by the rate Γ M and a width σ, Γ out (x) ∼ Γ M exp −x 2 /σ 2 (see Appendix D).Notice that, contrary to the case of a single dot, self-oscillations disappear for high bias voltages V SD (region (v)).This phenomenon is difficult to explain with the 1.0 0.5 model outlined above, since V SD simply determines the location of the leftmost Fermi level µ S and should not affect the stability of self-oscillations.However, according to Eq. ( 4), self-oscillations are maintained by a correlation between the charge of the dot n(t) and the instant velocity of the CNT, ẋ(t).This correlation could be lost due, for instance, to inelastic co-tunnelling, which is present for high values of V SD , see Appendix D.
For the double-dot case, which now considers interdot tunnelling and an inelastic contribution to the current proportional to V SD (see Appendix D), we obtain the current shown in Fig. 4(a) (orange line) and the stability diagram depicted in Fig. 4(b).The stability diagram shows good agreement with the measured current switches as a function of V SD in upwards (red) and downwards (blue) sweeps, although the boundaries of region (iii) are not precisely located due to its stochastic nature.We can also observe that the area with ∆E < 0 (light blue) appears close to zero A in region (iii) (inset).This would indicate that the rest position of the CNT is unstable and selfoscillations can be easily induced by fluctuations, a picture which is supported by the measurements of the onset of selfoscillations in region (iii) (see Appendix D).

III. CONCLUSION
In conclusion, we were able to construct stability diagrams that fully characterize the oscillations induced by electron tunneling in single and double quantum dot configurations.We achieve this by observing hysteresis cycles in the current flowing through the device, by using a novel protocol to probe, pump and kill these self-oscillations, and by developing a dynamical model.Our results reveal the subtleties in the coupling between mechanical motion and single electron transport, and open new venues to design autonomous motors and other types of energy transducers at the microscale.

ACKNOWLEDGMENTS
We thank Gerard Milburn for stimulating discussions on the subject of this research.This research was supported by grant number FQXi-IAF19-01 from the Foundational Questions Institute Fund, a donor advised fund of Silicon Valley Community Foundation.NA acknowledges the support from the Royal Society (URF-R1-191150), EPSRC Platform Grant (grant number EP/R029229/1) and from the European Research Council (ERC) under the European Union's Horizon 2020 research and innovation programme (grant agreement number 948932).AA acknowledges the support of the Foundational Questions Institute Fund (grant number FQXi-IAF19-05), the Templeton World Charity Foundation, Inc (grant number TWCF0338) and the ANR Research Collaborative Project "Qu-DICE" (grant number ANR-PRC-CES47).JA acknowledges support from EPSRC (grant number EP/R045577/1) and the Royal Society.JM acknowledges funding from the Vetenskapsrådet, Swedish VR (project number 2018-05061) and from the Knut and Alice Wallen-berg foundation through the fellowship program.JT-B and JMRP acknowledge financial support from Spanish Government through Grant FLUID (PID2020-113455GB-I00) Appendix A: Theoretical model In this section, we detail the theoretical model describing the coupling between the quantum dot and the mechanical oscillator.Similar models have been used in theoretical and experimental works [20,26,30].Since the quantum dot is strongly coupled to the right and left electrodes, the system behaves classically [45].The evolution of the vertical position x of the carbon nanotube (CNT) is described by an underdamped Langevin equation Here, Ω is the resonance frequency of the oscillator, Q the quality factor, m the oscillator mass, g m the coupling constant between the oscillator and n the electron occupation.The term − gm m n represents, precisely, the electrostatic force acting on the CNT due to the electron occupation.ξ represents the environmental thermal noise, following ⟨ξ(0)ξ(t)⟩ = 2γk B T δ(t), with T the environment temperature, k B the Boltzmann constant and γ = mΩ/Q the damping constant.Our experiments are performed at T = 60 mK.
The occupation n = 0, 1 evolves following a dichotomic stochastic process.During each time interval dt, the transition n = 0 → 1 takes place with probability Γ in dt, and n = 1 → 0 with probability Γ out dt. Γ in/out are the energy-dependent tunnel rates defined in the main text.The average occupation n evolves then following the master equation (A2)

Mean-field approximation
For a large quality factor Q and tunnel rates Γ in/out , the time required for the mechanics to thermalize becomes longer than the time required by n to equilibrate.In this case, both fluctuations in occupation and position are negligible, and we consider just the deterministic system of equations which corresponds to a mean-field approximation.In the main text, we wrote equations ( 2) and (3) under this approximation in order to explain the observed self-oscillations.

Mechanical energy
The CNT mechanical energy under the mean-field approximation is Taking the time derivative and using equations (A3) we obtain If x(t), n(t) is a given trajectory of Eqs.(A3), the change in mechanical energy in a certain time τ is Equation ( 4) in the main text uses this expression to evaluate ∆E along a single oscillation period τ = 2π/Ω.

Appendix B: Decay of self-oscillations
The switching between states in the bistability region finds its origin in the thermal and electric noise, see Eq. (A1).In order to clarify this, we expose the bistability in terms of the double-well problem [19,20,46] in the single-dot configuration.We now consider a long timescale T A ≫ 1/Ω where the amplitude of oscillation will vary with time.In this timescale, we consider an effective Langevin equation for mechanical energy, where f (E) = ∆E/τ is the average increment of energy per unit of time, and ∆E is given by Eq. ( 4) of the main text.Note that the value of ∆E depends on the amplitude of the oscillations, and therefore on the mechanical energy, E = mΩ 2 A 2 /2.τ = 2π/Ω is the period of the oscillation.ξ E is a noise containing both electrical and thermal fluctuations, and will depend on the actual value of the energy.Equation (B1) can be written as a Kramers equation considering an effective "potential" given by In Fig. 5, we represent the effective potential U (E) (red line), evaluated using harmonic trajectories x(t) = A cos Ωt.
Here, the problem is explicitly analogous to the typical double-well potential in energy space.In that figure, selfoscillations correspond to the right well, A, with a certain energy associated.The other well, C, corresponds with the motionless state.Switching between both states will happen when a fluctuation puts the system over the barrier B.
We use equation (A1) to simulate one trajectory in the selfoscillation state for 10 4 oscillation periods.The histogram of this trajectory is included in Fig. 5 (blue).This histogram is sharply peaked in the point A, and therefore a very large number of mechanical periods are required for the self-oscillations to vanish.According to this, the self-oscillation duration is given by the shape of the barrier.In the region (II) of figure 2 in the main text, the size of the barrier increases with the bias voltage V SD , generating longer and longer self-oscillations in that region.This increase in duration is observed in figure 2 of the main text along region (IIa).

Appendix C: Fitting of the Coulomb peak and tunneling rates
We fit the Coulomb diamond in the single dot configuration to estimate the tunnelling rates parameters Γ L,R and α L,R .The fitting procedure is similar to our previous work [14] with the addition of energy-dependent tunnelling rates, represented by the parameter α L,R , inspired from Ref. 12.The rates of charges tunnelling in and out between the left (L) or right (R) lead and the quantum dot are given by the expression: where, as in the main text, ϵ is the energy of the dot and f L,R are the Fermi distribution.
The overlap between the density of states of the quantum dot and left/right reservoirs is given by: Finally, the current across the quantum dot is: V SD = 0.9 mV V SD = 0.8 mV We fit the experimental the Coulomb diamonds in Fig. 6(a) by cutting it into multiple Coulomb peak with bias ranging from V SD = 0 mV to V SD = 0.9 mV (Fig. 6(b)).We fit these Coulomb peaks with a unique set of parameters to obtain Γ L = 400 GHz, Γ R = 15 GHz, α R /g m = 0 eV −1 and α L /g m = 400 eV −1 .
To take into consideration the voltage drop at the IV converter internal resistance (100 kΩ), we reduce the bias of the fit by V Corr SD (V G1 ) = V SD − I(V G3 )R s .Because the density of the point of the measurement is not sufficient, we calculate V Corr SD from a first fit of the Coulomb peaks.We then fit again the peaks considering the corrected bias voltage.The contribution of this correction on the final result is small.below a certain threshold, the measurement stops and the sequence restarts to the next column (saving significant measurement time).We ran 200 sequences for each figure (S4(a) and S4(c)), displaying a series of vertical lines from which we can extract the duration of the self-oscillations.We plot a histogram of the duration of the self-oscillations in Fig. S4(b) and S4(d).We note a pattern, particularly visible in Fig. S4(c), where long lived self-oscillations are followed by a rather repetitive number of short lived self-oscillations.In time domain, this correspond to an approximate 8 min period, in which the first 3 min shows a few long lived self-oscillations and the next 5 min a series of short lived self-oscillations.This regular variability of self-oscillations suggests an influence of external sources which could not be identified.

IV. TEMPERATURE DEPENDENCE
In yet another cool down run, we measure self-oscillations in the Coulomb blockade regime under various temperature conditions.In Fig. S5, we show the positive side of a Coulomb diamond, showing a region where self-oscillations are evident.The self-oscillations are first measured at 60 mK (base temperature of the cryostat) under both bias sweep directions (two leftmost panels in Fig S5 ) in order to highlight a bistability region.Then the same Coulomb diamond was measured at increasing cryostat temperatures up to a maximum of 820 mK (other panels in Fig. S5).Note that a 20 mK uncertainty in the device temperature should be considered to take into account the thermalization of the chip.We saw that the self-oscillation area decreases progressively, starting from the regions further away from the Coulomb diamond, as the temperature rises.
The damping of self-oscillations at high temperatures can be explained by the dispersion of the dynamics due to fluctuations.This magnitude can be estimated as σ therm x ∼ k B T /mΩ 2 , where k B is the Boltzmann constant and T is the cryostat temperature.For T = 60 mK, it takes the value σ therm x ≃ 3.7 × 10 −2 nm; since this value is orders of magnitude smaller  than the amplitude of the self-oscillations, its effect on the oscillations is negligible.However, for T = 1 K the dispersion is σ therm x ≃ 0.14 nm.In this regime, the dispersion in the dynamics due to temperature is of the size of the self-oscillations, preventing their appearance.

FIG. 1
FIG. 1. a,b) Schematic of the device in single dot (a) and double dot (b) configurations.A CNT is suspended between the source and drain electrodes.Five gate voltages, VG1-G5, are used to create either a single dot or a double dot.A bias voltage VSD drives a current I through the CNT.The chemical potentials of the source and drain electrodes are µS = eVSD and µD = 0, respectively.The right, left and interdot tunnel rates are indicated and the associated tunnel rates are labeled ΓR, ΓL and ΓM.c) Current measured in single-dot configuration as a function of the gate voltage VG3 and bias voltage VSD.White arrows point at features which indicate the presence of self-oscillations.The current traces in Fig. 2(a) were taken along the dashed vertical line.d) Current measured in the double-dot configuration by sweeping VG5 and stepping VG3 with VSD = 1.8 mV.White arrows point at current features which indicate the presence of selfoscillations.The white star indicates a triple point, where we observe the switch in current as a function of VSD plotted in Fig. 4(a).

6 FIG. 2
FIG. 2. a) Current switch hysteresis as a function of VSD measured in the single dot regime following the dashed line (VG3 = 1.8 V) in Fig. 1(c).The red and blue arrows indicate the direction of each current sweep.The orange dashed line is the current calculated from Eq. (5).The numerals (I, II, III) indicate regions of no oscillations, bistability and self-oscillations, respectively.b) Stability diagram: the dark (light) blue areas indicate ∆E > 0 (∆E < 0) for different values of A and VSD.Small arrows indicate the direction in which A would change in each area given the sign of ∆E.Blue and red arrows delimit regions (I), (II) and (III) and indicate a similar hysteresis cycle as that shown in panel a).We use Ω/2π = 270 MHz, gm = 0.01 eV/nm, m = 2 × 10 −22 kg, and Q = 1000, ΓL = 400 GHz, ΓR = 18 GHz, αL = 4 nm −1 , and αR = 0.

FIG. 3
FIG. 3. a) Sequence of bias voltages VSD applied (blue) and observed current during the protocol for VPROBE = 1.6 mV (IIa), 1.98 mV (IIb), and 2.5 mV (III), selected from panel b).The device is in the single-dot configuration (see Fig. 1(a)).VSD is initially set to VKILL = 0 mV to start the protocol with the CNT at rest.Then VSD is increased up to VPROBE for a given time.self-oscillations are pumped by setting VSD to VPUMP = 3 mV.Finally, the persistence of the self-oscillations is measured by setting VSD back to VPROBE.b) Observed current during the protocol for different values of VPROBE.We identify four regions: (I) absence of self-oscillations; (IIa) self-oscillations observed after the pumping step and spontaneously decaying after a random time; (IIb) self-oscillations spontaneously appearing at a bias potential VSD = VPROBE; and (III) stable self-oscillations.

FIG. 4 .
FIG. 4. a) Current switch hysteresis as function of VSD in the doubledot configuration, measured at the white star location in Fig. 1(d) (VG3 = 620 mV, VG5 = −120 mV).Red and blue arrows indicate the direction of each current sweep.The dashed orange line shows the current calculated from Eq. (5) with ΓM = Ω, σ = 0.1 nm and amplitude 1 nm (see Dfor details).The numerals (i, ii, iii, iv, v) indicate different regions in the stability diagram.b) Stability diagram: the dark (light) blue areas indicate ∆E > 0 (∆E < 0) for different values of A and VSD.Small arrows indicate the direction in which A would change in each area given the sign of ∆E.Blue and red dashed lines delimit regions (ii), (iii) and (iv).Blue and red arrows indicate a similar hysteresis cycle as that shown in panel a).Inset: zoom in on areas approaching zero amplitude.We use Ω/2π = 270 MHz, gm = 0.01 eV/nm, m = 2 × 10 −22 kg, and Q = 1000.

FIG. 5 .
FIG. 5. Numerical computation of the effective potential U (E) (red) and stochastic simulation of a single trajectory (blue) as a function of the mechanical energy for the one dot-configuration.The sharp peak in U at the point A is due to the end of the borders of the conduction region of the device.The letters A, B, C indicate the stable oscillation points (A, C) and the tipping point between both regimes (B).The parameters are the same as in the main text, Ω/2π = 270 MHz, gm = 0.01 eV/nm, m = 2×10 −22 kg, Q = 1000, ΓL = 400 GHz, ΓR = 18 GHz, αL = 4 nm −1 , and αR = 0.

FIG. 6 .
FIG. 6. Coulomb diamond and fit.a) Zoom on the Coulomb diamond from Fig. 1d of the main text.The white dashed rectangle represent the area from which are obtained the cut of panel b).b) Set of cuts of panel a) at different bias (circles) and fitting curves (lines).For each bias, the cuts and fit are vertically offset by 0.2 nA from the previous bias.The fitting parameters of the all set are: ΓL = 400 GHz, ΓR= 15 GHz, αL/g = 0 eV −1 , αR/g = 400 eV −1 and V0 = 1.826V.The lever arm 0.18 eV/V is obtained from panel a). of

FIG. 7
FIG. 7. a) Current measured in the double-dot configuration by repeating the sequence shown Fig.3 for different VPROBE.First, any selfoscillations are killed setting VSD to VKILL = 0 mV, then VSD is set to the target voltage VPROBE.About 1.5 s later, VSD is set to VPUMP = 1.8 mV to start self-oscillations for about 0.5 s.Finally, the voltage is set back to the target VSD.The different behaviours are classified in five areas: (i) absence of self-oscillations, (ii) self-oscillations observed after the pumping step, (iii) self-oscillations spontaneously appearing, (iv) self-oscillations observed after the pumping step and spontaneously decaying after a random time, (v) absence of self-oscillations.b) IV measurement with the two sweep direction in the configuration of panel b) showing the current switch hysteresis; same as Fig. 4(a) of the main text.
FIG. S4. a,b) 200 measurements in a row of the duration of self-oscillations for two different values of VSD.For each vertical line, the selfoscillations are first excited by setting the bias voltage to VSD = −4 mV for about 5 s then the bias is set to the target voltage VSD = −3 mV (a) or VSD = −3.2mV (c).Then we measure the current at a rate of about one point per second.We detect the stop of the self-oscillation when the current falls below a pre-determined threshold.The top of each row indicates the duration of the self-oscillations.b,d) Histograms of the duration of self-oscillation for (a) and (c), respectively.
FIG. S5.Measurements of self-oscillations in a Coulomb diamond for different temperatures.The leftmost panel was measured at base temperature, 60 mK, with two different voltage bias sweep directions indicated by white arrows.They reveal an area of self-oscillations and another of bistability showing self-oscillations in one of the sweep directions.In the other panels, the Coulomb diamond is always measured with the same bias sweep direction but with increasing temperatures, up to 820 mK.The self-oscillation area progressively vanishes as the temperature increases.