Generation and propagation of topological solitons in a chain of coupled parametric-micromechanical-resonator arrays
11 Generation and propagation of topological solitons in a chain of coupled parametric-micromechanical-resonator arrays
Hiroshi Yamaguchi and Samer Houri
NTT Basic Research Laboratories, NTT Corporation, Atsugi-shi, Kanagawa 243-0198, Japan ( [email protected] ) Using a coupled parametric-resonator array for generating and propagating a topological soliton in its rotating-frame phase space is theoretically and numerically investigated. In an analogy with the well-known π model, the existence of a soliton is topologically protected as the boundary of two different phase domains of parametric oscillation. Numerical simulation indicates that the propagation can be triggered by switching of the phase state of one specific resonator, and the effects of damping, collision, and the symmetry lifting by harmonic drive on the propagation dynamics are studied. The topological soliton can be implemented by using electromechanical resonators, which allow its propagation dynamics to be precisely electrically controlled and provide a fully controlled on-chip test bed for the study of a topological soliton. I. Introduction A parametric resonator is one of the simplest Floquet systems [1] -[3] that spontaneously breaks time translational symmetry [4] -[6] . Harmonic resonance forces a resonator to oscillate near the eigenfrequency, while parametric modulation at twice the eigenfrequency forces the system to maintain half-period time translational symmetry. Consequently, excited oscillations break the externally defined half-period symmetry and drop into one of the doubly degenerated phase states, each of which is transformed from the other through the half-period time translation. The parametric resonator is frequently described by a pseudo-energy Hamiltonian in its rotating frame, and the two oscillation states correspond to two local minima of a double-well-shaped potential [7] [8] . It can be implemented by using various physical systems, ranging from a macroscopic child swing, optical nonlinear systems [9] -[11] , and LCR circuits [12] [13] to state-of-the-art microscopic on-chip devices like micro- and nano-mechanical systems [14] -[17] and superconducting circuits [18] -[19] . As for these systems, oscillation states can be externally controlled by artificially βmanipulatingβ the double-well potential. Using the controllability, applications such as noise squeezing [14] [18] -[20] , computations [13] [17] [21] , and a symmetry-lifting detector for quantum measurement [22] -[24] have been proposed and experimentally demonstrated. Reflecting the recent progress of device-fabrication technology, interest in a βcoupled mechanical-resonator arrayβ has been increasing. Devices comprising a large number of micro- and nano-electromechanical systems are routinely fabricated, and peculiar collective behaviors, such as mode localization and synchronization, have been studied both theoretically and experimentally [25] -[40] . Among those behaviors, intrinsic localized modes (ILMs), also known as lattice solitons, have been studied [26] ,[32] -[40] . A parametric-resonator array has also been extensively studied and well-described in the rotating frame by using the so-called βparametrically driven nonlinear Schrodinger equationβ (PNSE) [32] -[40] , which was shown to have solitons as solutions. However, the study of topological solitons in a parametric-resonator array is limited, and, to the best of my knowledge, only the stability of a static dark soliton has been discussed [39] [40] . As in the case of other topological edge states, a topological soliton is a kind of topological defect between two different phase domains, and is a commonly observed excitation in one-dimensional systems with discontinuous symmetry breaking [41] . The existence of these solitons is robust in the sense that it is topologically protected against external perturbation, so it has the potential to be an important mode of robust device operation. Accordingly, in this paper, we study topological solitons in a one-dimensional array of parametric resonators. In particular, we numerically study the generation and following propagation dynamics of topological solitons, which can be directly demonstrated with real devices. Considering the results of this study, we discuss the feasibility of using a realistic coupled-electromechanical-resonator array to experimentally demonstrate the robustness of topological solitons. II.
Model The model used for this study is a coupled one-dimensional parametric-resonator array. The equation of motion is given by πΜ π = βπΎπΜ π β π [1 + 2π€ πππ (2π π‘)]π π β πΌπ π3 + ππ (π π+1 β π π )+ ππ (π πβ1 β π π ). (1) where π π is displacement of the n -th resonator ( n : integer) from its equilibrium position, π is resonant (angular) frequency, πΎ is the damping coefficient, π€ is parametric excitation amplitude, πΌ is the Duffing nonlinear coefficient, and π is the coupling constant between two adjacent resonators. For simplicity, no detuning is assumed, i.e., parametric excitation is at exactly twice the resonance frequency, and πΎ , π , π€ , πΌ , and π are assumed to be identical for all resonators. The implemented model based on an electromechanical parametric-resonator array is shown schematically in Fig. 1. Suspended parametric resonators are coupled through the overhang formed through the fabrication process [25] [26] [42] [43] . Parametric oscillation is electrically excited by applying an alternate voltage to the common electrode. A strong coupling regime, i.e. πΎ < π , is assumed, and rotating-frame approximation is applied by introducing slowly varying quadrature amplitudes, π π (π‘) and π π (t) , defined by π π (π‘) = π π[(π π (π‘) + ππ π (π‘))π ππ π‘ ] = π π (π‘) cos π π‘ β π π (π‘) sin π π‘. (2) The equations of motion for new variables become ππ π ππ‘ = β πΎ2 π π + π π€2 π π β 3πΌ8π π π (π π2 + π π2 ) + π2 (π π+1 β 2π π + π πβ1 ), ππ π ππ‘ = β πΎ2 π π + π π€2 π π + 3πΌ8π π π (π π2 + π π2 ) β π2 (π π+1 β 2π π + π πβ1 ). (3) Fig. 1 Schematic drawing of a coupled electromechanical parametric-resonator array and a travelling kink (topological soliton) as a boundary between two different oscillation phase domains. Parametric oscillations are simultaneously excited by applying a voltage to the common electrode (yellow) at twice the resonance frequency, . Each resonator is coupled with two adjacent resonators through an overhang formed by selective etching [25] [26] [42] [43] . Also shown are two quadrature amplitudes, π π β‘ π π /π and π π β‘ π π /π , for a travelling kink solution. Cosine component π π has maximum amplitude at the kink position, which is proportional to the traveling velocity of the soliton. The dynamics in the absence of damping ( πΎ~0 ) can be derived from the canonical equations of motion as ππ π ππ‘ = ππ» π ππ π , ππ π ππ‘ = β ππ» π ππ π , (4) where π» π is the rotating-frame Hamiltonian, given as π» π (π π , π π ) = β π» ππ + π» πππ‘ β‘ β [β π π€4 (π π2 β π π2 ) + 3πΌ32π (π π2 + π π2 ) ] π + π4 β[(π π+1 β π π ) + (π π+1 β π π ) ] π . (5) A three-dimensional plot of π» π for one specific resonator π is shown in Fig. 2, indicating the features of a double-well potential [8] . The two local minima (π π , π π ) = (Β±π , 0) , with π =2π β π€3πΌ , correspond to two steady-state phase states, π π (π‘) = βπ sin π π‘ . The parametrically amplified oscillation amplitude π π saturates at the value of Β±π according to Duffing nonlinearity. Fig. 2 Three-dimensional plot of rotating-frame Hamiltonian π― π , calculated with πΆ = π. ππ, π =π, π = π. ππ . The two fixed points corresponding to parametric-oscillation sates are shown by red markers at (π π , π π ) = (Β±π , π) with π = ππ βπͺ/ππΆ = π. ππ . The dashed closed orbit shows the trajectory of libration. To intuitively find the effect of coupling, the continuous variables (π π , π π ) are approximately replaced by discrete pseudo spin variables, π π β‘ π π /π , which have only two discrete values, π π = Β±1 , to obtain the effective Hamiltonian π» πππ for π π as π» πππ (π π ) β‘ π» π (π π π , 0) = β ππ π π π+1π + πΆπ π‘. (π π = Β±1), (6) where πΆπ π‘. is a π π -independent constant. Equation (6) is identical to a 1D Ising Hamiltonian [21] [43] . The analogy suggests that the phase boundary can exist at any finite temperature because the 1D Ising Hamiltonian does not show phase transition. In Section 4, it is more accurately shown that (i) the phase boundary of parametric-resonator oscillation states corresponds to a topological soliton in the phase space and (ii) the dynamics can be controlled by externally modifying the parameters of parametric excitation. Before the soliton solution is discussed, the weak linear excitation around the synchronized parametric oscillation is first described and used to find the stability and fundamental properties of coupled-parametric-resonator arrays. Without coupling π , each resonator amplitude (π π , π π ) sits at one of two steady-state fixed points (Β±π , 0) in the phase space (red makers in Fig. 2). When the resonator is forced to oscillate by an external perturbation, it departs from the fixed point to start a slow periodic oscillation around that point (dashed line) [44] . This rotating-frame orbit is referred to as βlibrationβ [34] , in an analogy with slow periodic astronomic motions. It is then considered that all parametrically oscillating resonators are synchronized through the coupling to the positive phase state, i.e., (π π , π π ) = (+π , 0) for all values of π , and small deviation βπ§ π is excited in its vicinity, i.e., βπ§ π β‘ π π + π(π π β π ) . It can be seen that the libration around the fixed point is also correlated and propagates as a traveling wave through the coupling. By linearizing Equation (3), the equation of motion for complex variable βπ§ π becomes πβπ§ π ππ‘ = ππ π€βπ§ π β π π2 (βπ§ π+1 β 2βπ§ π + βπ§ πβ1 ). (7) which is the wave equation for a coupled βlibratorβ array with on-site frequency of πΊ = π π€ , and coupling π . The libration waves (LWs) are stable and propagate along the resonator array if all the resonators are synchronized in one specific oscillation phase. Under the assumption of a travelling-wave solution, namely, βπ§ π ~ exp π(πΊ(π)π‘ β πππ₯π₯) , where π₯π₯ is adjacent resonator distance and k is wave number, the dispersion relation is given by πΊ(π) = πΊ β π(cos ππ₯π₯ β 1). (8) Group velocity can be obtained as π£ π = ππΊππ = π₯π₯π sin ππ₯π₯. (9) It therefore follows that π | β€ π₯π₯π β‘ π£ π,πππ₯ . In section IV, our numerical analysis confirms that as well as the soliton, a wave packet of travelling LWs, with velocity upper bound π£ π,πππ₯ , is also generated upon application of an external excitation, i.e. the phase state switching in one specific resonator. III. Continuous model and analytical solutions When considered along the s-axis, the double-well potential shown in Fig. 2 is similar to that of π model [41] . The Hamiltonian of π model is given by π» π = β« [12 π + (ππππ₯ ) + π(π)] ππ₯ , π(π) = 14 (π β 1) . (10) where π(π₯, π‘) is a field variable, and π(π₯, π‘) = ππ(π₯, π‘)/ππ‘ is its canonically conjugate momentum. Potential energy π(π) is double-well-shaped with the minima at π = Β±1 , and it is well known that the model has a topological soliton solution. If the boundary condition, lim π₯βΒ±β π(π₯) = Β±1 , is applied, a phase boundary, called a βkinkβ or βtopological defect,β which can propagate at arbitral velocity, must exist between the positive and negative domains. Strictly speaking, the phase boundary is not a soliton, and it is often referred to as βquasi-soliton,β because of the lack of waveform conservation when two excitations collide [41] . This is because upon the collision the energy is transferred to internal mode excitations. Hereafter, it is shown that a similar topological (quasi-) soliton exists as a phase boundary between two different oscillation states in a parametric-resonator array. The significant difference between the parametric-resonator array and the π model is that the phase boundary is in rotating-frame quadrature amplitudes in the former, whereas is in laboratory-frame coordinates in the latter. Therefore, while oscillation amplitude is finite over nearly the whole range of π₯ , the phase of the oscillation has a soliton-like behavior in parametric-resonator array. This type of soliton is therefore often called a βdarkβ or βgreyβ soliton, which plays a role also in nonlinear fiber optics [45] . A continuous model is used here to describe a static soliton solution. The Hamiltonian can be obtained from Equation (5) by replacing the variables as follows: π π β π (π₯)βπ₯π₯, π π β π(π₯)βπ₯π₯, π π+1 β π π β π₯π₯ ππ (π₯)ππ₯ , π π+1 β π π β π₯π₯ ππ(π₯)ππ₯ β π ππ β β« ππ₯π₯π₯ π(π₯). Hereafter, to simplify the equation, dimensionless time and length units, i.e., π = 1 and π₯π₯ = 1 , are used. It thus follows that π» π (π , π) = β« ππ₯ [β π€4 (π β π ) + 3πΌ32 (π + π ) ] + π4 β« ππ₯ [(ππ ππ₯) + (ππππ₯) ], (11) and the canonical field equations become ππ ππ‘ = πΏπ» π (π , π)πΏπ = π€2 π + 3πΌ8 (π + π )π β π2 (π πππ₯ ), (12) ππππ‘ = β πΏπ» π (π , π)πΏπ = π€2 π β 3πΌ8 (π + π )π + π2 (π π ππ₯ ). By introducing complex variables, π§(π₯, π‘) β‘ π(π₯, π‘) + ππ (π₯, π‘) , Equation (12) becomes a parametrically driven Nonlinear Schrodinger equation (PNSE) [32] -[40] π ππ§ππ‘ = π2 π π§ππ₯ β 3πΌ8 |π§| π§ β π€2 π§Μ . (13) The stability of bright and dark soliton solutions as the ILMs for several parameter ranges has been discussed [39] [40] . Zero-frequency dark mode solitons are focused on hereafter. The boundary condition of dark mode is given by The time-independent kink solution satisfying the boundary condition is given as [36] π (π₯, π‘) = Β±π ππππ (π₯) β‘ Β±π tanh (π₯ β π₯ π ) , c(π₯, π‘) = 0. (15) Solutions with positive and negative signs in Equation (15) are respectively referred to as a βkinkβ and βantikink.β They correspond to the phase boundary between two parametric-oscillation states (see Fig. 1) and a continuous version of a dark soliton obtained in [39] and [40] for small coupling lim π₯βΒ±β π (π₯, π‘) = Β±π , lim π₯βΒ±β π(π₯, π‘) = 0. (14) limit. Kink length π is given by π = β2ππ€ . (16) Equation (16) means that stronger coupling gives a longer kink. This finding is intuitively reasonable because stronger coupling makes the effect of the neighboring resonator stronger, leading to the wider phase-transition region. From Eq. (11), pseudo-energy density π π is given by π π = β π€4 (π β π ) + 3πΌ32 (π + π ) + π4 [(ππ ππ₯) + (ππππ₯) ] = π€
6πΌ [β2 tanh (π₯ β π₯ π ) + tanh (π₯ β π₯ π ) + 12 cosh π₯ β π₯ π ]. (17) Fig. 3 (a) Waveform π ππππ (π₯) and (b) pseudo-energy density Ξ΅ π (π₯) of the static kink solution as a function of normalized position π₯ , where the units used are π = 1 , π = 1 , and π€ /6πΌ = 1 . The waveform and pseudo-energy density of the kink are plotted in Fig. 3(a) and (b). The additional energy is located around the kink position. Due to the space-translational symmetry, the equation has solution of zero-frequency excitation, corresponding to a Nambu-Goldstone mode given by ππ ππππ (π₯)ππ₯ ~ sech ( π₯βπ₯ π ) ; i.e., π (π₯, π‘) = π ππππ (π₯) + π sech (π₯ β π₯ π ) , c(π₯, π‘) = 0 (18) satisfies Equation (12) within the linear approximation in π . In contrast to the π model, Equation (12) is not relativistically invariant. Therefore, an analytical solution with finite travelling velocity cannot be obtained by Lorentz transformation to a moving frame, so the solution is numerically calculated instead. Under the assumption of a travelling wave form as π (π₯, π‘) = π π(π§), π(π₯, π‘) = π πΆ(π§), π§ = π₯ β π£π‘, (19) Equation (12) becomes ππ β²β² (π§) + ππΆ β² (π§) + π(π§) β π(π§) β π(π§)πΆ(π§) = 0 ππΆ β²β² (π§) β ππ β² (π§) β πΆ(π§) β πΆ(π§) β πΆ(π§)π(π§) = 0. (20) where π = ππ€ = π (21) The boundary condition given by Equation (14) becomes lim π§βΒ±β π(π§) = Β±1, lim π§βΒ±β
πΆ(π§) = 0. (22) The solution for π = 1 is then obtained because it is straightforward to obtain a general solution by rescaling, i.e., π§ β π§/βπ and π β π/βπ . Equation (20) thus becomes π β²β² (π§) + ππΆ β² (π§) + π(π§) β π(π§) β π(π§)πΆ(π§) = 0, πΆ β²β² (π§) β ππ β² (π§) β πΆ(π§) β πΆ(π§) β πΆ(π§)π(π§) = 0. (23) If kink velocity π is small enough, the first-order perturbation approximation, given as π(π§) = tanh π§β2 + ππ (π§) + π(π ), πΆ(π§) = ππΆ (π§) + π(π ), (24) can be used with the boundary condition, namely, lim π§βΒ±β π (π§) = 0, lim π§βΒ±β πΆ (π§) = 0. (25) Equation (23) thus becomes π (π§) + (1 β 3 tanh π§β2) π (π§) = 0, πΆ (π§) β πΆ (π§) = tanh π§β2 πΆ (π§) + β12 sech π§β2. (26) The upper equation in (26) has solution π (π§) = π sech with arbitral coefficient π . This solution corresponds to the Nambu-Goldstone mode [Equation (18)] and simply provides a spatially displaced solution, so π (π§) = 0 can be set without loss of generality. The lower equation can be numerically solved, and the solution thereby obtained is shown in Fig. 4(a). The travelling kink solution has non-zero πΆ -quadrature, and its amplitude is proportional to normalized velocity π . The existence of the travelling-wave solution indicates that the waveform of the kink is preserved as the kink propagates. The kink and antikink can therefore be regarded as topological solitons. Fig. 4 (a) Numerically calculated πΆ (π§) . (b) Three-dimensional plot of S(z) and
πΆ(π§) for a travelling-kink solution (blue solid line). The dashed lines are guides for the eye. The solution has a clockwise rotation in phase space along the propagation axis. A three-dimensional plot of the calculated waveform is shown in Fig. 4(b). The solution has a clockwise rotation in phase space along the propagation axis. In contrast to the other topological solitons, where the phase rotation around the kink position is in the laboratory frame, the rotation is in the rotating-frame phase space. Because the oscillation phase at positive and negative infinities is constrained to zero and π , respectively, the boundary cannot disappear. This property protecting the existence of the kink plays a similar role as topological edge states. Number of rotations in phase space can thus be similarly defined as a topological charge, namely, π π β‘ 12π β« (π ππΆππ§ β πΆ ππππ§) +βπ§=ββ ππ§π , (27) which has a value of sgn(π)/2 for a kink solution. For an antikink solution, both π(π§) and
πΆ(π§) reverse signs, so the topological charge has the same sign as a kink. If the helicity, i.e., the number of windings in the propagation direction, is defined as β β‘ π π sgn(π), (28) β is always a positive half. This result simply reflects the fact that the rotation direction around a local-minimum fixed point is always anticlockwise in rotating-frame phase space. IV. Discrete model and numerical analysis To study the detailed propagation dynamics of topological solitons, the time evolution of two quadrature amplitudes was numerically calculated by using a discrete model. The stability of a static solution with no damping limit was studies by similar calculations [39] and [40] . Hereafter, to discuss the experimental feasibility of using a micro-electromechanical resonator array, a travelling wave with finite damping is studied in detail. This study focuses on the following three points of interests. The first point is the effect of damping on the soliton propagation. To experimentally demonstrate the soliton generation and study the propagation dynamics, the soliton should propagate long enough to be observed in a real resonator array with given damping. The second point is a scheme to excite the soliton. In real experiments, switching of an oscillating phase can be induced by electromechanical transduction means [23] [46] . It is confirmed hereafter that a travelling soliton can be generated by phase switching and can be driven by additional harmonic excitation. Such driving is also important to initialize all phase states. The third point is the collision dynamics of two solitons. This point is important because one of the required properties of a soliton is conservation of wave form after the collision process. As for the numerical study, the discrete model given by Eq. (3) is used, and the time evolution of two quadrature amplitudes is calculated by using a standard RungeβKutta method . As in the case of a continuous model, dimensionless amplitudes are introduced as π π = π π π , π π = π πΆ π . (29) A set of equations to be numerically solved is then obtained as ππΆ π ππ‘ = β πΎ2 πΆ π + π€2 π π [1 β πΆ π2 β π π2 ] + π2 [π π+1 β 2π π + π πβ1 ], ππ π ππ‘ = β πΎ2 π π + π€2 πΆ π [1 + πΆ π2 + π π2 ] β π2 [πΆ π+1 β 2πΆ π + πΆ πβ1 ]. (30) Stability of the static and propagating solitons is first confirmed. Time evolutions of π π and πΆ π when π = 0.02 and π€ = 0.01 are calculated using the initial conditions, π π (π‘ = 0) = tanh (ππ) , πΆ π (π‘ = 0) = β0.42π sech π1.2π. (31) These two initial waveforms are approximately obtained from continuous model analysis, i.e., Eq. (15), and the hyperbolic fitting curve applied to the numerical solution given in Fig. 4(a), respectively. The results of the calculation are shown in Fig. 5, which indicates that the kink propagates to the positive/negative direction when given π is positive/negative, as expected from the continuous model analysis in the previous section. It confirms that the propagation velocity π£ ππππ calculated from the numerically obtained π π (π‘) and πΆ π (π‘) is proportional to parameter π , as π£ ππππ ~0.0069π , as expected from the continuous-model analysis. Fig. 5
Time evolutions of π π and πΆ π calculated by using π = 0.02 , π€ = 0.01 , without damping, for (a) π = β0.2 , (b) π = 0 , and (c) π = 0.2 , under the initial condition given by Equation (15) and the fitting curve of Fig. 4(a). Next, the time evolution when the phase of one single resonator is switched is calculated as follows. The calculation starts from the initial condition as π π = β1, πΆ π = 0 for all values of π . Polarity of π is then reversed to +1 at π‘ = 0 , and the subsequent time evolution in πΆ π and π π is calculated. This scenario corresponds to a feasible experiment, where the phase states of one parametric resonator are switched by an external drive signal [23] [46] . The results when π = 0.02 and π€ = 0.01 with a finite damping πΎ are shown in Fig. 6. The calculation results for different values of quality factor π β‘ 1/πΎ are compared to find the effect of damping. First, it is confirmed that a kink in π π , i.e. the topological soliton, is generated by switching phase states at π = 0 resonator, and it propagates in the positive direction while keeping its shape. In addition, amplitude πΆ π has a peak at the kink position as expected from the continuous-model analysis. The straight lines indicated by blue arrows in Fig. 6 shows the generation of weak excitations, which are propagating at constant velocity. This wave packet corresponds to travelling LW discussed in the previous section and is also generated by phase switching. The wave packet has the velocity that equals maximum group velocity π£ π,πππ₯ . This is because the oscillation phase of only one resonator was reversed so that the shortest-wavelength component of LW was generated. The kink is initially generated at the π = 0 resonator and propagates with a slightly lower velocity than π£ π,πππ₯ , but it decelerates as it propagates. With increasing damping, the reduction of velocity becomes more significant, indicating that the two are correlated. Fig. 6
Time evolutions of π π and πΆ π calculated by using π = 0.02 and π€ = 0.01 , and three different values of π . Initially, all the resonators are set at the oscillation state π π = β1, πΆ π = 0 . At π‘ = 0 , variable π was switched to +1 and the subsequent time evolution was calculated. Blue arrows show the wave front of the LW packet. Calculated kink position π ππππ (π‘) is shown in Fig. 7(a) as a function of time π‘ . The kink position is fitted by an exponential function of time, π ππππ (π‘)~π β ππ βπ‘/π π , and the obtained relaxation time π π is shown in (b). Here, π π is nearly equal to the inverse of damping, Q . Because peak amplitude of πΆ π is proportional to velocity π (see Eq. (24)), the damping in πΆ π leads to deceleration of the soliton. The effect of damping is similar to that seen for sine-Gordon and π models [41] . In a topological soliton, the existence of a kink is guaranteed by the topological property, so the kinetic part of its energy should be responsible for the energy dissipation. Fig. 7 (a) Calculated time dependence of kink position for four different values of Q . The kink exponentially decelerates. (b) Relaxation time of kink velocity π π numerically obtained from (a) as a function of Q . Relaxation time π π is nearly equal to the quality factor. The unit-slope dashed line is a guide for the eyes. Hereafter, it is shown that propagation of the topological soliton can be driven by an external harmonic excitation. For a single parametric resonator, an additional harmonic drive breaks the discrete time translational symmetry with period Ο/Ο and lifts the symmetry between two phase states [22] [23] [24] [46] . This symmetry lifting leads to the preference for one specific phase state, and the injection of noise activate the transition. In the case of a topological soliton, the symmetry lifting activates the transition from the unpreferred to preferred phase states around the kink position without the help of noise. The subsequent transition drives the propagation of the topological soliton. For this study, we use the equation of motion with an additional harmonic drive, namely, πΜ π = βπΎπΜ π β π [1 + 2π€ πππ (2π π‘)]π π β πΌπ π3 + ππ (π π+1 β π π )+ ππ (π πβ1 β π π ) + π sin π π‘. (32) where π is harmonic drive amplitude. Equation (30) thus becomes ππΆ π ππ‘ = β πΎ2 πΆ π + π€2 π π [1 β πΆ π2 β π π2 ] + π2 [π π+1 β 2π π + π πβ1 ] β πΉ , ππ π ππ‘ = β πΎ2 π π + π€2 πΆ π [1 + πΆ π2 + π π2 ] β π2 [πΆ π+1 β 2πΆ π + πΆ πβ1 ]. (33) where πΉ = π /4π is dimensionless harmonic-drive amplitude. The harmonic drive increases propagation velocity of the soliton until it saturates by damping. Equation ( ) can be numerically solved as shown in Fig. 8. Applying the harmonic drive moves the kink position, whose propagation velocity saturates after damping time of ~π . As shown in Fig. 8(e), the saturation velocity is proportional to the harmonic drive amplitude and the quality factor. These results indicate that the topological soliton can be driven at constant propagation velocity under finite damping by simultaneously applying a harmonic drive to all the resonators. This symmetry lifting is practically feasible [23] [46] by applying an AC voltage at resonant frequency π to the same electrode as that to which the parametric drive voltage was applied (Fig. 1). The continuous drive can sweep out all the kinks from the focused area and is useful for initializing the oscillation states to one specific phase. Fig. 8 (a)-(d) Calculated time evolutions of π π with a harmonic drive. Both the harmonic drive amplitude ( F ) and quality factor ( Q ) are varied. (e) Saturated propagation velocity as a function of harmonic drive amplitude ( F ). The same parameters as those in Fig. 5(b) were used. Collision of two solitons is also numerically studied using the same model (Fig. 9). Two topological solitons, a kink and an anti-kink, are generated by phase switching at π‘ = 0 and at π = 0 and π = 100 resonators. They then propagate in opposite directions and collide at π‘~5,000 . When damping is small enough (Fig. 9(a)), the two waves maintain enough energy at the collision to survive afterwards. However, if damping is increased (Fig. 9(b)) their energy is reduced at the collision, and the two solitons disappear after the collision. In their place, a strong ILM, or a βbreatherβ mode [41] , at π~50 is generated by absorbing the energy associated with the solitons followed by outward propagation of LWs. Fig. 9
Calculated time evolutions of π π and πΆ π in the collision process of two solitons. Initially, all the resonators are set at oscillation state π π = β1, πΆ π = 0 . At π‘ = 0 , variables π and π were both switched to +1 , and their subsequent time evolution was calculated. The same parameters as those for Fig. 6 were used. Generation of a breather mode and LWs was observed also for π = 10,000 (Fig. 9(a)), where each soliton velocity is reduced by the collision. Some portion of the soliton energy is transferred to the local mode and LWs, leading to reduction of velocity of the LWs. These results indicate that the waveform is not maintained after the collision, as in the case of the π model, where the internal mode absorbs the energy to annihilate the two solitons. The topological soliton in this parametric-resonator array is therefore, rigorously speaking, not a soliton, but it can be referred to as a βquasi-solitonβ as in the case of the π model [41] . V. Discussion Although the topological solitons discussed in this paper are not pure solitons, the existence of an isolated soliton is topologically protected, and they show experimentally interesting features. Especially, the solitons are formed in rotating frames, which can be electrically well controlled in practical devices. Based on the similarity between statically bistable structures, such as a buckled beam [47] and dynamic bistability in parametric resonators, similar topological solitons can be formed by using an array of coupled buckled beams. Such a system would be more accurately analogous to the π model. However, because controlling static bistability is much more difficult than controlling dynamic bistability in parametric oscillators in real devices, our dynamical scheme allows more precise control and is promising for demonstrating on-chip test beds for propagation of topological solitons. Two examples using nanoelectromechanical devices for experimentally demonstrating topological solitons are given hereafter. Our numerical calculation suggests that the typical parameters required to observe the soliton propagation are π > 10 , π > 10 β2 , and Ξ~10 β2 for number of arrays of ~100 . Here, π can be estimated from frequency splitting between symmetric and antisymmetric vibration modes, and π€ can be estimated from the depth of frequency modulation Ξπ as Ξ~Ξπ /π . The first feasible example is given by the coupled doubly clamped beam array, already schematically shown in Fig. 1. In previous reports on two coupled parametric resonators [42] [48] , adjacent beams were coupled through an overhang formed by selective etching. Q is in the order of at room temperature and further increases up to by reducing the temperature. Frequency splitting between symmetric and antisymmetric modes is π₯π π΄βπ ~10 ππ»π§ as compared to resonance frequency f ~300 kHz , so the estimated coupling is π/π = π₯π π΄βπ /f ~0.03 , which satisfies the requirement. However, the modulation depth of present device is in the order of a few kHz, therefore Ξ~10 β3 . The numerical simulation with this parametric modulation depth shows that large-amplitude LWs are generated by phase switching, so it is difficult to clearly experimentally observe the generation of a soliton. Frequency-modulation depth should be increased by one order of magnitude by redesigning frequency-tuning efficiency. Alternatively, from the scale invariance of Eq. (3) under πΎ β π½πΎ, π€ β π½π€, π β π½π, π π β βπ½π π , π π β βπ½π π , and π‘ β π‘/π½ , it is clear that lowering π allows the observation of similar phenomena with longer time scale if Q can be further improved. The second example is given by a one-dimensional electromechanical phonon waveguide fabricated by selective etching of a sacrificial layer through periodically arranged etching holes [31] . Compared to coupled beam resonators, the waveguide structure is advantageous to have larger coupling ( π/π > 0.1 ) and modulation depth ( Ξ~Ξπ /π ~0.015 ) [49] , which are controllable through the device design. The demonstration of the topological soliton discussed in this paper is feasible using such an existing semiconductor-based electromechanical architecture. It is also interesting to extend this study to generate topological soliton in other systems. For example, the bistablity in a Duffing resonator array can be used to form a kink between two domains of high-amplitude and low-amplitude states when the nonlinear resonators are mutually coupled. Especially, the rotations in phase space between these two fixed points are opposite, and different signs of helicity can be formed. The extension to two-dimensional arrays and the effect of disorder and fluctuation are also interesting subjects, but they require further study and will need to be addressed elsewhere. VI. Conclusion Propagation of a topological soliton in a one-dimensional array of coupled parametric resonators was theoretically studied. Both in continuous and discrete models, a kink between two domains with different phase states can propagate stably. The travelling solution has a topological nature, showing a fixed chirality in phase space. The effect of damping on the propagation dynamics was studied, and the velocity of the kink was reduced by the damping. However, even with a finite damping, an additional harmonic excitation lifts the symmetry between two phase states and drives the soliton propagation at fixed velocity. Collision between the kink and antikink solutions is studied and it generates a localized breather mode, while the two solitons disappear if their energies are not large enough. These finding demonstrate the feasibility of using existing semiconductor-based electromechanical resonator array. ACKNOWLEDGMENTS We thank Dr. H. Okamoto, Dr. D. Hatanaka, and Dr. M. Kusoru for their support and helpful comments on the feasibility using devices. We are also grateful to Dr. R. Ohta and Dr. M. Asano for fruitful discussions. References [1]
M. Faraday, βOn a peculiar class of acoustical figures; and on certain forms assumed by groups of particles upon vibrating elastic surfacesβ, Phil. Trans. Roy. Soc. London 121, 299 (1831). [2]
L. D. Landau and E. M. Lifshitz, Mechanics (Pergamon, Oxford, United Kingdom, 1976), 3rd ed. [3]
M. Dykman, βFluctuating nonlinear oscillators: From nanomechanics to quantum superconducting circuitsβ (Oxford University Press, Oxford, United Kingdom, 2012). [4]
A. Leuch, L. Papariello, O. Zilberberg, C. L. Degen, R. Chitra, and A. Eichler, βParametric symmetry breaking in a nonlinear resonatorβ, Phys. Rev. Lett. 117, 214101 (2016). [5]
M. I. Dykman, C. Bruder, N. LΓΆrch, and Y. Zhang, βInteraction-induced time-symmetry breaking in driven quantum oscillatorsβ, Phys. Rev. B98, 195444 (2018). [6]
N. Y. Yao, C. Nayak, L. Balents, and M. P. Zaletel, βClassical discrete time crystalsβ Nature Phys. 16, 438β447 (2020). [7]
J. Woo and R. Landauer, βFluctuations in a parametrically excited subharmonic oscillatorβ, IEEE Quantum Electron. 7, 435 (1971). [8]
M. Marthaler and M. I. Dykman, βSwitching via quantum activation: A parametrically modulated oscillatorβ, Phys. Rev. A 73, 042108 (2006). [9]
J. A. Giordmaine and R. C. Miller, βTunable coherent parametric oscillation in LiNbO3 at optical frequenciesβ, Phys. Rev. Lett. 14, 973 (1965). [10]
S. Brosnan and R. Byer, βOptical parametric oscillator threshold and linewidth studiesβ, IEEE J. Quantum Electron. 15, 415 (1979). [11]
A. Dutt, K. Luke, S. Manipatruni, A. L. Gaeta, P. Nussenzveig, and M. Lipson, βOn-chip optical squeezingβ, Phys. Rev. Appl. 3, 044005 (2015). [12]
FitzGerald and F. George F. βOn the driving of electromagnetic vibrations by electro-magnetic and electrostatic enginesβ, The Electrician, 28, 329 (1892) [13]
E. Goto, βThe parametron, a digital computing element which utilizes parametric oscillationβ, Proc. IRE. 47, 1304 (1959). [14]
D. Rugar and P. GrΓΌtter, βMechanical parametric amplification and thermomechanical noise squeezingβ, Phys. Rev. Lett. 67, 699 (1991). [15]
K. L. Turner, S. A. Miller, P. G. Hartwell, N. C. MacDonald, S. H. Strogatz, and S. G. Adams, βFive parametric resonances in a microelectromechanical systemβ, Nature 396, 149 (1998). [16]
A. Dana, F. Ho, and Y. Yamamoto, βMechanical parametric amplification in piezoresistive gallium arsenide microcantileversβ, Appl. Phys. Lett. 72, 1152 (1998). [17]
I. Mahboob and H. Yamaguchi, βBit storage and bit flip operations in an electromechanical oscillatorβ, Nature Nanotechnol. (2008). [18]
B. Yurke, P. G. Kaminsky, R. E. Miller, E. A. Whittaker, A. D. Smith, A. H. Silver, and R. W. Simon, βObservation of 4.2-K equilibrium-noise squeezing via a Josephson parametric amplifierβ, Phys. Rev. Lett. 60, 764 (1988). [19]
M. A. Castellanos-Beltran, K. D. Irwin, G. C. Hilton, L. R. Vale, and K. W. Lehnert, βAmplification and squeezing of quantum noise with a tunable Josephson metamaterialβ, Nature Phys. 4, 929 (2008). [20]
J. Abadie et al., βA gravitational wave observatory operating beyond the quantum shot-noise limitβ, Nature Phys. 7, 962 (2011). [21]
T. Inagaki et al., βA coherent Ising machine for 2000-node optimization problemsβ, Science 354, 603 (2016). [22]
D. Ryvkine and M. I. Dykman, βResonant symmetry lifting in a parametrically modulated oscillatorβ, Phys. Rev. E 74, 061118 (2006). [23]
I. Mahboob, C. Froitier, and H. Yamaguchi, βA symmetry-breaking electromechanical detectorβ, Appl. Phys. Lett. 96, 213103 (2010). [24]
Z. R. Lin, K. Inomata, K. Koshino, W.D. Oliver, Y. Nakamura, J.S. Tsai, and T. Yamamoto, βJosephson parametric phase-locked oscillator and its application to dispersive readout of superconducting qubitsβ, Nature Communications 5, 4480 (2014). [25] E. Buks and M. L. Roukes, βElectrically tunable collective response in a coupled micromechanical arrayβ, J. Microelectromech. Syst. 11, 802 (2002). [26]
M. Sato, B.E. Jubbard, and A. J. Sievers, βColloquium: Nonlinear energy localization and its manipulation in micromechanical oscillator arraysβ, Rev. Mod. Phys. 78, 137 (2006) [27]
M. Spletzer, A. Raman, H. Sumali, and J. P. Sullivan, βHighly sensitive mass detection and identification using vibration localization in coupled microcantilever arraysβ, Appl. Phys. Lett. 92, 114102 (2008). [28]
E. A. Martens, S. Thutupalli, A. FourriΓ¨re, and O. Hallatschek, βChimera states in mechanical oscillator networksβ, PNAS, 110, 10563 (2013). [29]
M. Ludwig and F. Marquardt, βQuantum many-body dynamics in optomechanical arraysβ, Phys. Rev. Lett. 111, 073603 (2013). [30]
M. H. Matheny, J. Emenheiser, W. Fon, A. Chapman, A. Salova, M. Rohden, J. Li, M. H. de Badyn, M. PΓ³sfai, L. D.-Osorio, M. Mesbahi, J. P. Crutchfield, M. C. Cross, R. M. DβSouza, and M. L. Roukes, βExotic states in a simple network of nanoelectromechanical oscillatorsβ, Science 363, 1057 (2019). [31]
D. Hatanaka, I. Mahboob, K. Onomitsu, and H. Yamaguchi, "Phonon waveguides for electromechanical circuitsβ, Nature Nanotechnology, 9, 520-524 (2014). [32]
I. V. Barashenkov, M. M. Bogdan, and V. I. Korobov, βStability diagram of the phase-locked solitons in the parametrically driven damped nonlinear Schrodinger equationβ, Europhys. Lett. 15, 113 (1991). [33]
S. Trillo, M. Haelterman, and A. Sheppard, βStable topological spatial solitons in optical parametric oscillatorsβ, Opt. Lett. 22, 970 (1997). [34]
N. V. Alexeeva, I. V. Barashenkov, and G. P. Tsironis, βImpurity-induced stabilization of solitons in arrays of parametrically driven nonlinear oscillatorsβ Phys. Rev. Lett. 84, 3053 (2000). [35]
Ron Lifshitz and M. C. Cross, βResponse of parametrically driven nonlinear coupled oscillators with application to micromechanical and nanomechanical resonator arraysβ, Phys. Rev. B 67, 134302 (2003) [36]
I.V. Barashenkov, S. R.Woodford, and E.V. Zemlyanaya, βParametrically driven dark solitonsβ, Phys. Rev. Lett. 90, 054103 (2003). [37]
H. Susanto, Q. E. Hoq, and P. G. Kevrekidis, βStability of discrete solitons in the presence of parametric drivingβ, Phys. Rev. E 74, 067601 (2006). [38]
E. Kenig, B. A. Malomed, M. C. Cross, and R. Lifshitz, βIntrinsic localized modes in parametrically driven arrays of nonlinear resonatorsβ, Phys. Rev. E 80, 046202 (2009). [39]
M. Syafwan, H. Susanto, and S. M. Cox, βDiscrete solitons in electromechanical resonatorsβ, Phys. Rev. E 81, 026207 (2010). [40]
O. P. Swami, V. Kumar, B. Suthar, and A. K. Nagar, βStability of intersite dark solitons in a parametrically driven discrete nonlinear Schrodinger equationβ, Nanosystems: Physics, Chemistry, Mathematics, 10, 391 (2019) [41] T. Dauxois and M. Peyrard, βPhysics of Solitonsβ (Cambridge University Press, Cambridge, United Kingdom, 2006). [42]
H. Okamoto, A. Gourgout, C.-Y. Chang, K. Onomitsu, I. Mahboob, E. Y. Chang, and H. Yamaguchi, βCoherent phonon manipulation in coupled mechanical resonatorsβ Nature Phys. 9, 480 (2013). [43]
I. Mahboob, H. Okamoto, and H. Yamaguchi, βAn electromechanical Ising Hamiltonianβ, Sci. Adv. 2, e1600236 (2016). [44]
J. S. Huber, G. Rastelli, M. J. Seitner, J. KΓΆlbl, W. Belzig, M. I. Dykman, and E. M. Weig, βSpectral evidence of squeezing of a weakly damped driven nanomechanical modeβ, Phys. Rev. X 10, 021066 (2020). [45]
Y. S.
Kivshar and B. Luther-Davies, βDark optical solitons: physics and applicationsβ, Phys. Rep. 298, 81 (1998) and reference therein. [46]
H. B. Chan, M. I. Dykman, and C. Stambaugh, βSwitching-path distribution in multidimensional systemsβ, Phys. Rev. E, 78, 051109 (2008). [47]
S. M. Carr, W. E. Lawrence, and M. N. Wybourne, βStatic buckling and actuation of free-standing mesoscale beamsβ, IEEE Trans. on Nanotech. 4, 655 (2005). [48]