Suppressing birhythmicity by parametrically modulating nonlinearity in limit cycle oscillators
aa r X i v : . [ n li n . AO ] J u l Suppressing birhythmicity by parametrically modulating nonlinearity in limit cycleoscillators
Sandip Saha, ∗ Sagar Chakraborty, and Gautam Gangopadhyay † S. N. Bose National Centre for Basic Sciences,Block-JD, Sector-III, Salt Lake, Kolkata-700106, India Department of Physics, Indian Institute of Technology Kanpur, Kanpur, Uttar Pradesh 208016, India ‡ Multirhythmicity, a form of multistability, in an oscillator is an intriguing phenomenon foundacross many branches of science. From an application point of view, while the multirhythmicityis sometimes desirable as it presents us with many possible coexisting stable oscillatory states totap into, it can be also be a nuisance because a random perturbation may make the system settleonto an unwanted stable state. Consequently, it is not surprising that there are many naturaland artificial mechanisms available that can control the multirhythmicity. What we propose in thispaper is a hitherto unexplored mechanism of controlling birhythmicity—the simplest nontrivial formof the multirhythmicity. Our main idea is to incorporate parametric (periodic) modulation of thenonlinear damping in the limit cycle oscillators with a view to exciting resonance and antiresonanceresponses at particular angular driving frequencies, and controlling the resulting birhythmicity bychanging the amplitude of the modulation. To this end, we employ analytical (perturbative) andnumerical techniques on the van der Pol oscillator—a paradigmatic limit cycle system—havingadditional position dependent time delay term and its modified autonomous birhythmic version.We also bring the fact to the fore that introduction of delay—a commonly adopted method ofcontrolling multirhythmicity—in such a system can sometimes bring forth unwanted birhythmicity;and interestingly, our method of controlling birhythmicity through periodic modulation can suppresssuch a delay induced birhythmic response.
Keywords: Multistability; limit cycle; delay; perturbative methods; van der Pol oscillator
I. INTRODUCTION
Since Faraday’s observation [1] of parametric oscilla-tions as surface waves in a wine glass tapped rhyth-mically, almost two centuries have passed and over theyears, it has been realized that the phenomenon of para-metric oscillations is literally omnipresent [2, 3] in physi-cal, chemical, biological, and engineering systems. Para-metric oscillations are essentially effected by periodicallyvarying a parameter of an oscillator which, thus, is aptlycalled a parametric oscillator. The simplest textbook ex-ample with wide range of practical applications is theMathieu oscillator [4] where the natural frequency of asimple harmonic oscillator is varied sinusoidally and theinteresting phenomenon of parametric resonance [5] is ob-served. The effect of additional nonlinearity in the Math-ieu oscillator has also been extensively investigated, e.g,in Mathieu–Duffing [6, 7], Mahtieu–van-der-Pol [8, 9, 10],and Mathieu–van-der-Pol–Duffing [11, 12, 13, 14] oscil-lators. However, only rather recently, the effect of pe-riodically modulating the nonlinearity in a limit cyclesystem, viz., van der Pol oscillator has been investi-gated [15]. The resulting parametric oscillator, termedPENVO ( p arametrically e xcited n onlinearity in the v ander Pol o scillator), along with the standard phenomenonof resonance, exhibits the phenomenon of antiresonance ∗ [email protected] † [email protected] ‡ [email protected] that is said to have occurred if there is a decrease in theamplitude of the limit cycle at a certain frequency of theparametrical drive (cf. [16, 17, 18]).In the context of the limit cycle oscillations [19], oneis readily reminded of the limit cycle systems possess-ing more than one stable limit cycle. A plethora ofsuch multicycle systems are manifested in biochemicalprocesses [20, 21, 22, 23, 24, 25, 26, 27]; one of thesimplest of them being a multicycle version of the vander Pol oscillator [28, 29, 30] modelling some biochem-ical enzymatic reactions. This oscillator has two stablelimit cycles (and an unstable limit cycle between them inthe corresponding two dimensional phase space) owingto the state dependent damping coefficient that has upto sextic order terms. Consequently, it shows birhyth-mic behaviour wherein depending on the initial condi-tions, the long term asymptotic solution of the oscilla-tor corresponds to one of the stable limit cycles thathave, in general, different frequencies and amplitudes.Needless to say, birhythmicity is a widely found phe-nomenon across disciplines—and not only in biochemi-cal processes—because so are the ubiquitous limit cycleoscillations.Since different initial conditions lead to different solu-tions for a birhythmic oscillator, the inherent uncertaintyin the amplitude and the frequency in the eventually real-ized stable oscillations can be gotten rid of if the oscillatoris somehow made monorhythmic. It is known [31, 32, 33]from the studies on the Hénon map and rate equationsof laser that while a small change of one of the systemparameters of a birhythmic oscillator may not in generalconvert it to a monorhythmic system, an external controlin the form of a slow periodic parameter modulation canannihilate one of the coexisting attractors resulting ina monostable oscillatory system. Technically speaking,birhythmicity is a simple type of multistability which, inother words, mean coexistence of different attractors atfixed parameter values in the system. The existence ofmultistability in diverse systems and the need to controlit are elaborately discussed in a review article [34] whichalso reviews various control strategies including their ex-perimental realizations.Interestingly, time delay is known to have significanteffect on the attractors of a nonlinear system and can alsobrings forth new ones. For example, even in a relativelysimple system like the Rössler oscillator, time delayedfeedback control [35] induces a large variety of regimes,like tori and new chaotic attractors, nonexistent in theoriginal system; furthermore, the delay modifies the pe-riods and the stabilities of the limit cycles in the systemdepending on the strength of the feedback and the magni-tude of the delay. As another example, we may point outthat the direct delayed optoelectronic feedback can sup-press hysteresis and bistability in a directly modulatedsemiconductor laser [36]. The coexistence of two stablelimit cycles with different frequencies in the presence ofdelayed feedback has been discussed in detail [37] for thevan der Pol oscillator and its variants. Mutlicycle van derPol oscillator has also been investigated from the pointof view of control of birhythmicity using some differentforms of time delay [30, 38, 39, 40, 41, 42, 43].However, to the best of our knowledge, there hasbeen no investigation into the control of multistabilityin a parametric oscillator whose parameter, determin-ing the strength of the nonlinear term, is varied. Itshould be noted that periodic variation of such a pa-rameter is not inconceivable [44]; in fact, it can resultin parametric spatiotemporal instability leading to in-teresting time-periodic stationary patterns in reaction-diffusion systems. In view of the above, it is imperativethat an investigation of the PENVO and its relevant ex-tension be carried out and the interplay, if any, betweenthe time-delayed feedback and the parametric forcing berevealed.To this end, in this paper, we first discuss in Sec. IIhow presence of time delayed feedback affects the res-onance and the antiresonance in the PENVO. Further-more, we discuss how the resulting birhythmicity thereinis supressed by tuning the strength of the period modu-lation. Subsequently, in Sec. III, we consider multicyclePENVO—multicycle van der Pol oscillator whose nonlin-earity is sinusoidally varying—and argue in detail thatit is possible to control birhythmicity in this system aswell. Finally, we reiterate the main results of this paperin Sec. IV. FIG. 1: Limit cycles in PENVO with delay haveoscillating amplitudes.
We time-evolve Eq. (1) with γ = 1 . , K = µ = 0 . , τ = 0 . for Ω = 2 (black) and 4 (red) to arrive at the correspondingtime-series plots (subplot a), x vs. t , and phase spaceplots (subplot b), ˙ x vs. x . II. PENVO WITH DELAY
Even a simple harmonic oscillator with its quadraticpotential modified so as to have a term that is time de-layed, exhibits nontrivial dynamics. The resulting solu-tions, including the oscillatory ones, in the weak non-linear limit can be iteratively extracted using perturba-tive methods based on the concept of renormalizationgroup [45, 46]. An extended version of the delayed sim-ple harmonic oscillator, that possesses limit cycle, hasalso been analyzed [47] using the Krylov–Bogoliubovmethod [48, 49]. Motivated by these results, we nowconsider the PENVO with a time delay term as follows: ¨ x + µ [1 + γ cos(Ω t )]( x −
1) ˙ x + x − Kx ( t − τ ) = 0 , (1)where < K, µ, τ ≪ ; γ ∈ R ; and Ω ∈ R + .Note that for K = γ = 0 , we get back the van derPol oscillator that in weak nonlinear limit shows stablelimit cycle oscillations with amplitude 2. For appropri-ate non-zero values of γ ( K still zero), we arrive at theequation for the PENVO [15] that is known to show an-tiresonance (oscillations with amplitude smaller than 2)and resonance (oscillations with amplitude greater than2) at Ω = 2 and
Ω = 4 respectively.
Our specific goal inthis section is to find out what happens to the resonanceand the antiresonance states once the time delay is intro-duced (i.e., when
K, γ = 0 and Ω = 2 , ), and to explorethe possible existence of birhythmicity and its control inthe system. To begin with we have extensively searched for numer-ical solutions of Eq. (1) at different parameter values. InFig. 1, we present two particular oscillatory solutions forthe cases
Ω = 2 and
Ω = 4 . We note that the limitcycles have oscillating amplitudes. In order to under-stand the origin of oscillating amplitude and to discoverbirhythimicity in the course of our investigation, we em-ploy the Krylov–Bogoliubov method on Eq. (1). We,thus, make an ansatz: x ( t ) = r ( t ) cos( t + φ ( t )) where wehave adopted polar coordinate, ( r, φ ) = ( √ x + ˙ x , − t + t0.61.21.82.4 r (a) γ 〈 r 〉 (c) FIG. 2:
Anti-resonant responses with oscillating amplitudes in PENVO with delay.
This figure panel has beengenerated by time-evolving Eq. (1) with γ ∈ [0 , , K = µ = 0 . , τ = 0 . ; and Ω = 2 (black) and 4 (red) . Thetime-series, r vs. t , (subplot a) depicts oscillating limit cycles in the PENVO with delay and the reason behind theoscillations is best understood as the corresponding non-circular limit cycle attractors in the p - q plane (subplot b).While for subplots (a) and (b), γ = 1 . , subplot (c) showcases the variation of the averaged amplitudes with γ , thus,highlighting the presence of antiresonances ∀ γ ∈ [0 , .FIG. 3: Strength of periodic modulation of nonlinear damping controls delay-induced birhythmicity.
This figure panelof streamline plots depicts repellers [unstable focus (red dot) and saddle (orange dot)] and attractors [stable focus(blue dot) and stable limit cycle (around red dot; not explicitly shown)] in p - q space of the PENVO with delay at γ = 1 . , . , and . K = µ = 0 . τ = 0 . ; and Ω = 2 . The stable foci on(approximately) principle diagonal of the figures have same p p + q -value, and so is the case with the stable focion (approximately) anti-diagonal of the figures. Note how with change in γ -value, the number of attractors changesfrom one (limit cycle) to four (foci that have only two distinct p p + q -value). tan − ( − ˙ x/x )) . r and φ are very slowly varying functionof time since we are working under the assumption that < µ ≪ ; we set r ( t ) = r + O ( µ ) and φ ( t ) = φ + O ( µ ) .Here, we have used the definition that average of a func-tion, f ( x, ˙ x ) (say), over a period π is conveniently de-noted as f ( t ) = (1 / π ) R π f ( s ) ds . Furthermore, Taylor-expanding r ( t − τ ) as r ( t − τ ) = r ( t ) − τ ˙ r ( t ) = r ( t )+ O ( µ ) (since ˙ r ( t ) ∼ O ( µ ) ), one finally obtains ˙ r = − r (cid:0) K sin τ + µ (cid:0) r − (cid:1)(cid:1) A Ω ( r, φ ; γ ) + O ( µ ) , (2a) ˙ φ = − K cos τ B Ω ( r, φ ; γ ) + O ( µ ) , (2b)where, O ( µ ) terms can be neglected and A Ω and B Ω denote the γ dependent parts. It is interesting that thesetwo functions’ denominators blow up at Ω equal to and . We, thus, resort to the L’Hôspitals’ rule to find thefunctions at Ω = 2 , : A ( r, φ ; γ ) = − γµr cos(2 φ ) , (3a) B ( r, φ ; γ ) = − γµ sin(2 φ ) (cid:0) r − (cid:1) ; (3b) A ( r, φ ; γ ) = 116 γr µ cos(4 φ ) , (3c) B ( r, φ ; γ ) = − γr µ sin(4 φ ) . (3d)Here the subscripts specify the value of Ω at which A Ω and B Ω have been determined.As an illustration, in Fig. 2(a), we present r as a func-tion of t for both Ω = 2 and
Ω = 4 after fixing γ = 1 . , τ = 0 . , and K = µ = 0 . . The solutions are oscilla-tory in sharp contrast to the case of the weakly nonlinearvan der Pol oscillator for which the plot of r vs. t wouldbe a horizontal straight line passing through r = 2 at large times. Obviously, it is a little ambiguous to definethe resonance and the antiresonance states in terms ofthe magnitude of the oscillations’ amplitude because theamplitude itself is oscillating. Hence for the sake of con-sistency, to define the resonance and the antiresonancestates, we henceforth use the average of the oscillatingamplitude. Consequently, in Fig. 2(c), we plot averageof r i.e. h r i t (after removing enough transients) with γ to note that at both Ω = 2 and
Ω = 4 the system showsantiresonance. Note that one of the interesting effects ofthe delay is to suppress the uncontrolled growth of oscil-lations (for
Ω = 4 and as γ → ) present in the absenceof delay.The oscillations in the amplitudes of the limit cyclesis best explained by recasting the equations for r and φ in ( p, q ) -plane where ( p, q ) = (cid:0) r cos φ, r sin φ (cid:1) or conse-quently, ( r, φ ) = ( p p + q , tan − ( q/p )) . Substitutingthese relations in equations (2), one arrive at the follow-ing dynamical flow equations: ˙ p | = − Kp sin τ Kq cos τ − µp − γµp µp γµpq − µpq , (4a) ˙ q | = − Kp cos τ − Kq sin τ − γµp q − µp q − µq γµq µq (4b) ˙ p | = − Kp sin τ Kq cos τ γµp − µp µp − γµpq − µpq , (5a) ˙ q | = − Kp cos τ − Kq sin τ − γµp q − µp q + 116 γµq − µq µq . (5b)Here again subscripts and refer respectively to thecases corresponding to Ω = 2 and
Ω = 4 . Fig. 2(b)exhibits the limit cycles that are not perfect circles aboutthe origin in p - q plane. Thus, it is clear that for either ofthe cases, the slow variation of the limit cycle amplitudeis manifested through the slow variation of the distanceof the phase point on the closed trajectory from the originin p - q plane.Now, we ask the question if the system allows forbirhythmicity. We realize that a convenient way to searchfor it is to look for stable fixed points (except the one atthe origin) and stable limit cycles in the corresponding p - q plane. A closer look at Eqs. (4) and (5) reveals that (0 , is a common fixed point and, additionally, we haveseen that they possess limit cycles. Straightforward lin-ear stability analysis about the fixed point for the case Ω = 4 yields (cid:0) µ ± iKe ± iτ (cid:1) / as the eigenvalues thatclearly has real negative part and there is no local bi-furcation possible with change in γ . In fact, detailed nu-merical study suggests that, for the appropriately fixedparameters and Ω = 4 , no changes occur except that theoscillation in the amplitude of the limit cycle becomesless perceptible with increase in γ . Naturally, one ex- pects only monorhythmicity in the system.The case of Ω = 2 is, however, very interesting:The linear stability about (0 , yields the eigenvalues ( ± p γ µ − K cos(2 τ ) − K − K sin τ + 2 µ ) / andthus the character of the fixed point can change with thevalue of γ , e.g., it is quite clear that for small values of γ (other parameters being appropriately fixed) the originshould be a focus and for larger values it should be asaddle. The full study of Eq. (4) being analytically quitecumbersome, we present a numerical illustration of howbirhythmicity is generated by varying γ .In this respect, please see Fig. 3 where we have depictedthe vector plots corresponding to Eq. (4) for γ = 1 . , γ = 2 . and γ = 3 . . We have fixed Ω = 2 , τ = 0 . ,and K = µ = 0 . . Careful study reveals that, as γ is in-creased, after γ ≈ . the origin becomes a saddle froman unstable focus. The saddle however is born along withtwo stable foci (say, F − and F +1 ) at which the stablemanifolds of the saddle terminate; two other stable fociare also born (say, F − and F +2 ) and the limit cycle, thatexists around the origin for γ . . , is annihilated. Oneobserves that at a given γ , the value of p + q is same for F − and F +1 , and also for F − and F +2 , meaning that onlyFIG. 4: Birhythmic response of PENVO with delay.
The time series plot (a) and the phase space plot (b) forEq. (1) with γ = 3 . , K = µ = 0 . , τ = 0 . , and Ω = 2 . The blue and the black lines correspond to twodifferent initial conditions.two (and not four) different limit cycles can be observedin the PENVO with delay when γ & . . We verify thisconclusion by numerically solving Eq. (1) for two differ-ent initial conditions but at the same set of parametervalues and as shown in Fig. 4, we observe birhythmic os-cillations. To conclude what we have shown is that bychanging γ we can induce birhythmicity or conversely,one can say that if the system is already birhythmic, wecan make the system monorhythmic by using γ as a con-trol parameter. III. MULTICYCLE PENVO
Up to now we have seen how a delay term added in thePENVO modifies the antiresonance and the resonance at
Ω = 2 and
Ω = 4 respectively, and furthermore, givesrise to birhythmicity that in turn can be controlled bythe strength of the periodically modulated nonlinearityin PENVO. Another natural modification of the van derPol oscillator with multiple limit cycles is a variant ofthe van der Pol oscillator—originally proposed [28, 29]to model enzyme reaction in biochemical system—witha sextic order polynomial as damping coefficient: ¨ x + µ ( − x − αx + βx ) ˙ x + x = 0 . (6)Here, < µ ≪ and α, β > . We call it Kaiser os-cillator. It has three concentric limit cycles surround-ing an unstable focus at the origin: two of them arestable and the unstable one acts as the boundary sep-arating the basins of attractions of the two stable cy-cles. However, whether there are two stable limit cycles(birhythmicity) or only one (monorhythmicity) strictlydepends on values of α and β . Under the assumptionthat µ ≪ , straightforward application of the Krylov–Bogoliubov method helps to demarcate the regions ofbirhythmicity and monorhythmicity in α − β parameterspace (see Fig. 8 in Appendix A). In the context of thispaper, it is of immediate curiosity to ponder upon theimportant questions like ‘can one find resonance and an-tiresonance in the Kaiser oscillator’, ‘would periodically 〈 r 〉 (c) 〈 r 〉 (d) FIG. 5:
Resonant and antiresonant responses inmulticycle PENVO.
Presented are time series plots(subplot a and b) corresponding to both small (solidline) and large (dotted line) cycles for
Ω = 2 (black) , , .Furthermore, subplots (c) and (d) depict how theaveraged amplitudes of the responses change with γ ∈ [0 , . It is depicted that the smaller limit cycleshows resonances for the case Ω = 4 , and butantiresonance for the case Ω = 2 ; the larger limit cycleadmits resonance for
Ω = 8 but antiresonance for thecase
Ω = 2 , and . The values of the parameters usedto numerically solve Eq. (8) for the purpose of thefigure are α = 0 . , β = 0 . , µ = 0 . and γ = 1 . (insubplot a and b).modulating the nonlinearity control the inherent birhyth-micity in the Kaiser oscillator’, etc.The addition of the periodic modulation of nonlinearityin the Kaiser oscillator get us the following equation: ¨ x + µ [1 + γ cos(Ω t )] ( − x − αx + βx ) ˙ x + x = 0 , (7)where γ > . For obvious reasons, henceforth we aptlycall this system: multicycle PENVO. Again, the Krylov–Bogoliubov method yields, ˙ r = 1128 rµ (cid:0) − βr + 8 αr − r + 64 (cid:1) + A Ω ( r, φ ; γ ) , (8a) ˙ φ = B Ω ( r, φ ; γ ) + O ( µ ) . (8b)Here the symbols are in their usual meaning as detailedin Sec. II. The subscripts specify the value of Ω at which A Ω and B Ω have to be determined; the functions havesingularities at Ω = 2 , , and , and their limiting val-ues at these Ω -values are respectively, A = − γrµ cos(2 φ ) (cid:0) βr − αr + 16 (cid:1) , (9a) B = − γµ sin( φ ) cos( φ ) (cid:0) βr − αr + 16 r − (cid:1) ; (9b) A = 164 γr µ cos(4 φ ) (cid:0) βr − αr + 4 (cid:1) , (9c) B = − γr µ sin(4 φ ) (cid:0) βr − αr + 8 (cid:1) ; (9d) A = − γr µ cos(6 φ ) (cid:0) α − βr (cid:1) , (9e) B = 1128 γr µ sin(6 φ ) (cid:0) α − βr (cid:1) ; (9f) A = 1256 βγr µ cos(8 φ ) , (9g) B = − βγr µ sin(8 φ ) . (9h)As before, we go on to p - q plane to recast set of equa-tions (8) for all four Ω -values in terms of p and q vari-ables (see Appendix B) in order to understand the dy-namics conveniently. For all the four values of Ω , theorigin— p, q =(0,0)—is a fixed point that on doing linearstability analysis, turns out to be unstable for all valuesof γ . Since now the corresponding equations of motionare much more cumbersome to handle analytically, weresort to a numerical investigation of the systems. Firsthowever we need to pick appropriate value of α and β .We choose α = 0 . and β = 0 . that would allowthe Kaiser oscillator (multicycle PENVO with γ = 0 )to exhibit birhythmicity (see Appendix A); the ampli-tudes of the limit cycles that are concentric circles about ( x, ˙ x ) = (0 , in the limit µ → are approximately . and . respectively. In what follows, we work with µ = 0 . .We now turn on the periodic modulation of the nonlin-ear term, i.e., we work with the multicycle PENVO withnonzero γ . We scan the system for various values of γ andpresent the results for γ up to in Fig. 5. For illustrativepurpose, consider γ = 1 . . We note that the amplitudeof the smaller limit cycle of the Kaiser oscillator increasesfor the case Ω = 4 , and (resonances) but decreasesfor the case Ω = 2 (antiresonance). Similarly, while theamplitude of the larger limit cycle of the Kaiser oscillatorincreases for the case
Ω = 8 (resonance), but it decreasesfor the case
Ω = 2 , and (antiresonances). As an aside,for the case Ω = 6 , we also note that the amplitudes ofboth the cycles themselves oscillate and the response cor-responding to the outer limit cycle changes from antires-onance to resonance as γ increases (see Fig. 5d).More interesting, however, is the fact that the res-onance and the antiresonance, manifested as limit cy-cles with oscillating amplitudes, for Ω = 6 merge—asimplicitly shown in Fig. 6(a)—for a range of γ -values: γ ∈ ( γ c , γ c ) ≈ (0 . , . . This means that γ is yet again acting as a control parameter in bringing about monorhythmicity by suppressing the birhythmic-ity. To understand the phase dynamics of control ofthe aforementioned birhythmicity, we consider the sys-tem (8) in ( p, q ) plane at three representative valuesof γ , viz., γ = 0 . (Fig. 6b), γ = 1 . (Fig. 6c), and γ = 1 . (Fig. 6d). For γ = 0 . < γ c , a case of birhyth-micity, there are twelve stable nodes—the only attractorsin the phase space—that can be classified into two groupssuch that one group of nodes has p p + q ≈ . andthe other group has p p + q ≈ . . This correspondsto the fact that there are two distinct limit-cycles in the x - ˙ x plane, and their radii are . and . ; in otherwords, the system is birhythmic. In the monorhythmiccase of γ = 1 . ∈ ( γ c , γ c ) , we note that the attractorsnow are twelve limit cycles whose centers (unstable focus)lie on a circle of radius . (approximately). Thus, thesystem has now become monorhythmic and the limit cy-cle in the x - ˙ x plane has periodically oscillating amplitude.The bifurcation leading to the creation of the twelve sym-metrically placed limit cycles takes place at γ = γ c whenthe stable nodes and the unstable saddles (present at γ < γ c ) merge appropriately to give rise to the limitcycles (seen at γ > γ c ). Finally, For γ = 1 . > γ c ,the system showcases birhythmic behaviour yet again:the six symmetrically placed asymptotically stable nodesin the corresponding p - q plane have identical values for p p + q , viz. , . that corresponds to the amplitude ofthe limit cycle of the multicycle PENVO.We note that the birhythmicity present at other res-onance and antiresonance conditions, i.e., for Ω =2 , , and , could not be controlled to monorhythmicityby the variation in γ . However, recalling that in Sec. IIthe combination of γ and delay could effect control ofbirhythmicity, one is tempted to add delay term, viz.,‘ − K ( t − τ ) ’ in the left hand side of Eq. (7) with a hopeto effect control of birhythmicity for Ω = 2 , , and .The introduction on such a delay term in the Kaiseroscillator shifts the region of birhythmicity in the α - β plane (see Appendix. A). In the simultaneous presenceof non-zero γ and K , the multicycle PENVO’s responseat Ω = 2 , , , and can be analyzed using the Krylov–Bogoliubov method just as has been done in detail forEq. (1) and Eq. (7). We omit the repetitive details andrather present the summary of the analyses in Fig. 7(a-b).We note that the delay does indeed suppress birhythmic-ity; and interestingly in the case of Ω = 8 , γ can be seento be a control parameter even in the presence of delay. IV. DISCUSSION AND CONCLUSIONS
How to control birhythmicity in an oscillator is an in-teresting question. In this paper we have illustrated thatthe birhythmicity seen in the delayed van der Pol os-cillator and the van der Pol oscillator modified to havehigher order nonlinear damping (the Kaiser oscillator)can be suppressed if the nonlinear terms of the oscillatorsare periodically modulated. This periodic modulation γ 〈 r 〉 (a) γ c γ c FIG. 6:
Strength of periodic modulation of nonlinear damping controls birhythmicity in multicycle PENVO.
Subplot(a) presents the observation that the average amplitudes of the periodic responses—the smaller limit cycle (solidblue line) and the larger limit cycle (dotted blue line)—merge for an intermediate range of γ between γ c ≈ . to γ c ≈ . resulting in monorhythmicity. Streamplots (b)-(d) depict repellers [unstable node (black dot) , unstablefocus (red dot) and saddle (orange dot)] and attractors [stable node (green dot) and stable limit cycle (around eachred dot; not explicitly shown)] in p - q space of the multicycle PENVO at γ = 0 . , . , and . , respectively. Otherparameter values have been fixed at α = 0 . , β = 0 . , µ = 0 . and Ω = 6 . In subplot (b), there are two sets ofstable foci with two distinct values of p p + q (hence birhythmicity), while in subplot (c) only attractor (and hencemonorhythmicity) is a limit cycle—a circle that passes through all the unstable foci with same p p + q -values andcentred at origin. In subplot (c), in addition to this limit cycle, another set of stable foci appear with same p p + q -value (hence birhythmicity). 〈 r 〉 (a) 〈 r 〉 (b) FIG. 7:
Controlling birhythmicity via delay in multicyclePENVO.
Subplots (a) and (b) exhibit how the averagedamplitudes change with γ ∈ [0 , corresponding to bothsmall (solid line) and large (dotted line) cycles for Ω = 2 (black) , , . Thevalues of the relevant parameters used in the figure are α = 0 . , β = 0 . , µ = 0 . and τ = 0 . . We notethat the responses are mostly monorhythmic.of the nonlinear damping also brings about resonanceand antiresonance responses in the aforementioned oscil-lators. In order to characterize the responses, we havepresented perturbative calculations using the Krylov–Bogoliubov method and supplemented them with am-ple numerical solutions for the systems of ordinary dif-ferential equations under consideration. We have alsodiscussed in detail how to understand the bifurcationsleading to monorhythmicity from birhythmicity (and viceversa) from the relevant phase space trajectories obtainedvia the perturbative technique.We recall that the introduction of delay is one of thepopularly known method of controlling birhythmicity. However, as we have seen in Sec. II, delay can intro-duce birhythmicity as well. It is interesting to realizein such cases periodically modifying the nonlinear termscan change the birhythmic behaviour to monorhythmic.A comparison of responses due to delay and parametricexcitation in a limit cycle system provides an extra tool-kit for controlling birhythmicity when one alone may notbe fruitful. We may point out that the delay term wehave used in this paper is completely position dependentas opposed to the more commonly investigated velocitydependent delay terms [40, 41, 42] in the literature.We strongly believe that the proposed idea of control-ling multirhythmicity by invoking periodic modulation ofnonlinear terms could be useful in plethora of limit cyclesystems. It is also worth pondering if such a mechanismof suppressing multirhythmicity is present in nature be-cause, after all, there is no dearth of the limit cycle os-cillations [19] in nature. However, we do not believe thatbuilding a general universal mechanism behind this phe-nomenon can be proposed easily; each system has to beanalysed on a case-by-case basis. ACKNOWLEDGMENT
SS acknowledges RGNF, UGC, India for the partialfinancial support. SS is grateful to Rohitashwa Chat-topadhyay for his enormous support during a visit to IITKanpur, and Pratik Tarafdar for some help with Mathe-matica. SC is thankful to Anindya Chatterjee (IIT Kan-pur) for insightful discussions.
Appendix A: Birhythmicity in the Kaiser Oscillator:Effect of Delay
Consider the Kaiser model in presence of a positiondependent delay: ¨ x + µ (cid:0) − x − αx + βx (cid:1) ˙ x + x − Kx ( t − τ ) = 0 , (A1) (0 < ǫ, τ ≪ ). When K = 0 , the system is eithermonorhythmic or birhythmic depending on the values of α and β as depicted in Fig. (8). It is expected that forsmall values of K and τ , the behaviour of the Kaiser oscil-lator should be qualitatively similar, although the regionin the α - β plane where the birhythmic behaviour is seenwould be shifted slightly. This is shown in Fig. (8) thathas been obtained by employing the Krylov–Bogoliubovmethod to write the equations for the amplitude as wellas the phase of the system’s response as ˙ r = − r (cid:0) K sin τ + µ (cid:0) βr − αr + 16 r − (cid:1)(cid:1) , (A2a) ˙ φ = − K cos τ, (A2b)respectively. Here higher order terms have been ne-glected. It is clear from the existence of non-overlappingregions of birhythmicity that introducing delay may in-duce monorhythmicity in birhythmic cases or vice versa. α . . β FIG. 8:
Delay changes rhythmicity.
This figureshowcases for what values of α and β , systems (6) and(A1) are birhythmic—the green and the red zonesrespectively. In other words, the changes in thebirhythmic zone in α - β parameter space in the presenceof the time delay ( K = 0 . and τ = 0 . ) have beenexhibited. The systems are monorhythmic when notbirhythmic. Here, µ = 0 . . Appendix B: Flow Equations: Multicycle PENVOwith Delay
On imposing parametric excitation to the nonlinearityin Eq. (A1), we can write, ¨ x + µ [1 + γ cos(Ω t )] (cid:0) − x − αx + βx (cid:1) ˙ x + x − Kx ( t − τ ) = 0 . (B1)The corresponding amplitude and phase equations are ˙ r = − r (cid:0) K sin τ + µ (cid:0) βr − αr + 16 r − (cid:1)(cid:1) + A Ω ( r, φ ; γ ) + O ( µ ); (B2a) ˙ φ = − K cos τ + B Ω ( r, φ ; γ ) + O ( µ ) , (B2b)where higher order terms have been neglected, and A Ω and B Ω are functions with singularities at Ω = 2 , , and . One may resort to the L’Hôspitals’ rule and go to p - q plane to rewrite the amplitude and the phase equationsin terms of the coordinate of the plane: ˙ p = − Kp sin τ Kq cos τ − βγµp − βµp + 164 αγµp + 116 αµp + 332 βγµp q − βµp q − µp βγµp q − βµp q − αγµp q + 18 αµp q − γµp µp βγµpq − βµpq − αγµpq + 116 αµpq + 14 γµpq − µpq , ˙ q = − Kp cos τ − Kq sin τ − βγµp q − βµp q − βγµp q − βµp q + 1164 αγµp q + 116 αµp q − βγµp q − βµp q + 532 αγµp q + 18 αµp q − γµp q − µp q + 164 βγµq − βµq − αγµq + 116 αµq − µq γµq µq p = − Kp sin τ Kq cos τ βγµp − βµp − αγµp + 116 αµp + 964 βγµp q − βµp q + 116 γµp − µp − βγµp q − βµp q − αγµp q + 18 αµp q + µp − βγµpq − βµpq + 732 αγµpq + 116 αµpq − γµpq − µpq , ˙ q = − Kp cos τ − Kq sin τ − βγµp q − βµp q − βγµp q − βµp q + 732 αγµp q + 116 αµp q + 964 βγµp q − βµp q − αγµp q + 18 αµp q − γµp q − µp q + 164 βγµq − βµq − αγµq + 116 αµq + 116 γµq − µq µq p = − Kp sin τ Kq cos τ βγµp − βµp − αγµp + 116 αµp − βγµp q − βµp q − µp − βγµp q − βµp q + 532 αγµp q + 18 αµp q + µp βγµpq − βµpq − αγµpq + 116 αµpq − µpq , ˙ q = − Kp cos τ − Kq sin τ − βγµp q − βµp q + 1564 βγµp q − βµp q + 564 αγµp q + 116 αµp q + 332 βγµp q − βµp q − αγµp q + 18 αµp q − µp q − βγµq − βµq + 164 αγµq + 116 αµq − µq µq p = − Kp sin τ Kq cos τ βγµp − βµp + 116 αµp − βγµp q − βµp q − µp
8+ 35256 βγµp q − βµp q + 18 αµp q + µp − βγµpq − βµpq + 116 αµpq − µpq , ˙ q = − Kp cos τ − Kq sin τ − βγµp q − βµp q + 35256 βγµp q − βµp q + 116 αµp q − βγµp q − βµp q + 18 αµp q − µp q + 1256 βγµq − βµq + 116 αµq − µq µq . The subscript indicates the value of Ω in Eq. (B1) for which the pair of above first order equations are writtenin ( p , q ) coordinates. [1] M. Faraday, Philosophical transactions of the Royal So-ciety of London , 299 (1831)[2] A. H. Nayfeh and D. T. Mook, Nonlinear oscillations (Wiley Classic Library ed., New York, 1995)[3] M. E. Marhic,
Fiber optical parametric amplifiers, oscil-lators and related devices (Cambridge university press,UK, 2008)[4] É. Mathieu, J. Math. Pures Appl , 137 (1868)[5] L. D. Landau and E. Lifshitz, Mechanics , 121 (1988)[6] E. Esmailzadeh and G. Nakhaie-jazar,International Journal of Non-Linear Mechanics , 905 (1997)[7] D. J. Braun, Phys. Rev. Lett. , 044102 (2016)[8] F. Lakrad, A. Azouani, N. Abouhazim, and M. Belhaq,Chaos, Solitons & Fractals , 813 (2005)[9] M. Momeni, I. Kourakis, M. Moslehi-Fard, and P. K. Shukla,Journal of Physics A: Mathematical and Theoretical , F473 (2007)[10] F. Veerman and F. Verhulst,Journal of Sound and Vibration , 314 (2009)[11] J. Li, W. Xu, X. Yang, and Z. Sun,Journal of Sound and Vibration , 330 (2008)[12] M. Belhaq and A. Fahsi, Nonlinear Dynamics , 139 (2008)[13] M. Pandey, R. H. Rand, and A. T. Zehnder, NonlinearDynamics , 3 (2008)[14] L. Lu and X. Li, J. Math. Phys. , 122703 (2015)[15] S. Chakraborty and A. Sarkar,Physica D , 24 (2013)[16] D. J. Ewins, Modal testing: theory and practice , Vol. 15(Research studies press Letchworth, 1984)[17] D. B. Saakian, Phys. Rev. E , 016126 (2005)[18] P. Sarkar and D. S. Ray, Phys. Rev. E , 052221 (2019)[19] A. Jenkins, Phys. Rep. , 167 (2013)[20] F. Morán and A. Goldbeter,Biophys. Chem. , 149 (1984)[21] E. E. Sel’kov, Euro. J. Biochem. , 79 (1968)[22] M. Morita, K. Iwamoto, and M. Sen,Phys. Rev. A , 6592 (1989)[23] J.-C. Leloup and A. Goldbeter,J. Theo. Bio. , 445 (1999)[24] I. D. la Fuente, Biosystems , 83 (1999)[25] M. Stich, M. Ipsen, and A. S. Mikhailov,Phys. Rev. Lett. , 4406 (2001)[26] M. Stich, M. Ipsen, and A. S. Mikhailov, Physica D , 19 (2002)[27] S. Kar and D. S. Ray,Phys. Rev. Lett. , 238102 (2003)[28] F. Kaiser, “Theory of resonant ef-fects of rf and mw energy,” in Biological Effects and Dosimetry of Nonionizing Radiation ,edited by M. Grandolfo, S. M. Michaelson, and A. Rindi(Springer US, Boston, MA, 1983) pp. 251–282[29] F. Kaiser and C. Eichwald,Int. J. Bifurc. Chaos , 485 (1991)[30] R. Yamapi, B. R. Nana Nbendjo, and H. G. E. Kadji,Int. J. Bifurc. Chaos , 1343 (2007)[31] A. N. Pisarchik and B. K. Goswami,Phys. Rev. Lett. , 1423 (2000)[32] A. N. Pisarchik, Y. O. Barmenkov, and A. V. Kir’yanov,Phys. Rev. E , 066211 (2003)[33] B. K. Goswami and A. N. Pisarchik,Int. J. Bifurc. Chaos , 1645 (2008)[34] A. N. Pisarchik and U. Feudel,Phys. Rep. , 167 (2014), control of multistabil-ity[35] A. G. Balanov, N. B. Janson, and E. Schöll,Phys. Rev. E , 016222 (2005)[36] S. Rajesh and V. Nandakumaran, Physica D , 113 (2006)[37] T. Erneux and J. Grasman,Phys. Rev. E , 026209 (2008)[38] H. G. E. Kadji, R. Yamapi, and J. B. Chabi Orou,Chaos , 033113 (2007)[39] H. G. E. Kadji, J. B. Chabi Orou, R. Yamapi, andP. Woafo, Chaos, Solitons & Fractals , 862 (2007)[40] P. Ghosh, S. Sen, S. S. Riaz, and D. S. Ray,Phys. Rev. E , 036205 (2011)[41] D. Biswas, T. Banerjee, and J. Kurths,Phys. Rev. E , 042226 (2016)[42] D. Biswas, T. Banerjee, and J. Kurths,Chaos , 063110 (2017)[43] D. Biswas, T. Banerjee, and J. Kurths,Phys. Rev. E , 062210 (2019)[44] S. Ghosh and D. S. Ray, Phys. Rev. E , 032209 (2016)[45] S.-i. Goto, Prog. Theo. Phys. , 211 (2007)[46] A. Sarkar, J. K. Bhattacharjee, S. Chakraborty, andD. B. Banerjee, Euro. Phys. J. D , 479 (2011)[47] S. Saha and G. Gangopadhyay,J. Math. Chem. , 750 (2019)[48] N. M. Krylov and N. N. Bogolyubov, Introduction to non-linear mechanics (Princeton Univ. Press, 1947)[49] J. K. Bhattacharjee, A. K. Malik, and S. Chakraborty,Indian J. Phys.81