Energy Transfer into Period-Tripled States in Coupled Electromechanical Modes at Internal Resonance

Efficient energy transfer often occurs between oscillation modes in a resonator when they are tuned to internal resonance. We design the eigenfrequencies of two vibrational modes of an electromechanical resonator to be close to a ratio of 3:1 and demonstrate that the energy supplied to the upper mode can be controllably transferred to the lower mode. With the lower mode vibrating with a period tripled that of the upper mode, the discrete time-translation symmetry imposed by the periodic drive is broken. The lower mode settles into one of three stable period-tripled states with different phases. This channel for energy transfer from the upper mode can be turned on or off without changing system parameters. When the upper mode itself becomes multistable under strong resonant or parametric drive, additional sets of coexisting period-tripled states emerge in the lower mode. In the latter case, we measure a total of 6 coexisting vibration states with identical amplitude but phases differing by $\pi$/3. Excitation of coexisting states with three different phases could open new opportunities in designing mechanical memory based on ternary logic. Coupled resonators with period-tripled states can also be used to model complex interacting systems with spin equals one.


I. INTRODUCTION
Nonlinear mode coupling enables energy exchange between oscillatory modes with significantly different resonant frequencies. When the resonant frequencies are incommensurate or far from each other, application of parametric modulation is often necessary to enhance the coupling. The enhanced interaction allowed mechanical resonators to be controlled through their coupling with optical, microwave or phononic cavities [1,2]. When the frequency of one mode becomes commensurate with an integer multiple of another mode, nonlinear coupling of the modes becomes strong, giving rise to efficient energy transfer and internal resonance, which profoundly affects the response of both modes. Internal resonances in macroscopic systems with mode frequencies in the ratios of 3:1 or 2:1 exhibit a wealth of phenomena including nonlinear oscillations, amplitude-modulated oscillations, multi-stable behavior, dissipative nonlinearities and chaotic motions [3][4][5][6][7][8]. Recently, internal resonance has been demonstrated to play an important role in micro-and nano-mechanical resonators [9][10][11][12][13][14][15][16]. When the lower mode is excited by a periodic drive, energy exchange with the upper mode leads to novel behaviors in the steady state response as well as the decay of vibrations when the energy supply is removed. For example, the frequency fluctuations in the self-sustained oscillations of the lower mode of a micromechanical device is dramatically reduced at internal resonance with a higher mode [9]. In resonators made of nanotube or two-dimensional materials, complex lineshapes and multiple hysteresis are observed [10,11].
The energy exchange between the two modes at internal resonance also leads to anomalous decay of vibrations with rates that changes with energy in both micromechanical structures [15] and suspended graphene [16].
While much progress has been made by pumping energy into the lower mode, the regime in which energy is supplied to the upper mode for systems at internal resonance is much less explored. When the nonlinear coupling opens up a channel of energy transfer from the upper mode to the lower mode, stable vibrations are induced in the latter with a period that is an integer multiple of the period of the drive. As a result, the discrete time-translation symmetry [17,18] imposed by the periodic drive is broken. The breaking of discrete time translation symmetry is well-known in single mode parametric oscillators [19][20][21][22][23][24][25][26][27][28][29], where vibrations are period-doubled relative to the modulation. The two coexisting stable oscillation states with identical amplitude but opposite phase can be used to represent a classical bit [30] or an Ising spin [31][32][33][34][35] for logic or computation applications. For parametrically-driven quantum oscillators, the breaking of time-translation symmetry imposed by the periodic drive, together with many-body interactions and disorder, can lead to the formation of time crystals [36][37][38]. While period-doubled oscillations in a single mode rely on modulations of the parabolic term in the confinement potential, it has been predicted [39][40][41] that if the potential contains a cubic term that is modulated at a frequency ω f close to 3 times the resonant frequency ω o , period-tripled oscillations can be excited to break the discrete timetranslation symmetry imposed by the modulation. The simplest model for a single-mode system with stable period-tripled oscillations is given by the equation of motion: where q is the displacement of the mode, m is the effective mass, γ denotes the Duffing nonlinearity and F p is the strength of the parametric pump at frequency ω f . The system settles into one of three stable vibrations states with phases differing by 2π/3 or the zero- Apart from direct parametric modulation of the mode through a term of q 2 cos(3ω 1 t), period-tripled oscillations can also be excited by supplying energy to a higher mode that is at 3:1 internal resonance with the first mode with resonant frequency ω 1 . Through the interaction energy γ 12 q 3 1 q 2 (q 1,2 represent the displacements of modes 1 and 2 respectively), oscillations in mode 2 at amplitude A 2 at its eigenfrequency ω 2 generates an effective parametric modulation of the form 3γ 12 q 2 1 [A 2 cos (3ω 1 t) ] in the equation of motion of mode 1. Internal resonance facilitates the transfer of energy from mode 2 to mode 1, opening the possibility to excite oscillations in the latter with a period that is three times the periodic drive.
Whether stable period-tripled states can be excited in mode 1 depends on the dynamics of both modes because oscillations in mode 1 also affect mode 2. Enhancement of the drive by a second mode to generate period-tripled oscillations was demonstrated in superconducting resonators [42]. The signal includes contributions from the period-tripled states and the zero-amplitude state as the system switches among them due to noise. Maintaining the system in a period-tripled state with practically no interstate switching may be beneficial for a number of applications such as memory and frequency conversion. To our knowledge,  Here, we show that a properly designed electromechanical resonator with two vibrational modes at frequencies close to 3:1 internal resonance exhibits stable period-tripled vibrations in the lower mode (mode 1) that breaks the discrete time-translation symmetry of a periodic drive applied to the higher mode (mode 2) near its eigenfrequency. Since the zero-amplitude state of mode 1 remains stable for all drive frequencies and amplitudes, energy transfer from mode 2 to create the period-tripled states requires activation. Mode 1 needs to be perturbed by a secondary drive, which is subsequently removed. By choosing the phase of the perturbation, mode 1 can be controlled to settle into one of the three period-tripled states that differ in phase by 2π/3 when the periodic drive on mode 2 is above a certain threshold value, but remains relatively weak. At certain system parameters, a set of 3 coexisting limit cycles emerges in the rotating frame. When the drive amplitude on mode 2 is increased beyond another threshold, a second set of period-tripled vibrations at a different amplitude can be excited for mode 1, yielding a total of seven coexisting states for a range of drive frequencies. If the resonant drive on mode 2 is replaced by parametric modulation of the confinement potential near twice its eigenfrequency, we observe six coexisting vibrational states in mode 1 with the same amplitude but with phases differing by integer multiples of π/3. These states are period-sextupled relative to the parametric modulation on mode 2.
Apart from fundamental interests relevant to the breaking of discrete time symmetry and classical time crystals [43], the ability to prepare the system in a specific period-tripled state opens new opportunities in using mechanical resonators in ternary memory and logic [44], as well as in simulations of spin systems with spin of 1 [45] in a context similar to Ising machines for spin of 1/2 [31][32][33][34][35].
We emphasize that the stable period-three vibrations described in this paper are qualitatively different from the multiple-period oscillations often observed in strongly nonlinear systems where the nonlinear part of the energy is of the same order of magnitude as the harmonic part. The latter are best known in the context of the period-multiplication cascades of limit cycles leading to dynamical chaos [46][47][48]. They are generally profoundly nonsinusoidal, with multiple Fourier components. In contrast, the period-tripling phenomena discussed in this paper rely on much weaker nonlinearity. The magnitude of the nonlinear terms in the potential energy of the modes is a factor of ∼10 −7 smaller than the magnitude of the harmonic terms. Vibrations in both modes are near perfectly sinusoidal. As we show, the very occurrence of the onset of these vibrations is a consequence of the interplay of the internal resonance between the modes and the weak damping. Therefore a comparatively weak drive is sufficient to excite them. The results presented in this paper are measured well below the threshold at which chaos occur. There are two electrodes underneath each plate, as shown in the schematic in Fig. 1(d). By adjusting the dc voltages V dc1 and V dc2 applied to these electrodes, ω 1,2 can be fine-tuned for the system to go into and out of internal resonance (see Appendix A). Previous studies of internal resonance in micro-and nano-mechanical systems were often hampered by the difficulty to independently excite and simultaneously detect the vibrations of the two modes.
Our design circumvents such problems. For each plate, an ac voltage applied to the electrode on the left in Fig. 1 We first characterize the system by applying ac voltages to the electrodes to periodically drive the two modes. Both modes have softening nonlinearity, with Duffing constants of −8.1 × 10 −7 kg m 2 s −2 and −6.1 × 10 −8 kg m 2 s −2 respectively. There is strong dispersive coupling between the two modes, with coupling energy 1 2 γθ 2 1 θ 2 2 , where θ i is the rotation angle of the plate i for mode i (i = 1, 2). As shown in Figs. 1(e) and 1(f), the resonant frequency of one mode decreases by an amount proportional to the square vibration amplitude of the other mode. The linear fits yield the constant γ = −1.02 × 10 −7 kg m 2 s −2 . While dispersive coupling is not required to generate period-tripled oscillations, we find that it must be included to obtain agreement between our measurement and theory. For the rest of the paper, energy is only supplied to mode 2. Unless otherwise stated, the periodic drive on mode 1 is turned off.
Another relevant term in the interaction energy is of the form γq 3 1 q 2 , which becomes important for small mismatch 1 = ω 2 /3 -ω 1 when the two modes are near internal resonance.
The equations of motions are given by: where I 1,2 are the effective moment of inertia of the two modes and τ d is the amplitude of the periodic torque applied to mode 2. The last and second last terms on the left side of Eqs. (2) originate from the dispersive coupling and the interaction energy γq 3 1 q 2 respectively. With the driving frequency ω d on mode 2 close to ω 2 , we change from θ(t) andθ(t) to complex amplitudes θ 1 (t) = u 1 (t) exp [i (ω d /3) t] + c.c.., θ 2 (t) = u 2 (t) exp [iω d t] + c.c.. Under the rotating wave approximation in which the fast-oscillating terms are dropped, u 1,2 (t) evolves according to: where = ω d −ω 2 and 1 = ω 2 /3−ω 1 . Parameters τ 2 , γ 11 , γ 22 , γ 12,21 and γ 12,21 are determined by τ d , γ 1 , γ 2 , γ and γ in Eq. (2) respectively, taking also into account contributions from the renormalization due to the coupling between the two modes. We note that the rotating wave approximation is justified as the nonlinear terms in the potential energy are much smaller (∼ 10 −7 ) than the harmonic terms.

III. ENERGY TRANSFER INTO PERIOD-TRIPLED STATES
We focus on the steady state oscillations, withu 1 =u 2 = 0. The zero-amplitude state of mode 1 is found to be stable for all τ 2 , and 1 . Measurements on the device yield results that are consistent with this theoretical finding: if the vibration amplitude of mode 1 is initially zero, finite-amplitude vibrations of this mode cannot be excited when ω d is swept up or down. In this scenario, the response of mode 2 is identical to that of a single mode Duffing oscillator, as illustrated by the black data in Fig. 2(a) for pump voltage of 6 mV.
V dc1 and V dc2 are chosen to be 1.475 V and 1.975 V respectively for the system to be near internal resonance, with 1 /2π = -16.4 Hz. To achieve stable period-tripled oscillations, it is necessary to apply a proper perturbation to the system to excite mode 1 out of the zeroamplitude state. We choose to perturb the system by applying a secondary ac drive (with amplitude 23 mV) to mode 1 to produce a torque of the form τ 1 cos( ω d 3 t + φ). τ 1 and φ can be adjusted to excite vibrations in mode 1 with different amplitude and phase, as shown in Fig. 2 Upon removal of the secondary drive, the system either settles into one of the three period-tripled states, represented by the red, blue and green squares, or the zero-amplitude state represented by the black square, depending on τ 1 and φ. Until the perturbation is applied again, mode 1 remains in one of these 4 stable states.
We emphasize that Fig. 2(b) plots the behavior of mode 1 only. In analyzing the dynamics of the system, mode 2 also needs to be included. Hence the phase space is four dimensional (4D) in the rotating frame. To settle into a specific state, the system must start inside the basin of attraction of the corresponding state in 4D phase space. Trajectories from different initial points do not intersect with each other in the 4D space. Figure 2  (shown in red) is disconnected from the zero-amplitude branch. The corresponding isolated branch of mode 2 is shown in red in Fig. 2(a). When is increased beyond ω H , mode 1 jumps to the zero-amplitude state. Correspondingly, mode 2 jumps to the response of an ordinary Duffing oscillator. The system also undergoes a jump if is decreased beyond ω L .
Once a jump to the zero-amplitude state takes place, getting back to any one of the three period-tripled states in mode 1 requires activation, as discussed earlier. The amplitude for the unstable oscillation states (the saddle points) is calculated from Eq. (3) and is shown as dashed lines (Supplementary Information). We find that a saddle-node bifurcation occurs at ω H , where one stable state merges with one saddle point and disappear. At ω L , the system undergoes a Hopf bifurcation where limit cycles develop in the rotating frame, as we will discuss later. In mode 1, the branches of stable and unstable states form an isolated loop.
For mode 2, the stable and unstable branches attain peak amplitude at /2π = -5.5 Hz and -6.8 Hz respectively. If the amplitude of the pump is reduced, the loops decrease in size and eventually shrink to zero at pump amplitude of 5.5 mV. Below this threshold value of pump amplitude, only the zero-amplitude state is stable and the period-tripled states cannot be excited even with activation.
The energy transfer from mode 2 to mode 1 can be controllably turned on and off. 4.5 a result, a new, isolated branch of stable vibrations emerges at lower amplitude, plotted in light red in Fig. 3(a). The isolated branch of unstable states, calculated using Eq. (3), also becomes wider and higher. It is connected to the stable branch as shown in the inset of Fig. 3(a). For the upper branch of stable vibrations, there are saddle-node bifurcations at both ends (ω L1 and ω H1 ). Starting from /2π = -5.1 Hz [red square in Fig 3(a)], if the drive frequency is swept up beyond ω H1 , as indicated by the thin dotted arrows pointing to the right, mode 2 jumps up to the Duffing response of a single mode (black). If, instead, the driving frequency is swept down, mode 2 jumps down to the Duffing response at ω L1 .For the lower branch of stable vibrations, there is a saddle-node bifurcation on the right end at ω H2 . The left end at ω L2 is a Hopf bifurcation where limit cycles develop in the rotating frame, as we will discuss later. Starting from a stable vibration state in the lower branch, sweeping the drive frequency up or down eventually lead to jumping of mode 2 down to the Duffing response, as indicated by the two light red arrows. Once mode 2 jumps to the Duffing response from either the upper or lower branch, merely reversing the direction of frequency sweeping does not bring the system back to the isolated branch. Returning to the isolated branch requires activation, similar to the case of a single set of period-tripled states discussed earlier.
Each stable state in the isolated loop of mode 2 is associated with one set of three stable period-tripled oscillation states of mode 1. Figure 3(b) shows the amplitude of the periodtripled states in mode 1. The zero-amplitude state, despite lying out of the range of the vertical axis of the plot, remains stable. For driving frequencies between ω L2 and ω H2 there are a total of 7 stable states in mode 1. The inset of Fig. 3(b) shows the seven coexisting states in phase space for /2π = -9.72 Hz. We note that no switching between the states is observed. Thermal noise is not strong enough to induce transitions in the duration of our measurement. To place mode 1 in the low-amplitude branch, one method is to start from a smaller drive amplitude on mode 2 that yields only one single branch of period-tripled states, as in Fig. 2. When the drive voltage is gradually increased, the system settles to the low-amplitude branch in mode 1. Placing mode 1 in the high-amplitude branch requires slight modifications to the scheme because the frequency range for stable period-tripled states (between ω L2 and ω H2 ) does not overlap with that at smaller drive amplitude of Fig.   2 (between ω L and ω H ). The drive amplitude is increased in small steps. Between successive steps, the drive detune frequency is decreased by an appropriate amount so that mode 1 settles in the high-amplitude branch at the end of the process. More details of preparing mode 1 in the high-amplitude state is described in Appendix D.
Coexisting stable states of vibrations of mode 2 can also be excited by parametric modulating its resonant frequency. We modulate the gradient dτ /dθ of the electrostatic torque exerted by the electrode on the small plate at frequency ω p near 2ω 2 (see Appendix A). When there is no energy transfer to mode 1, mode 2 follows the response of a single-mode parametric resonator, as shown in black in Fig. 3(c). The oscillation states are period-doubled with respect to the modulation. They occurs in pairs that are out of phase with each other.
For certain range of parameters, internal resonance enables the energy transfer into mode 1, leading to stable vibrations in mode 1 with period 12π ωp . Figure 3(d) plots the dependence of the amplitude of mode 1 on the modulation frequency detuning 2 = ω p /2 -ω 2 . We measure six coexisting stable vibration states of mode 1, as shown in the inset of Fig. 3(d) where The six stable states are measured to have phase offset of π/3 from one another. They are period-tripled compared to oscillations in mode 2, and period-sextupled compared to the parametric modulation on mode 2. Unlike Fig. 2(b) and the inset of Fig. 3(b), the inset of Fig. 3(d) displays six-fold symmetry. When energy is transferred to mode 1, mode 2 follows the red response in Fig. 3(c) that is slightly shifted from the black branch. In calculating the vibrations amplitudes and phase, Equations (2) and (3) need to be modified (See Appendix F).

V. COEXISTING LIMIT CYCLES
Interesting phenomena such as limit cycles and chaos have been demonstrated in coupled micro-and nanomechanical resonators [49][50][51][52][53]. In particular, limit cycles are observed for two coupled modes in 3:1 internal resonance when energy is supplied to the lower mode [52].
Our system, in contrast, is under periodic drive on the higher mode and shows qualitatively different behavior. Figure 4(a) plots the amplitude of mode 1 when mode 2 is driven at pump amplitude of 9 mV. It is a zoom-in of the high-amplitude branch of Fig. 3(b). At /2π = -9.93 Hz, a Hopf bifurcation takes place and limit cycles emerge. Each of the three period-tripled states evolves into a limit cycle in the rotating frame, with a three-fold symmetry illustrated in Fig. 4(b). As is further decreased, the limit cycles become larger and eventually disappears at /2π -10.07 Hz. Without energy transfer from mode 2, the vibration amplitude of mode 1 drops to zero. The grey band in Fig. 4(a) indicates the frequency range for the occurrence of limit cycles both in measurement and in calculations.
Limit cycles only appear in a narrow frequency range of ∼ 0.14 Hz that is small compared to the frequency range at which stable period-tripled vibrations take place in mode 1 (∼ 12 Hz). For smaller drive amplitudes in Fig. 2(c), limit cycles also develop at the left side of the isolated branch of period-tripled vibrations. However, the frequency range at which limit cycles occur is too small to be visible.

VI. DISCUSSIONS AND OUTLOOK
Our experiment shows that by using an electromechanical resonator with two modes at internal resonance, the nonlinear coupling leads to an effective parametric modulation on the lower mode that is resonantly enhanced by the higher mode. Upon activation, energy can be transferred from the higher mode to the lower mode to excite oscillations in the latter with a period that is tripled the period of pumping. The system can be prepared in any one of the three period-tripled states of mode 1 and remains stable in the duration of the experiment, opening the possibility of using parametrically driven mechanical modes for ternary memory and logic. More efficient schemes for writing bits in period-doubled states in "parametrons" [20,22,29,54] can be modified to be used in the period-tripled system.
Demonstration of such schemes in our two-mode system is in progress.
Away from the internal resonance with the lower mode, the higher mode has multiple stable states with different amplitude and/or phase due to its own nonlinearity. Each of these states can enhance the effective parametric modulation on the lower mode when the two modes are tuned into internal resonance. If the coupling is sufficiently strong, each state will be accompanied by a set of period-tripled states in the lower mode. We demonstrate the excitation of two sets of period-tripled states when the higher mode is under a strong resonant drive. Alternatively, we apply parametric modulation to the higher mode, resulting in two sets of states in mode 1 that are period-tripled with respect to vibrations in mode 2 and period-sextupled with respect to the parametric modulation applied. Mode 2 can also be designed and excited to have even larger number of coexisting stable states. For example, if the Duffing (third order) and fifth order nonlinear coefficients in the restoring force are of opposite signs, the higher mode can possess three stable states at different amplitude and phase [55]. If each of these states generates a set of period-tripled states in the lower mode, there will be a total of ten coexisting states in the latter including the zero-amplitude state. A single resonator can therefore, in principle, be used to represent a decimal digit.
Alternatively, if mode 2 is subjected to parametric modulation described by Eq. (1), three coexisting states can be excited. Internal resonance with the lower mode can potentially lead to nine states in the lower mode with finite vibration amplitude but phase shifted by 2π/9. With a period that is 9 times the parametric modulation, the reduction of discrete time-translation symmetry will be even more pronounced than the data presented in this paper.
Near the bifurcation points where stable period-tripled states disappear, effects of fluctuations are expected to be strong. For a single mode system experiencing period-tripling described by Eq. (1), the basins of attraction of the period-tripled states are isolated from each other in the 2D phase space. As shown in Fig. 1(b), each basin is surrounded by the basin of the zero-amplitude state. In the presence of thermal or external noise, the system fluctuates about the stable state. If the system crosses the separatrix, it switches to the zero-amplitude state. For weak noise, the system does not switch to the other two periodtripled states. Interestingly, it has been predicted that the behavior is qualitatively different in quantum oscillators, where direct switching between the period-tripled states can take place [56]. In the two-mode classical system of this paper, the phase space is 4D, making it difficult to plot and visualize the basin of attraction. By injecting noise into the drive to increase the effective temperature, the system can be induced to switch out of a stable period-tripled state. Preliminary results indicate that when the system switches out of one of the period-tripled states, it settles in the zero-amplitude state for the system parameters chosen in the current experiment. The results are consistent with the notion that the basins of the period-tripled states are isolated from each other in the 4D phase space, similar to the much simpler case of the single mode system described by Eq. (1). Whether the behavior changes for other system parameter values, especially those that exhibit limit cycles in phase space, remains an open question that warrants further studies.
The switching from the isolated branch of period-tripled states to zero amplitude in both directions can be used for detecting small perturbations of the eigenfrequencies of the modes. At the saddle-node bifurcations discussed in this paper, the activation barrier for switching is expected to scale with the frequency detuning from the bifurcation point with an exponent of 3/2 [57][58][59]. Unlike traditional bifurcation amplifiers [60] that are only sensitive to unipolar change in parameters, here perturbations in both directions lead to jumps in amplitude that can be easily detected. High sensitivity can be achieved by tuning the amplitude of the pump on the mode 2 to produce a narrow frequency range at which the isolated period-tripled branch in the lower mode is stable.
With the capability to excite and detect vibrations of the two modes independently, our system provides a versatile platform for studying nontrivial nonlinear effects including limit cycles and dynamical chaos that are absent from single mode systems with only two dynamic variables described by Eq.
(1). We demonstrated that a set of 3 coexisting limit cycles develops over a narrow range of driving frequencies when mode 2 is under resonant drive. The limit cycles in the rotating frame correspond to non-sinusoidal vibrations that could lead to the development of frequency combs [61][62][63]. Just before the limit cycles disappear, fluctuations are expected to become strong. Future efforts will explore whether period doubling cascades and/or chaotic motion arise when the driving amplitude on mode 2 is further increased.
While the period-tripled vibrations reported here are observed in two vibrational modes in a micromechanical resonator designed to be near internal resonance, the findings are generic and applicable to other resonators, such as coupled superconducting coplanar waveguides [42], optical cavities [64], membranes of two-dimensional materials [11,16] and carbon nanotubes [10]. The main requirements are that the ratio of the frequencies of the two nonlinearly coupled modes are near 3:1 and the higher mode is subjected to a resonant drive sufficiently strong but well-below the threshold of driving the system into dynamical chaos. In fact, the coexistence of period-tripled oscillation states and the zero-amplitude state was demonstrated in superconducting coplanar waveguide resonators [42]. However, the detected signal includes contributions from all four states as the system switches among them due to noise.
The interplay of fluctuations and dynamics has novel, and yet fairly general features in systems that display period tripling. Some of these features arise because the zero-amplitude state remains dynamically stable independent of the drive frequency and strength. Activation is required to excite the system into one of the three period-tripled states. The activation can come from thermal noise [65] or, as recently suggested [66], from quantum fluctuations.
A nontrivial aspect of the activation is that, as the driving amplitude is increased, the unstable states in the phase space [hollow squares in Fig. 1(b)] move toward the origin (the zero-amplitude state). The activation barrier for switching out of the zero-amplitude state approaches zero, but never disappears completely. It has been predicted that in this limit the switching dynamics is qualitatively different from the escape from metastable states near bifurcation points: in contrast to conventional wisdom, there is no detailed balance and the switching is not controlled by a soft mode. Even though the theoretical results [65,66] refer to the case of a single mode, the analysis can likely be extended to the coupled modes. It is possible that such generic behavior in switching could be revealed in resonator systems with modes at internal resonance [10,11,16,42,64], including the one measured in this paper.
Another important feature of period tripling, which makes it qualitatively different from the seemingly similar and extensively studied period doubling, is the nontrivial topology.
This topology has been found in the analysis of quantum tunneling between period-three states of a single mode and were related to the geometric phase between the period-three states [39,67,68]. We demonstrated in a coupled-

TION OF VIBRATIONS
For each plate, a voltage V d applied to the left electrode in Fig. 1(d) generates an electrostatic torque: where C is the capacitance between the movable plate and the electrode. The torque can be written as a Taylor expansion about the equilibrium angle θ 0 : where C , C , C (3) and C (4) represent the first, second, third and fourth derivative of C with respect to θ respectively. Higher order terms are neglected. The voltage V d consists of a sum of a dc component V dc and an ac component with amplitude V ac and frequency ω f : The time-independent component of τ given by Eqs. (A-2) and (A-3) changes the system parameters. First, the constant term produces a shift in the equilibrium position of θ to θ 0 . Second, electrostatic spring softening is induced by the term proportional to θ − θ 0 . The restoring torque on the plate is modified by an amount 1 2 C (θ 0 )V 2 dc (θ − θ 0 ). Our experiment relies on such electrostatic tuning of the resonant frequencies to control the frequency mismatch 1 of the two modes by adjusting V dc1 and V dc2 . Third, contributions from the nonlinear restoring torques of 1 4 C (3) (θ 0 )V 2 dc (θ − θ 0 ) 2 and 1 12 C (4) (θ 0 )V 2 dc (θ − θ 0 ) 3 dominate the Duffing nonlinearity.
The ac voltage leads to both an additive torque C (θ 0 )V dc V ac cos(ω f t) and a torque C (θ 0 )V dc V ac cos(ω f t)(θ − θ 0 ) that modulates the effective spring constant. In the main text, when the drive is applied to mode 2 at a frequency is close to ω 2 , ω f is replaced by ω d to denote resonant driving. When the frequency is close to 2ω 2 , the main effect is parametric modulation of the spring constant. ω f is replaced by ω p to distinguish the parametric modulation from the case of resonant driving.
Vibrations of each plate are detected by capacitive detection. The application of dc voltages to the right electrodes in Fig. 1(d) leads to build-up of charges between the top plates and the electrodes. Rotations of the plates lead to changes in their capacitance.
Charge flowing out of the two plates are measured simultaneously and independently by two separate charge sensitive amplifiers, the outputs of which are fed into the two channels of a lock-in amplifier (Zurich Instruments HF2LI).

POINTS
The stable states of vibration and the saddle points are calculated by assuming that the complex amplitudes u 1 and u 2 are independent of time. Equations (3) reduce to two algebraic equations, the solutions of which give the stationary states: where u 1,2 = |u 1,2 |e iϑ 1,2 . The stability of the states is found by linearizing Eq.  The eigenfrequency ω 1 and damping coefficient Γ 1 of mode 1 are determined by fitting to the response when a small resonant drive is applied to mode 1 only but not to mode 2.
Subsequently, the driving amplitude is increased so that the response becomes nonlinear.
Dispersive coupling energy of the form 1 2 γθ 2 1 θ 2 2 leads to a shift in the resonant frequency of one mode by an amount proportional to the square vibration amplitude of the other mode. In Figs. 1(e) and 1(f), the slopes of the linear fits are given by γ 12 4ω 1 and γ 21 4ω 2 respectively, yielding γ 12 = −1.92 × 10 11 s −2 and γ 21 ≈ −1.34 × 10 12 s −2 . Measurement of γ 12 and γ 21 also allows the ratio of the effective moment of inertia of the two modes to be determined, through the relation I 2 I 1 = γ 12 γ 21 . In mode 2, the rotation of the large plate is negligible. By taking I 2 to be the moment of inertia of the small plate about the rotation axis (7.61 × 10 −20 kg m 2 ), I 1 is calculated to be 5.34 × 10 −19 kg m 2 .
All the system parameters in Eq. (3) are thus individually characterized. The only exceptions are γ 12 and γ 21 . Both originate from the coupling energy γq 3 1 q 2 . Their ratio is known, given by γ 12 γ 21 = I 2 I 1 , similar to the case for dispersion coupling discussed above. By settinġ u 1 = 0 in the first equation in Eq. (3), the steady state vibration amplitudes satisfy: γ 12 can be written in terms of the other system parameters and u 1,2 : respectively. As the drive amplitude is increased, the frequency range for stable vibrations widens. The peaks become higher and tilt towards low frequencies. At the largest driving amplitude of 10 mV, the peak tips over and a second, isolated stable branch appears at low amplitudes. Figure 6(b) plots the amplitudes for the corresponding vibrations of mode 1 with period tripled that of mode 2.
Starting from the red curve for a driving amplitude of 6 mV in Fig. 6(a), increasing the drive amplitude while keeping frequency detune fixed brings mode 2 to the upper branch for driving amplitude of 10 mV. To place mode 2 in the lower branch, it is necessary to decrease the driving frequency as the driving amplitude is increased in small steps. One choice for the driving frequencies is illustrated by the black squares in Fig. 6(a). The frequency mismatch 1 = ω 2 /3 − ω 1 is fixed for all results presented in the main text, with 1 /2π = −16.4 Hz. We chose 1 to be non-zero because excitation of period-tripled vibrations in our system requires both modes to vibrate at moderately large amplitudes so that the shifts in eigenfrequencies due to the Duffing nonlinearities are significant. While the Duffing coefficients of both modes are negative, the shift in the frequency of mode 1 is much larger (∼ 5 times) than mode 2. It is therefore necessary to pick ω 1 exceeding ω 2 /3 so that at the larger amplitudes, the ratio of the shifted eigenfrequencies becomes close to 3 to satisfy the conditions of internal resonance. For 1 > 0, we do not observe any stable period-tripled vibrations in mode 1 in our system.
The range of 1 that yields stable period-tripled vibrations can be calculated using the procedure described in Appendix B. mV. When the pump amplitude is increased to 9 mV, it is possible to excite a second branch of period-tripled states, as illustrated by the light red region in Fig.7 When the resonant drive on mode 2 in Eq. (2) is replaced by parametric modulation of its eigenfrequency, the equation of motions are modified to:   θ 1 + 2Γ 1θ1 + ω 2 1 θ 1 + γ 1 I 1 θ 3 1 + 3 γ I 1 θ 2 1 θ 2 + γ I 1 θ 1 θ 2 2 = 0 where k e = −C V dc2 V ac2 is the amplitude of the modulation of the torsional spring constant, C is the second derivative of the capacitance between plate 2 and the driving electrode with respect to θ 2 and V ac2 (80 mV) is the amplitude of the periodic ac voltage applied to the driving electrode. With the modulation frequency ω p close to 2ω 2 , we change from θ(t) andθ(t) to complex amplitudes θ 1 (t) = u 1 (t) exp [i (ω p /6) t] + c.c.., θ 2 (t) = u 2 (t) exp [i (ω p /2) t] + c.c..
Under the rotating wave approximation, u 1,2 (t) evolves according to: The stationary states and their stability is determined using the procedure described in Appendix B.  Fig. 3(a)] is omitted for clarity. Figure 8(  The vibrational modes of our device are modelled using finite element analysis (COMSOL). Figure   S. 1 shows the results for the first four modes. The second and fourth modes are torsional modes. They correspond to modes 1 and 2 in the main text respectively. Their eigenfrequencies are designed to be close to a ratio of 1:3 so that they can be tuned into internal resonance by applying voltages to the electrodes.
The other two modes can also be observed in measurement. Figures S. 1(b) and (c)  When the periodic drives to mode 1 and mode 2 are removed, the amplitude of vibrations decays.  Due to the Duffing nonlinearity, the eigenfrequency of each mode increases as the vibration amplitude decreases. In Fig. S. 2, the frequency mismatch 1 = 2 1 is chosen to ensure that the ratio of the eigenfrequencies never gets close to 3 for internal resonance and there is no energy exchange between two modes. For other choices of 1 that allow the system to enter into the regime of internal resonance during the ringdown, the interaction energy 1 3 2 can give rise to decay coefficients dependent on the vibration amplitude. Such nonlinear damping has been observed by a number of groups [s1, s2].  contains contributions from both modes. The strains for beams 2 and 3 can also be calculated using a similar procedure.
The elastic potential energy is then calculated with where cij are the second order stiffness tensor. Contributions due to material nonlinearities are assumed to be negligible.