Citation |

- Permanent Link:
- http://digital.auraria.edu/AA00007209/00001
## Material Information- Title:
- Modeling the response of Rayleigh - Van der Pol oscillators to stochastic excitation near the Hopf bifurcation
- Creator:
- Lu, Bingshen
- Place of Publication:
- Denver, CO
- Publisher:
- University of Colorado Denver
- Publication Date:
- 2019
- Language:
- English
## Thesis/Dissertation Information- Degree:
- Master's ( Master of integrated sciences)
- Degree Grantor:
- University of Colorado Denver
- Degree Divisions:
- College of Liberal Arts and Sciences, CU Denver
- Degree Disciplines:
- Integrated sciences
- Committee Chair:
- Asadi-Zeydabadi, Masoud
- Committee Members:
- Tagg, Randall
Simon, Burt
## Notes- Abstract:
- In this thesis we investigate the Rayleigh-van der Pol oscillator with additive noise. This dynamical system is a fundamental model that can be tested experimentally by the Wien-bridge electrical circuit. The analysis of the Rayleigh-van der Pol oscillator in the presence of noise provides a good understanding of the stochastic nonlinear dynamical system. The Rayleigh-van der Pol oscillator is a two-variable nonlinear system with one control parameter. As the control parameter varies at the bifurcation (transition) point the state of system changes from steady to self-sustained oscillation. This transition is known as Hopf bifurcation. The focus of this project is on the behavior of the system near bifurcation point in the presence of noise. We found that the bifurcation point widens to a â€œbifurcation regionâ€ and the â€bifurcation regionâ€ shifts as the level of noise increases.
## Record Information- Source Institution:
- University of Colorado Denver
- Holding Location:
- Auraria Library
- Rights Management:
- Copyright Bingshen Lu. Permission granted to University of Colorado Denver to digitize and display this item for non-profit research and educational purposes. Any reuse of this item in excess of fair use or other copyright exemptions requires permission of the copyright holder.
## Auraria Membership |

Downloads |

## This item has the following downloads: |

Full Text |

MODELING THE RESPONSE OF RAYLEIGH - VANDER POL OSCILLATORS TO STOCHASTIC EXCITATION NEAR THE HOPF BIFURCATION by BINGSHEN LU B.S., University of Colorado Denver, 2015 A thesis submitted to the Faculty of the Graduate School of the University of Colorado in partial fulfillment of the requirements for the degree of Master of Integrated Sciences Integrated Sciences 2019 This thesis for the Master of Integrated Sciences degree by Bingshen Lu has been approved for the Integrated Sciences Program by Masoud Asadi-Zeydabadi, Chair Randall Tagg Burt Simon Date: May 18, 2019 m Lu, Bingshen (M.S., Integrated Sciences) Modeling the Response of Rayleigh - Van der Pol Oscillators to Stochastic Excitation Near the Hopf Bifurcation Thesis directed by Associate Professor Masoud Asadi-Zeydabadi ABSTRACT In this thesis we investigate the Rayleigh-van der Pol oscillator with additive noise. This dynamical system is a fundamental model that can be tested experimentally by the Wien-bridge electrical circuit. The analysis of the Rayleigh-van der Pol oscillator in the presence of noise provides a good understanding of the stochastic nonlinear dynamical system. The Rayleigh-van der Pol oscillator is a two-variable nonlinear system with one control parameter. As the control parameter varies at the bifurcation (transition) point the state of system changes from steady to self-sustained oscillation. This transition is known as Hopf bifurcation. The focus of this project is on the behavior of the system near bifurcation point in the presence of noise. We found that the bifurcation point widens to a â€œbifurcation regionâ€ and the â€™â€™bifurcation regionâ€ shifts as the level of noise increases. The form and content of this abstract are approved. I recommend its publication. Approved: Masoud Asadi-Zeydabadi IV TABLE OF CONTENTS I INTRODUCTION......................................................... 1 II ESSENTIAL IDEAS FROM NONLINEAR DYNAMICS.............................. 7 11.1 Linear Systems................................................... 7 11.11 Simple Harmonic Motion with Damping.............................. 9 11.111 Limit Cycle and Hopf Bifurcation................................ 11 II.IV A Nonlinear System: the Classic Van der Pol Model............... 11 II.V Modified Van der Pol Model...................................... 18 II. VI Rayleigh-van der Pol Model...................................... 21 III ESSENTIAL IDEAS FROM STOCHASTIC DIFFERENTIAL EQUATION (SDE).............................................................. 27 111.1 Classical Example of Stochastic Process: Brownian Motion........ 27 111.11 Generating a Stochastic Process Numerically..................... 27 111.111 Embedding Stochastic Processes in Differential Equations....... 28 III. IV Numerical Methods for SDE...................................... 28 III.IV.I Euler-Maruyama method......................................... 29 III.V Solving Stochastic Integrals..................................... 29 III.V.I Itoâ€™s Rule.................................................... 29 III.V.II Itoâ€™s Chain Rule.............................................. 30 III. V.Ill Analytical examples......................................... 30 IV STOCHASTIC DYNAMICS................................................ 37 IV. I Stochastic Van der Pol Model.................................... 37 IV. II Stochastic Rayleigh-van der Pol................................. 37 V BEHAVIOR NEAR THE BIFURCATION...................................... 42 V. I Numerical Simulation............................................. 42 V.I.I Poincare Section.............................................. 42 v V.I.II Behavior of Mean and Variance with Bifurcation Parameter for Fixed Noise Level................................ 46 VI CONCLUSION AND FUTURE WORK............................. 60 VI.I Future Work.......................................... 60 VI.I.I Coupled System..................................... 60 REFERENCES................................................ 62 APPENDIX A. ANALYSIS AND SIMULATION FOR COUPLED SYSTEM............. 64 vi CHAPTER I INTRODUCTION Bifurcation, a sudden change of the dynamical system, lies in the center of nonlinear dynamics[29][19]. It has plenty of applications in daily activity and fields of science such as the quick change of medical condition from stable to severe, financial market change from balanced to crisis[ 19], or neural network change from stationary to generates motion. The needs to do a proper prediction of such phenomena rise as time goes by. However, reality has uncertainties so it is important to understand how randomness effects the nonlinear systems. The general theory of random excitation in physical systems was first studied by Albert Einstein in 1905 and was generalized by Marian Smoluchowski in 1916[8]. Later on, physicists and mathematicians Andrey Kolmogorov, Adriaan Fokker, Max Planck, Leonard Ornstein, George Eugene Uhlenbeck, Norbert Wiener, Joseph L. Doob and Kiyoshi Ito laid the foundation for solving linear SDE between 1930s and 1950s[8], However, those earlier works are for solving linear SDE. Difficulties arise for nonlinear system because of the complexity nature of nonlinear systems[8]. Back in 1957, Garstens started studying noise in nonlinear oscillators[12], Garstens used a nonlinear Van der Pol equation with driving noise term to model nonlinear oscillators that â€™â€™the output is fed back into the system thus modifying in a complicated manner the effective noise inputâ€[12] here the nonlinear terms in such model are multiplicative (i.e. Xi X2, where ^ =X2). Garstens presented a method 1 for estimating the effect of nonlinear noise in an oscillator at low levels of oscillation. This required an approximate solution of a non-homogeneous van der Pol type equation with a noisy driven force by reducing the equation to a linear form. The power spectrum of the output was used to study the effect of noise in the nonlinear system. Caugheys paper in 1959[6] showed that there was interest in the how nonlinear systems response to random excitation, especially in the systems with multiplicative nonlinear terms. He studied the self-excited system like Van der Pol Oscillator and its response to external stochastic terms acting as external random excitation. Caughey found that the solution has a stochastic/noise term because of the stochastic term in the self-excitation system, and that the strength of the stochasticity decreases as the amplitude of oscillation increase. Interestingly enough, so does the band width of the stochasticity/noise which prove the early work done by Garstens incorrect. In Caugheys paper in 1967[7], the use of a Fokker-Planck equation had been introduced to study the stationary response of self-excited oscillators. Using such technique, Caughey obtained the exact solution for the stationary joint probability density function, as well as the spectral density for oscillator with weak damping and weak excitation. Following his previous work, Caughey developed a nonlinear theory of random vibrations[8]. As Caughey pointed out, the fieldâ€™s importance had rise rapidly because of the developments in fields like high speed flight and the needs of dealing the random fluctuation cause by turbulence back at that time. Caughey summarized and developed some solution techniques for random vibrations in nonlinear systems[8]. During the similar time in Japan (1974), Kawakubo and Kabashima also used the Fokker-Planck approach to study stochastic processes in self-excited oscillation using a Van der Pol stochastic nonlinear system in a more experimental way[17], Kawakubo and Kabashima focused on how the fluctuations, namely the amplitude and phase 2 fluctuations in the system, behave in the experiment of Wien-Bridge oscillators. They showed how the voltage fluctuates with the feedback factor and how the variance fluctuates with the parameter controlling the nonlinearity. Kawakubo and Kabashima also included the probability distribution function obtained with the Fokker-Planck equation, examining how it changed with the parameter controlling the nonlinearity. In agreement with Caughey and others find. Spanos showed that in 1978 such problems were drawing more and more attention[28]. This included several more applications of stochastic phenomena in physics and engineering, such as control, aerospace, natural disaster, ocean wave and transportation engineering. Spanos pointed out that even though linear random vibration has be well studied, the theory of non-linear random vibrations was not. This is still the case: there are still a lot of development needed in this topic. The Fokker-Planck approach was used in Spanos work to obtain solutions, Spanos also showed the application of Fokker-Planck approach in Van der Pol and Rayleigh oscillators. Later on in 1995, Leung did a systematic study of the behavior of Van der Pol oscillator with additive and multiplicative stochastic terms using different nonlinearity parameters[20], comparing the bifurcation behaviors of the system with different stochastic terms. By such methods Leung stated that we could get a more realistic understanding of the bifurcation of a nonlinear system. Leung found that all level of stochastic terms / noise chosen disturb the usual behavior of the Van der Pol oscillator. Moreover, he found that noise could postpone the bifurcation between fixed point and limit cycle and that noise could induce a transition from a stable fixed point to a stable limit cycle, and vice versa. Leung suggested that noisy nonlinear systems could have more and different phenomena than the deterministic one, such as stochastic resonance. 3 In Leungâ€™s later work in 1998[21], he pointed out studies on stochastic bifurcation focus on the shift of bifurcation point, where a shifting of bifurcation point is an indication of the noise effect on the two competing mechanisms, namely self-sustained oscillation and nonlinear dissipation. Yuan and others made further developments in this direction in 2013 and found that zero-mean multiplicative noise can induce the bifurcation points and sometimes even change the topology of the system[31], or in another words, the structure does not correspond to the deterministic counterpart, which make the analysis of such system more difficult. They also used the steady-state distribution of the process to approach such case, in where correspond means the probability density function does not change over time and have a maximum at the stable attractors, which later on define as the steady state distribution multiplying negative one being a Lyapunov function of the deterministic counterpart of the stochastic dynamics. Yuans group apply a new stochastic interpretation rather than the common Ito and Stratonovichs integration to obtain a correspondence between stochastic and deterministic dynamics. Just last year (2018), Doan and others were trying to merge bifurcation theory for stochastic dynamical systems[10], Doan and others work shows that the earlier concepts developed in the 1990s, namely theory of phenomenological (P) bifurcation and dynamical (D) bifurcation, cannot comprehensively capture the complexity of bifurcation in random dynamical system. Doan points out that despite its relevance for many applications of topical interest, a bifurcation theory of random dynamical systems is still in its infancy[10], With the developments made in mid 20th century on time dependent probability theory[2][16][26] and stochastic differential equations[25], and we began to have tools to solve the stochastic system in the form of stochastic differential equation[13][25]. Yet The bifurcation of stochastic chaotic system is still something we do not sufficiently understand yet. But the importance of bifurcation in such system is becoming 4 increasingly unavoidable[l][10]. In recent years, neuroscience (as well as machine learning and artificial intelligent) draws more and more attention[22] and funding, so do the nonlinear and stochastic components of this field[l][22][15]. There are some breakthroughs on the applications of such topics in crucial areas such as neuroscience[14][l 1][17], However, most of the developments about stochastic differential equations are in the field of stock market forecasting and related topics in quantitative/financial engineering[25], i.e. the early detection of economic bubbles[9]. How to analyze the stochastic nonlinear dynamical model is still a difficult task which been rarely discussed in the context of stochastic differential equation but mostly done by deterministic approximation, especially for nonlinear system with multiplicative terms in neural networks (i.e. System contains Xi X2, where ^ = X2). Most of the analytical and numerical methods cannot be applied to such terms[5][18][25]. Yet the interest in self-excited oscillation and its application in neural network pushes us forward to study such system[3][17], especially the response of a class of oscillators to stochastic excitation[28]. Now we might have tools and methods to study it in the content of stochastic differential equation instead of deterministic approximation methods. In this thesis we are building a modified Rayleigh Van der Pol model with additive stochastic components. We hope to test if the system behaves as others observed, when we put the system in the proper context of stochastic differential equation instead of deterministic methods. We also hope to lay a foundation for future work and potential applications of the Rayleigh Van der Pol model. The applications include neural networks for the Rayleigh Van der Pol[6][7] based central pattern generator (CPG) model[3], as well as other potential applications in physics, neuroscience, finance and other related fields. This work is also offered as a guide for some methods of how to simulate and analyze models of nonlinear system with stochastic components in the context of 5 stochastic differential equations. This is primarily focused on how to solve the system with terms containing random variable products Xi â€¢ X2, where = X2, which we do not have many good options to solve yet. 6 CHAPTER II ESSENTIAL IDEAS FROM NONLINEAR DYNAMICS This chapter organizes a conceptual framework for analyzing two-dimensional nonlinear dynamics. This provides a broad and powerful foundation for analyzing many applications especially after stochastic forcing is added in later chapters. II.I Linear Systems In general a two-dimensional dynamical linear system has the form of X-\ â€” dy\X-\ + d-\2X2 (II. la) X2 = 921-^1 + 922-^2 (II.lb) where x 1 = â€”, x 2 = â€”, and an, 9i2, 9i3, 3u are constant parameters. dt dt We can write the system in matrix form as where A x = Ax [?n ai2D ^1 â–¡ and x = â–¡ â–¡ a2i a22 â–¡ â–¡ x2 We seek solution of the form [29], (II 2) 7 x=eAv. (II. 3) Substituting (2.3) into (2.2), we get Av = Av, (II 4) We say the solution in (2.3) is an eigensolution if v is an eigenvector of A with corresponding eigenvalue A [29], The characteristic equation of matrix A is given by Equation (2.5) yields di 1 ~A a-i2 det C&21 ^22 â€”A o, (II. 5) where A2 -tA + A = 0 (II. 6) t = trace(A) = a-n+ 322 A = det(A) = 311 322 -312321 (II. 7a) (II. 7b) Solving equation (2.6), we get M = A2 â€” V y+ 7-2 _4A V 2 T - T2 - 4A 2 (II. 8 a) (II. 8b) where A1 and ^2 are the eigenvalues of matrix A. We find the eigenvectors by substituting eigenvalues into one of the equations in (2.4). We get two linearly independent eigenvectors Vi and V2. Now we can get the general solution for x, x=c-\eMtv-\ +C2eA2tv2. (II 9) 8 II.II Simple Harmonic Motion with Damping One good example for the linear system is simple harmonic motion with damping. The F in this case is the sum of the spring and damping forces: dx F=-kx-b~ (II. 10) where k and b are the spring and damping constants. Applying Newtonâ€™s second law for (2.10), we get (fix dx m_____= -kx - bâ€” dt2 dt We can rewrite equation (2.11) as or where (3=^, U)0 = (fix dx m_____+ b____+ kx = 0 dt2 dt x + /3x + (fiQx = 0 The solution (2.12) and (2.13) has the form x = Ae*( Substituting (2.14) into (2.12), then we have (11.11) (11.12) (11.13) (11.14) (,mA2 + bA + k)AeAt = 0 Equation (2.15) for the nontrivial solution of (2.12) implies >2 b i k A2 +___A +____= o. m m Equation (2.16) is the characteristic equation of (2.12). (11.15) (11.16) 9 The roots of equation (2.16), eigenvalues of (2.12), for underdamped or oscillatory case, (|L < J-), is given by b k b2 A = â€” Â±i - 9 O m m A (11.17) The general solution of (2.12) or (2.13) for weak damping, 4hs< is given by b j K x = Aie^te m - b2 f A - b~ f Am- -h A2e2nl m Am~ (11.18) or V V x = A e^e' -fPXtt a gjSfg-/ w2-p2/4t (11.19) 1 2 For real valued solutions, A-\ and A2 will be complex conjugates. We call x = 0, x = 0 a fixed point solution. For (3 < 0 all solutions tend to this fixed point solution and it is therefore said to be stable, for /3 > 0 the fixed point is unstable. Figure (2.1) shows time series and phase plot for damped harmonic motion. Figure II. 1: Underdamped Harmonic Motion 10 II.III Limit Cycle and Hopf Bifurcation A limit cycle occurs when the phase plane (x v.s. contains an isolated closed trajectory. Isolated means the neighboring trajectories either spiral towards or away from the limit cycle. [29] A Hopf bifurcation happens where a fixed point of a system changes to a limit cycle. II.IV A Nonlinear System: the Classic Van der Pol Model The Van der Pol equation has been used as a model for neuron networks with feedback controls. It is named after Dutch physicist Balthasar van der Pol[29][30], who applied the equation to electronic circuit oscillators. CpX 2 Cfx 2 dt2 ~(P )Wo dt + W0X = (11.20) Provided that [J /= 0 we can write: CpX ( 2 cfx 2 dt2 ~^1 ~fjx )'wÂ° dt + = 0 (II.21) Using the rescaled variables, t =U)ot and x = ^ equation (2.21) becomes: cfx dt2 ~/J( 1 -x ) dx 2 â€” dt + x = 0 (11.22) This is the standard, or classical form of the van der Pol equation. However, the scale factor ^becomes singular at jt/ = 0 at the same point that the actual output, x, has an amplitude that is just about to rise from a value of zero[30]. This is one of the reasons that we modified the classical Van der Pol model as shown in the later section.(see section 2.3 and 2.4) An alternate dimensionless form of the dynamics is the single scaled nonlinear 11 equation in x: 12 2 dx cPx dp-^-x)d;+x=0' (11.23) Here, the scaled time and dynamic variable are t = U)ot and x = - where x* =V -ipor notational clarity, we will drop the tildeâ€™s in the following but assume the equations are still in the scaled form. We can applying method of averaging to the Van der Pol equation to find approximate solutions to the nonlinear system. The dynamics, equation (2.23), can be re-written as two first-order differential equations dx = w dt dw 2 (II. 24a) dt=~x + Qj -x )w. (II. 24b) We are going to assume that the the solution will evolve from the simple sinsoidal solution when the nonlinear term is turned on. Thus we think of time variations in amplitude and phase around the limit cycle such that the phase-space trajectory can be written in the form: x(0 = /\(0cos[f+ w(t) = -A(t) sin[f+
Our goal is to find equations for the amplitude function A(t) and the phase function In order to do so, we insert (2.25a) and (2.25b) for x and w into the dynamical equations (2.24a) and (2.24b) to obtain dA(t) dw(t}
13
Collect terms and use the fact that cos2[f + + sin2[f + |

Full Text |

PAGE 1 MODELING THE RESPONSE OF RAYLEIGH VAN DER POL OSCILLATORS TO STOCHASTIC EXCITATION NEAR THE HOPF BIFURCATION by BINGSHEN LU B.S., University of Colorado Denver, 2015 A thesis submitted to the Faculty of the Graduate School of the University of Colorado in partial fulfillment of the requirements for the degree of Master of Integrated Sciences Integrated Sciences 2019 PAGE 2 iii This th esis for the Master of Integrated Sciences degree by Bingshen Lu has been approved for the Integrated Sciences Program by Masoud Asadi Zeydabadi, Chair Randall Tagg Burt Simon Date: May 18, 2019 PAGE 3 iv Lu, Bingshen (M.S., Integrated Sciences) Modeling the Response of Rayleigh Van der Pol Oscillators to Stochastic Excitation Near the Hopf Bifurcation Thesis directed by Associate Professor Masoud Asadi Zeydabadi ABSTRACT In this thesi s we investigate the Rayleigh van der Pol oscillator with additive noise. This dynamical system is a fundamental model that can be tested experimentally by the Wien bridge electrical circuit. The analysis of the Rayleigh van der Pol oscillator in the presence of noise provides a good understanding of the stochastic nonlinear dynamical system. The Rayleigh van der Pol oscillator is a two variable nonlinear system with one control parameter. As the control parameter varies at the bifurc ation (transition) point the state of system changes from steady to self sustained oscillation. This transition is known as Hopf bifurcation. The focus of this project is on the behavior of the system near bifurcation point in the presence of noise. We fou nd that the bifurcation point widens to a shifts as the level of noise increases. The form and content of this abstract are approved. I recommend its publication. Approved: Masoud Asadi Zeydabadi PAGE 4 TABLEOFCONTENTS IINTRODUCTION...............................1 IIESSENTIALIDEASFROMNONLINEARDYNAMICS..........7 II.ILinearSystems...............................7 II.IISimpleHarmonicMotionwithDamping.................9 II.IIILimitCycleandHopfBifurcation.....................11 II.IVANonlinearSystem:theClassicVanderPolModel...........11 II.VModiedVanderPolModel........................18 II.VIRayleigh-vanderPolModel........................21 IIIESSENTIALIDEASFROMSTOCHASTICDIFFERENTIALEQUATION SDE......................................27 III.IClassicalExampleofStochasticProcess:BrownianMotion.......27 III.IIGeneratingaStochasticProcessNumerically...............27 III.IIIEmbeddingStochasticProcessesinDierentialEquations........28 III.IVNumericalMethodsforSDE........................28 III.IV.IEuler{Maruyamamethod........................29 III.VSolvingStochasticIntegrals........................29 III.V.IIt^o'sRule.................................29 III.V.IIIt^o'sChainRule.............................30 III.V.IIIAnalyticalexamples...........................30 IVSTOCHASTICDYNAMICS..........................37 IV.IStochasticVanderPolModel.......................37 IV.IIStochasticRayleigh-vanderPol......................37 VBEHAVIORNEARTHEBIFURCATION..................42 V.INumericalSimulation............................42 V.I.IPoincareSection.............................42 v PAGE 5 V.I.IIBehaviorofMeanandVariancewithBifurcationParameterforFixed NoiseLevel..........................46 VICONCLUSIONANDFUTUREWORK...................60 VI.IFutureWork.................................60 VI.I.ICoupledSystem.............................60 REFERENCES ...................................62 APPENDIX A.ANALYSISANDSIMULATIONFORCOUPLEDSYSTEM........64 vi PAGE 6 1 dX 1 CHAPTER I INTRODUCTION Bifurcation, a sudden change of the dynamical system, lies in the center of nonlinear dynamics[29][19]. It has plenty of applications in daily activity and fields of science such as the quick change of medical condition from stable to severe, financial mar ket change from balanced to crisis[19], or neural network change from stationary to generates motion. The needs to do a proper prediction of such phenomena rise as time goes by. However, reality has uncertainties so it is important to understand how ran domness effects the nonlinear systems. The general theory of random excitation in physical systems was first studied by Albert Einstein in 1905 and was generalized by Marian Smoluchowski in 1916[8]. Later on, physicists and mathematicians Andrey Kolmogorov, Adriaan Fokker, Max Planck, Leonard Ornstein, George Eugene Uhlenbeck, Norbert Wiener, Joseph L. D o ob and Ki y oshi I t o laid the foundation for solving linear SDE b e t w een 1930s and 1950s[8]. However, those earlier works are for solving linear SDE. Difficulties arise for nonlinear system because of the complexity nature of nonlinear systems[8]. Back in 1957, Garstens started studying noise in nonlinear oscillators[12]. Garstens used a nonlinear Van der Pol equation with driving noise term to model nonlinear s in such model are multiplicative (i.e. X 1 Â· X 2 , where dt = X 2 ). Garstens presented a method PAGE 7 2 for estimating the effect of nonlinear noise in an oscillator at low levels of oscillation. This required an approximate solution of a non homogeneous van der Pol type equa tion with a noisy driven force by reducing the equation to a linear form. The powe r spectrum of the output was used to study the effect of noise in the nonlinear system. Caugheys paper in 1959[6] showed that there was interest in the how nonlinear systems response to random excitation, especially in the systems with multiplicative nonli near terms. He studied the self excited system like Van der Pol Oscillator and its response to external stochastic terms acting as external random excitation. Caughey found that the solution has a stochastic/noise term because of the stochastic term i n the self excitation system, and that the strength of the stochasticity decreases as the amplitude of oscillation increase. Interestingly enough, so does the band width of the stochasticity/noise which prove the early work done by Garstens incorrect. In C augheys paper in 1967[7], the use of a Fokker Planck equation had been introduced to study the stationary response of self excited oscillators. Using such technique, Caughey obtained the exact solution for the stationary joint probability density function, as well as the spectral density for oscillator with weak damping and weak excitation. Following his previous work, Caughey developed a nonlinear theory of random cause of t he developments in fields like high speed flight and the needs of dealing the random fluctuation cause by turbulence back at that time. Caughey summarized and developed some solution techniques for random vibrations in nonlinear systems[8]. During the simi lar time in Japan (1974), Kawakubo and Kabashima also used the Fokker Planck approach to study stochastic processes in self excited oscillation using a Van der Pol stochastic nonlinear system in a more experimental way[17]. Kawakubo and Kabashima focused on how the fluctuations, namely the amplitude and phase PAGE 8 fluctuations in the system, behave in the experiment of Wien Bridge oscillators. They showed how the voltage fluctuates with the feedback factor and how the variance fluctuates with the parameter controlling the nonlinearity. Kawakubo and Kabashima also included the probability distribution function obtained with the Fokker Planck equation, examining how it changed with the parameter controlling the nonlinearity. In agreement with Caughey and others find. Spanos showed that in 1978 such problems were drawing more and more attention[28]. This included several more applications of stochastic phenomena i n physics and engineering, such as control, aerospace, natural disaster, ocean wave and transportation engineering. Spanos pointed out that even though linear random vibration has be well studied, the theory of non linear random vibrations was not. This is still the case: there are still a lot of development needed in this topic. The Fokker Planck approach was used in Spanos work to obtain solutions, Spanos also showed the application of Fokker Planck approach in Van der Pol and Rayleigh oscillators. Later on in 1995, Leung did a systematic study of the behavior of Van der Pol oscillator with additive and multiplicative stochastic terms using different nonlinear ity parameters[20], comparing the bifurcation behaviors of the system with di fferent stochastic terms. By such methods Leung stated that we could get a more realistic understanding of the bifurcation of a nonlinear system. Leung found that all level of stochastic terms / noise chosen disturb the usual behavior of the Van der Pol o scillator. Moreover, he found that noise could postpone the bifurcation between fixed point and limit cycle and that noise could induce a transition from a stable fixed point to a stable limit cycle, and vice versa. Leung suggested that noisy nonlinear sys tems could have more and different phenomena than the deterministic one, such as stochastic resonance. 3 PAGE 9 In later work in 1998[21], he pointed out studies on stochastic bifurcation focus on the shift of bifurcation point, where a s hifting of bifurcation point is an indication of the noise effect on the two competing mechanisms, namely self sustained oscillation and nonlinear dissipation. Yuan and others made further developments in this direction in 2013 and found that zero mean multiplicative noise can induce the bifurcation points and sometimes even change the topology of the system[31], or in another words, the structure does not correspond to the deterministic counterpart, which make the analysis of such system more difficult. They also used the steady state distribution of the process to approach such case, in where correspond means the probability density function does not change over time and have a maximum at the stable attractors, which later on define as t he steady state distribution multiplying negative one being a Lyapunov function of the deterministic counterpart of the stochastic dynamics. Yuans group apply a new stochastic interpretation rather than the common Ito and Stratonovichs integration to ob tain a correspondence between stochastic and deterministic dynamics. Just last year (2018), Doan and others were trying to merge bifurcation theory for stochastic dynamical systems[10]. Doan and others work shows that the earlier concepts developed in t he 1990s, namely theory of phenomenological (P) bifurcation and dynamical (D) bifurcation, cannot comprehensively capture the complexity of bifurcation in random dynamical system. Doan points out that despite its relevance for many applications of topical interest, a bifurcation theory of random dynamical systems is still in its infancy[10]. With the developments made in mid 20th century on time dependent probability theory[2][16][26] and stochastic differential equations[25], and we began to have tools to solve the stochastic system in the form of stochastic differential equation[13][25]. Yet The bifurcation of stochastic chaotic system is still something we do not suffi ciently understand yet. But the importance of bifurcation in such system is becoming 4 PAGE 10 5 dX 1 increasingly unavoidable[1][10]. In recent years, neuroscience (as well as machine learning and artificial intelligent) draws more and more attention[22] and funding, so do the nonlinear and stochastic components of this field[1][22][15]. There are some breakthroughs on the applications of such topics in crucial areas such as neuroscience[14][11][17]. However, most of the developments about stochastic differential equations are in the field of stock market forecasting and related topics in quantitative/financial engineering[25], i.e. the early detection of economic bubbles[9]. How to analyze the stochastic nonlinear dyn amical model is still a difficult task which been rarely discussed in the context of stochastic differential equation but mostly done by deterministic approximation, especially for nonlinear system with multiplicative terms in neural networks (i.e. System contains X 1 Â· X 2 , where dt = X 2 ). Most of the analytical and numerical methods cannot be applied to such terms[5][18][25]. Yet the interest in self excited oscillation and its application in neural network pushes us forward to study such system[3][17], especially the response of a class of oscillators to stochastic excitation[28]. Now we might have tools and methods to study it in the content of stochastic di fferential equation instead of deterministic approximation methods. In this thesis we are building a modified Rayleigh Van der Pol model with additive stochastic components. We hope to test if the system behaves as others observed, when we put the system in the proper context of stochastic differential equation in stead of deterministic methods. We also hope to lay a foundation for future work and potential applications of the Rayleigh Van der Pol model. The applications include neural networks for the Rayleigh Van der Pol[6][7] based central pattern generator (CPG) model[3], as well as other potential applications i n physics, neuroscience, fi nance and other related fields. This work is also offered as a guide for some methods of how to simulate and analyze models of nonlinear system with stochastic components in the context of PAGE 11 6 dX 1 stochastic differential eq uations. This is primarily focused on how to solve the system with terms containing random variable products X 1 Â· X 2 , where dt = X 2 , which we do not have many good options to solve yet. PAGE 12 7 CHAPTER II ESSENTIAL IDEAS FROM NONLINEAR DYNAMICS This chapter organizes a conceptual framework for analyzing two dimensional non linear dynamics. This provides a broad and powerful foundation for analyzing many applications especially after stochastic forcing is added in later chapters. II.I Linear Systems In general a two dimensional dynamical linear system has the form of x 1 = a 11 x 1 + a 12 x 2 (II.1a) x 2 = a 21 x 1 + a 22 x 2 (II.1b) where x 1 = dx , x 2 = dx 2 , and a 11 , a 12 , a 13 , a 14 are constant parameters. dt dt We can write the system in matrix form as x = A x (II.2) where A = a 11 a 1 2 and x = x 1 . a 21 a 22 x 2 We seek solution of the form [29], PAGE 13 8 x = e v. (II.3) Substituting (2.3) into (2.2), we get Av = (II.4) We say the solution in (2.3) is an eigensolution if v is an eigenvector of A with corresponding eigenvalue [29]. The characteristic equation of matrix A is given by det a 11 a 12 = 0 , (II.5) Equation (2.5) yields a 21 a 22 where 2 0 (II.6) = trace ( A ) = a 11 + a 22 (II.7a) = det ( A ) = a 11 a 22 a 12 a 21 (II.7b) Solving equation (2.6), we get 1 = 2 = + 2 , (II.8a) 2 2 . (II.8b) 2 where 1 and 2 are the eigenvalues of matrix A . We find the eigenvectors by substituting eigenvalues into one of the equations in (2.4). We get two linearly inde pendent eigenvectors v 1 and v 2 . Now we can get the general solution for x , x = c 1 e 1 t v 1 + c 2 e 2 t v 2 . (II.9) PAGE 14 9 0 m m II.II Simple Harmonic Motion with Damping One good example for the linear system is simple harmonic motion with damping. The F in this case is the sum of the spring and damping forces: dx F = kx b dt , (II.10) where k and b are the spring and damping for (2.10), we get d 2 x dx m = kx b (II.11) dt 2 dt We can rewrite equation (2.11) as d 2 x dx m + b + kx = 0 (II.12) dt 2 dt or x Â¨ + x + 2 x = 0 ( I I.13) where = b 0 = l k . The solution (2.12) and (2.13) has the form x = Ae (II.14) Substituting (2.14) into (2.12), then we have ( 2 + + k ) Ae = 0 (II.15) Equation (2.15) for the nontrivial solution of (2.12) implies 2 + b + k = 0 . (II.16) m m Equation (2.16) is the characteristic equation of (2.12). PAGE 15 10 4 m 2 m = 2 m Â± i m 4 m 2 (II.17) 4 m 2 m x = A 1 e 2 m te m 4 m 2 + A 2 e 2 m m 4 m 2 The roots of equation (2.16), eigenvalues of (2.12), for underdamped or oscillatory case, ( b 2 < k ), is given by b I k b 2 The general solution of (2.12) or (2.13) for weak damping, b 2 < k , is given by b i I k b 2 t b te i I k b 2 t or x = A e t e i w 2 2 / 4 t + A e t e i w 2 2 / 4 t ( I I.19) 1 0 2 0 For real valued solutions, A 1 and A 2 will be complex conjugates. We call x = 0 , x = 0 a fixed point solution. For 0 all solutions tend to this fixed point solution and it is therefore said to be stable, for 0 the fixed point is unstable. Figure (2.1) shows time series and phase plot for damped harmonic motion. Figure II.1: Underdamped Harmonic Motion (II.18) PAGE 16 11 dt Âµ 2 This is the standard, or classical form of the van der Pol equation. However, the Âµ II.III Limit Cycle and Hopf Bifurcation A limit cycle occurs when the phase plane ( x v.s. dx ) contains an isolated closed trajectory. Isolated means the neighboring trajectories either spiral towards or away from the limit cycle.[29] A Hopf bifurcation h appens where a fixed point of a system changes to a limit cycle. II.IV A Nonlinear System: the Classic Van der Pol Model The Van der Pol equation has been used as a model for neuron networks with feedback controls. It is named after Dutch physicist Balthasar va n der Pol[29][30], who applied the equation to electronic circuit oscillators. d 2 x 2 dx 2 dt 2 ( Âµ ) 0 dt + 0 x = (II.20) Provided that Âµ / = 0 we can write: d 2 x 2 dx 2 dt 2 Âµ (1 Âµ x ) 0 dt + 0 x = 0 (II.21) Using the rescaled variables, t 0 t and x l x , equation (2.21) becomes: d 2 x d t 2 dx Âµ (1 x ) d t + x = 0 (II.22) scale factor l becomes singular at Âµ = 0 at the same point that the actual output, x , has an amplitude that is just about to rise from a value of zero[30]. This is one of the reasons that we modified the classical Van der Pol model as shown in the later section.(see section 2.3 and 2.4) An alternate dimensionless form of the dynamics is the single scaled nonlinear PAGE 17 12 equation in x : PAGE 18 13 x 1 l dt l d 2 x 2 dx dt 2 ( Âµ x ) dt + x = 0 . (II.23) Here, the scaled time and dynamic variable are t 0 t and x x where x . For notational clarity, we the equations are still in the scaled form. We can applying method of averaging to the Van der Pol equation to find ap proximate solutions to the nonlinear system. The dynamics, equation (2.23), can be re written as two first order differential equations dx = w (II.24a) dt dw 2 dt = x + ( Âµ x ) w. (II.24b) We are going to assume that the the solution will evolve from the simple sinsoidal solution when the nonlinear term is turned on. Thus we think of time variations in amplitude and phase around the limit cycle such that the phase space trajectory can be written in the form: x ( t ) = A ( t ) cos[ t + ( t )] (II.25a) w ( t ) = A ( t ) sin[ t + ( t )] (II.25b) Our goal is to find equations for the amplitude function A ( t ) and the phase function ( t ). In order to do so, we insert (2.25a) and (2.25b) for x and w into the dynamical equations (2.24a) and (2.24b) to obtain dA ( t ) ( t ) cos[ t + ( t )] A ( t ) 1 + sin[ t + ( t )] = A ( t ) sin[ t + ( t )] (II.26) dA ( t ) ( t ) sin[ t + ( t )] A ( t ) 1 + cos[ t + ( t )] = dt dt A ( t ) cos[ t + ( t )] Âµ A 2 ( t ) cos 2 [ t + ( t )] 1 A ( t ) sin[ t + ( t )] . (II.27) dt PAGE 19 14 2 2 2 2 2 2 2 1 2 2 1 2 dA ( t ) ( t ) dt dt 2 2 2 2 3 This simplifies to dA ( t ) cos[ t + ( t )] A ( t ) dt ( t ) sin[ t + ( t )] = 0 (II.28) dt dA ( t ) sin[ t + ( t )] A ( t ) dt ( t ) cos[ t + ( t )] = dt Âµ A ( t ) cos [ t + ( t )] 1 A ( t ) sin[ t + ( t )] . (II.29) Use (2.28) to write ( t ) dt dA ( t ) cos[ t + ( t )] = dt . (II.30) A ( t ) sin[ t + ( t )] Substitute (2.30) into (2.29) and multiply through by sin[ t + ( t )]: dA ( t ) 2 dA ( t ) 2 sin [ t + ( t )] dt cos [ t + ( t )] = dt E A ( t ) cos [ t + ( t )] 1 A ( t ) sin [ t + ( t )] . (II.31) Collect terms and use the fact that cos 2 [ t + ( t )] + sin 2 [ t + ( t )] = 1 to obtain dA ( t ) dt = Âµ A ( t ) cos [ t + ( t )] A ( t ) sin [ t + ( t )] . (II.32) We can also get ( t ) dt = Âµ A ( t ) cos [ t + ( t )] cos[ t + ( t )]sin[ t + ( t )] . (II.33) So far this has been an exact change of coordinates. We now make an approxi mation. Anticipating that A 2 ( t ) O ( Âµ ), we take both time derivatives and to likewise be O ( Âµ ) and thus the amplitude and phase are slowly varying deviations from the circular lim it cycle. We will then treat them as effectively constant over one cycle and replace the right hand side with averages. Taking t + ( t ), dA ( t ) 1 2 dt 2 Âµ A ( t ) cos ( t ) sin (II.34a) 1 2 1 2 0 0 0 2 2 PAGE 20 15 2 8 ÂµA ( t ) 2 sin A ( t ) 2 cos sin (II.34b) 1 1 3 = ÂµA ( t ) A ( t ) . (II.34c) PAGE 21 16 1 So we have dA ( t ) 1 2 dt 8 4 Âµ A ( t ) A ( t ) . (II.35) This is an example of an *amplitude equation*. It may be formally integrated starting with to arrive at dA 8 = dt (II.36) [4 Âµ A 2 ] A A + 2 Âµ Âµt Âµ = e . (II.37) At long times A 2 (2 Âµ A ) A 2 Âµ (II.38) (justifying thies assumption that A 2 O ( Âµ )), while the next order correction is Âµ ( 2 e Âµt ) . (II.39) Returning to dimensioned coordinates, x ( 2 e 0 t ) x Âµ cos( 0 t ) , (II.40) and the predicted final steady oscillation is given by x 2 x Âµ cos( 0 t ) . (II.41) where we recall that x . For future purpose, we re write the classcial Van der Pol system as two first order equations, dx = y (II.42) dt PAGE 22 17 dy 2 dt = Âµ Â· (1 x ) Â· y x (II.43) As we mentioned the Hopf bifurcation happens where a fixed point changes to a limited cycle. The fixed point of classic Van der Pol model is ( x , y ) = (0 , 0). In this system, the Hopf bifurcation happens at Âµ = 0, when the system from steady state figure (2.3) turn into a self sustain system (2.2). Figure II.2: Classical Van der Pol model, Stable Fixed Point, X versus t for Âµ = 0.1. PAGE 23 18 Figure II.3: Classical Van der Pol model, Stable Fixed Poi nt, Y versus t for Âµ = 0.1. Figure II.4: Classical Van der Pol model, Stable Fixed Point, Phase Portrait for Âµ = 0.1. We can see from both the figures of X and Y versus t that the model started oscillating and damping to zero after some time. The phase portrait shows the red dot, which is our initial point, spiral in and become a fixed point, that indicates the model reaches a fixed point. PAGE 24 19 Figure II.5: Classical Van der Pol model, Limit Cycle, X versus t for Âµ = 0.1. Figure II.6: Classical Van der Pol model, Limit Cycle, Y versus t for Âµ = 0.1. PAGE 25 20 Figure II.7: Classical Van der Pol model, Limit Cycle, Phase Portrait for Âµ = 0.1. We can see from both the figures of X and Y versus t that the model oscillates at an fixed amplitude after some time. The phase portrait shows the red dot, which is our initial point, spiral out and become a circular trajectory, that also indicate the model reaches a limit cycle. II.V Modified Van der Pol Model The classical Van der Pol is a good simple model for the current goal. However, during the numerical simulation, the simulation breaks when Âµ is close to zero, which agree with the discussion in section 2.2. Therefore, we are using one version of modified Van der Pol model, dx = y (II.44) dt dy 2 dt = ( Âµ x ) Â· y x (II.45) Comparing the second equation in classical with modified Van der Pol equations, (comparing equations (2.43) with and (2.45)) we find Âµ is outside the parenthesis in (2.43) but is inside in (2.45).That makes the nonlinear term more smooth (with no singularity) and more close to real world phenomena. For this reason it i s easier to PAGE 26 21 more consistent and does not break during simulation when Âµ is close to zero. For modified Van der Pol Model, the bifurcation still happens at Âµ = 0, as figures 2.8 2.13 indicate. Figure II.8: Modified Van der Pol model, Stable Fixed Point, X versus t for Âµ = 0.1. Figure II.9: Modified Van der Pol model, Stable Fixed Point, Y versus t for Âµ = 0.1. PAGE 27 22 Figure II.10: Modified Van der Pol model, Stable Fixed Point, Phase Portrait for Âµ = 0.1. Given the same initial condition and parameters, we can see that the modified Van der Pol behaves almost the same with the classical counterpart. However, the modified Van der Pol model reaches fixed point faster than the classical model because the nonlinear term Âµ is less destructive. Figure II.11: Modified Van der Pol model, Limit Cycle, X versus t for Âµ = 0.1. PAGE 28 23 Figure II.12: Modified Van der Pol model, Limit Cycle, Y versus t for Âµ = 0.1. Figure II.13: Modified Van der Pol model, Limit Cycle, Phase Portrait for Âµ = 0.1. Given the same initial condition and parameters, we can see as well that the mod ified Van der Pol behaves almost the same with the classical counterpart. However, the modified Van der Pol model reaches limit cycle faster than the classical model because the nonlinear term Âµ is less destructive. II.VI Rayleigh van der Pol Model We use (Modified) Rayleigh van der Pol equations as our standard system. One of the characteristics worth noting about Rayleigh van der Pol equations is that its PAGE 29 24 / phase portrait is circular, even for larger value of Âµ which make our future work on the phase plane a lot easier. dx = y (II.46) dt dy 2 2 dt = ( Âµ y x ) Â· y x (II.47) Rayleigh van der Pol is also a good representation of the Wein bridge circuit which can be conduct into experiment for testing [30]. Now we change the system to polar coordinate, x = r cos (II.48a) y = r sin (II.48b) r = x 2 + y 2 (II.49a) y = arctan x (II.49b) PAGE 30 25 / dt dt / r r 2 2 2 2 dr 1 dx dy = (2 x + 2 y ) (II.50a) dt 2 x 2 + y 2 = x 2 + y 2 ( x dt + y dt ) (II.50b) 1 2 2 = r xy + y [( Âµ y x ) y x ] (II.50c) 1 = ( xy + Âµy 2 y 4 x 2 y 2 xy ) (II.50d) 1 = ( Âµy 2 y 4 x 2 y 2 ) (II.50e) y 2 = ( Âµ x 2 r y ) (II.50f) r 2 sin 2 = [ Âµ ( x 2 r + y 2 )] (II.50g) = ( Âµ r ) r sin (II.50h) 1 dy dx = ( x y ) (II.51a) dt r 2 dt dt 1 2 2 2 = r 2 x [( Âµ y x ) y x ] y (II.51b) 1 2 2 2 2 = r 2 [( Âµ y x ) xy ( x + y )] (II.51c) 1 2 2 2 = r 2 [( Âµ r ) r sin cos r ] (II.51d) = ( Âµ r ) sin cos 1 (II.51e) The constant amplitude solution is r = Âµ (II.52a) = 1 . (II.52b) It predicts a pure sinusoidal oscillation at all values of Âµ , that is, one that contains no harmonics unlike the van der Pol model. Moreover the limit cycle is a circle of fixed radius proportional to Âµ in ( x, y ) phase space. Figures figures 2.6 and 2.7 1 dx dy PAGE 31 26 show the bifurcation at Âµ = 0. PAGE 32 27 Figure II.14: Rayleigh van der Pol model, Stable Fixed Point, X versus t for Âµ = 0.1. Figure II.15: Rayleigh van der Pol model, Stable Fixed Point, Y versus t for Âµ = 0.1. PAGE 33 28 Figure II.16: Rayleigh van der Pol model, Stable Fixed Point, Phase Portrait for Âµ = 0.1. Given the same initial condition and parameters, the figures shows that the Rayleigh van der Pol model reaches the fixed point as fast as the modified Van der Po l model with a different oscillation trajectory and a more circular phase portrait even for the case of fixed point. Figure II.17: Rayleigh van der Pol model, Limit Cycle, X versus t for Âµ = 0.1. PAGE 34 29 Figure II.18: Rayleigh van der Pol model, Limit Cycle, Y versus t for Âµ = 0.1. Figure II.19: Rayleigh van der Pol model, Limit Cycle, Phase Portrait for Âµ = 0.1. Given the same initial condition and parameters, the figures shows that the Rayleigh van der Pol model reaches the limit cycle as fast as the modified Van der Pol model with a different oscillation trajectory and a more circular phase portrait. The black circle in the phase portrait indicate the constant amplitude solution where r = Âµ . PAGE 35 30 CHAPTER III ESSENTIAL IDEAS FROM STOCHASTIC DIFFERENTIAL EQUATION (SDE) III.I Classical Example of Stochastic Process: Brownian Motion A Brownian motion B t (Wiener process W t ), is a sequence of random variables satisfies the following conditions [25]: 1) B t (0) = 0 with probability 1, 2) for 0 s < t T the random variable given by difference B t B s is normally distributed with mean 0 and variance t s , or B t B s t sN (0 , 1), where N (0 , 1) is normal distribution with zero mean and unit variance. 3) in 0 s < t u < v T , the increments B t B s and B v B u are independent. III.II Generating a Stochastic Process Numerically For the purpose of the numerical simulation, we discretize Brownian Motion as dW t tN (0 , 1). A Sample path of Brownian motion is shown in figure 3.1, it is a random path start at 0 with mean 0 and variance dt. PAGE 36 31 Figure III.1: sample path of Brownian Motion With the help of measure theory, it has been proved that Brownian Motion like showed in figure 3.1 is nowhere differentiable [2], we can consider it as a special case of fractal to understood its non differen tiablity. III.III Embedding Stochastic Processes in Differential Equations A stochastic process [25] has a general form of dX t = f ( X t , t ) Â· dt + g ( X t , t ) Â· dB t (III.1) where B t stand for Brownian motion/Wiener process( W t ), f ( X t , t ) dt is the drift term which carry the deterministic characteristics and g ( X t , t ) B t is the diffusion term which carry the stochastic characteristics. It is worthy noted that because Brownian motion is nowhere differentiable [2], which lead to a result that stochastic process is nowhere differentiable as well, there fore we need write it in difference equation form before we try to solve it analytically and numerically. III.IV Numerical Methods for SDE III.IV.I Euler Maruyama method PAGE 37 32 t A stochastic differential equation in general in the integral form is given by X ( t ) = X 0 + t f ( X ( s )) ds + 0 t g ( X ( s )) dW ( s ) , 0 t T In general w e use I t o v ersion of i n tegral since most of the meth o ds in st o c hastic di ff ere n tial equation are b ased on I t o c hain rule n u merically and analyticall y . integral solution in differential form: dX ( t ) = f ( X ( t )) dt + g ( X ( t )) dW ( t ) , X (0) = X 0 , 0 t T To solve this equation numerically, we t = T/N where N is the total steps and X i = X ( t i ). The Euler Maruyama (EM) method is similar to method for ODE with additional stochastic term. X i = X i 1 + f ( X i 1 t + g ( X i 1 W, i = 1 , 2 ...N. For the case of additive noise, we write the diffusion terms W tN (0 , 1) (III.2a) = (0 , t ) (III.2b) where is a constant parameter that influence the strength of noise. The level of stochasticity is controlled by this parameter. III.V Solving Stochastic Integrals III.V.I I t o Rule I t rule pr o v e that for st o c hastic pr o cess dX t = f ( X t , t ) dt + g ( X t , t ) dB t , dt 2 = 0 , dtdB t = 0 (III.3a) dB t dt = 0 , dB 2 = dt (III.3b) This is the one of the major breakthrough leading the method to solve the stochas tic integrals. 0 PAGE 38 30 i i i i t t j 2 x x j t j k for i = 1,2. Note that g 2 t t 2 t III.V.II I t o Chain Rule III.V.II.I Single V ariable I t o Chain Rule I t Chain Rule s a ys, for Y t = g ( t, X t ), where X t is st o c hastic pr o cess, w e h a v e 1 2 g 2 dY t = ( t, X t ) dt + ( t, X t ) dX t + ( t, X t )( dX t ) 2 2 (III.4) III.V.II.II Multi v ariate I t o Chain Rule for Y i = g i ( t, X t ), where X t is stochastic process, we have 2 2 2 2 dY = ( t, X ) dt + ( t, X ) dX 1 g + ( t, X ) dX dX (III.5) 2 i x j x k = 0 j, k . III.V.III Analytical examples III.V.III.I Simple Brownian Motion We can see the main difference between Stochastic differential equation (SDE) and ordinary differential equation (ODE) by integrating over Brownian Motion B t : t B s dB s . (III.6) 0 first we guess the form of the solution from deterministic model, in this case, it is y = g ( t, x ) = 1 x 2 , then we try the form for Y t . 1 Y = g ( t, X ) = B 2 (III.7) t t 2 t Using It c hain rule w e h a v e dY = 0 Â· dt + B Â· dB 1 2 + Â· 1 Â· ( dB ) (III.8) 1 d ( B 2 ) = B dB 1 + dt (III.9) 2 t where ( dB t ) 2 = dt b y It Rule, k j =1 k =1 j j =1 k j =1 k =1 j j =1 t PAGE 39 31 t t 2 PAGE 40 32 t t 1 B 2 = 2 t 1 B s dB s + 0 2 t 1 ds. (III.10) 0 Which give us the solution, t B s dB s = 0 1 1 B 2 t (III.11) 2 2 III.V.III.II RLC circuit For a simple RLC circuit, we have system where F t = G t + t , 1 LQ t + RQ t + C Q t = F t (III.12) W t : White noise ( W t dt = dB t ), B t : one dimensional Brownian Motion, Q t : the change at time t at a fixed point in an electric circuit, F t : the potential source at time t, G t : the potential source at time t without noise, R : the inductance, L : the resistance, C : the resistance. We rename the variable as following: X = X = X 1 = Q t (III.13) We write equation (3.12) as a system of the first order SED: X 1 = X 2 (III.14a) 1 LX 2 = RX 2 C X 1 + G t + t (III.14b) Now we write the system of SED, equations (3.14a) and (.14b), in the from of difference equation as following Q t X 2 Q t X 2 t PAGE 41 33 L 1 2 x 2 dX 1 = X 2 dt (III.15a) R 1 1 dX 2 = L X 2 dt CL X 1 dt + L G t dt + L W t dt (III.15b) where W t dt = dB t . In order to solve this two dimensional SDE, we represent the system, (3.15a) and (3.15b), in the matrix form equation, dX = dX t = AXdt + Hdt + KdB t (III.16a) d X = dX 1 , A = 0 1 , H ( t ) = 0 , K = 0 (III.16b) dX 2 1 R 1 G t N o w solving the ma t rix equation (Note that e At = k =0 A K t K k ! converges.), we move AX t dt to the left side of equation (3.16a) and multiply both sides by e At . dX t AX t dt = H ( t ) dt + KdB t (III.17a) e At [ dX AXdt ] = e At [ H ( t ) dt + KdB t ] (III.17b) Now we want to show if we can write the left hand side of equation 3.14b as d ( e At X ). I n attempt to fi gure this out, w e try to apply m ulti v ariate It o c hain rule to g ( t, x , x ) = e At x 1 , dY = d ( e At X ) = Ae At Xdt + e At dX (III.18a) = e At [ dX AXdt ] . (III.18b) (Note that Ae At = e At A .) L CL L L PAGE 42 34 d ( e At X ) = e At [ H ( s ) dt + KdB t ] . (III.19) Integrate both sides, e At X X 0 = t e As H ( s ) ds + 0 t e As KdB s (III.20) 0 Recall integration by parts we have t f ( s ) dB s = f ( t ) B t 0 t d s ( f ( s )) B s (III.21) 0 Use integration by parts for the last term on the right hand side in equation (3.20), t e As KdB s = e At KB t 0 Substituting (3.22) into (3.20), we get t ( A ) e As KB s ds (III.22) 0 e At X t = X 0 + t e As H ( s ) ds + e At KB t + 0 t Ae As KB s ds (III.23) 0 Which give us the solution X t = e At X 0 + KB t + e At t e As [ H ( s ) + AKB s ] ds (III.24) 0 III.V.III.III Damped Harmonic Oscillator The equation of motion for Damped Harmonic Oscillator with Brownian Noise is given by PAGE 43 35 dx B t 2 X t kx b dt + a dt = m 2 (III.25) B t : Brownian noise ( W t dt = dB t where W t is white noise), k : spring constant, b : damping constant, a : constant, m : mass. Now let us to write it in the content of stochastic differential equation. 2 X t m 2 + b dX t + kX = a dt t dB t (III.26) dt Equation (3.26) can be written as a system of two first order SDE as following X 1 = X 2 (III.27a) mX 2 + bX 2 + kX 1 = aB t (III.27b) Writing the system of SDE in the from of difference equation, we have dX 1 = X 2 dt (III.28a) b k a dX 2 = m X 2 dt m X 1 dt + m dB t (III.28b) Now we rewrite the system in the form of matrix equation, X = X ( t ) = X 1 (III.29) X 2 X 2 PAGE 44 36 dX 2 k m b m a m dX = AXdt + CdB t (III.30a) d X = dX 1 , A = 0 1 , c = 0 (III.30b) get Moving AXdt to the left side of (3.30a) and multiplying both sides by e At , we dX AXdt = CdB t (III.31a) e At [ dX AXdt ] = e At CdB t (III.31b) We already know from the previous example that we can write the left hand side of equation 3.31b as d ( e At X ). Which give us d ( e At X ) = e At CdB t (III.32) Integrate both sides, e At X t X 0 = t e As CdB s (III.33) 0 Use integrate by parts for the last term in the right hand side in equation 3.33, t e As CdB s = e At CB t 0 Bring the term in we have e At X t = X 0 + e At CB t Which give us the solution t ( A ) e As CB s d s (III.34) 0 t Ae As CB s d s (III.35) 0 PAGE 45 37 4 m 2 m X t = e At X 0 + CB t + e At t Ae As CB s ds (III.36) 0 It is worth noting that the damped harmonic oscillator and RLC circuit in the previous section have similar form of solution. In the RLC circuit, in addition to the noise, the external source is also included in the system. The stochastic RLC in previous section mathematically is similar to the stochastic damped driven harmonic oscillator. Figure 3.2 show the result for the case of damped harmonic motion with weak noise ( = 0 . 01) for weak damping ( b 2 < k ), we can see that for small , the stochastic term did not influence the damping that such, but started to take over after. The phase portrait indicate that the stochastic term keep kicking the system out of the fix point. Figure III.2: 1, m = 5, b = 0.01, k = 1. PAGE 46 38 2 2 2 CHAPTER IV STOCHASTIC DYNAMICS IV.I Stochastic Van der Pol Model In general, the system always have some level of noise. We take the external additive white noise into consideration and treated it as a forcing term. The stochastic Van der Pol model with additive noise has the form dX = Y Â· dt + dW 1 ( t ) (IV.1a) dY = ( Âµ X ) Â· Y Â· dt X Â· dt + dW 2 ( t ) (IV.1b) IV.II Stochastic Rayleigh van der Pol Stochastic Rayleigh van der Pol equation is given by dX = Y Â· dt + dW 1 ( t ) (IV.2a) dY = ( Âµ Y X ) Â· Y Â· dt X Â· dt + dW 2 ( t ) (IV.2b) We found out that the analytical method for linear stochastic differential equa tions does not work for multiplicative variables. We cannot rewrite stochastic nonlin ear system into matrix equations when there are multiplicative variables. Moreover, we found no theory about analytical method for stochastic nonlinear systems with PAGE 47 39 2 2 multiplicative variables. However, we can solve such system using Euler Maruyama method. In the case of circuit, we set it up in a more realistic way such that dW 1 ( t ) = 0. dX = Y Â· dt (IV.3a) dY = ( Âµ Y X ) Â· Y Â· dt X Â· dt + dW ( t ) (IV.3b) In figures 4.1 to 4.6 the time sires and phase portraits are shown for Âµ = 0 . 1 and Âµ = 0 . 1, below and above deterministic bifurcation value Âµ = 0, for = 0 . 01. Figure IV.1: Rayleigh van der Pol model, Stable Fixed Point, X versus t for Âµ = 0.01. PAGE 48 40 Figure IV.2: Rayleigh van der Pol model, Stable Fixed Point, Y versus t for Âµ = Figure IV.3: Rayleigh van der Pol model, Stable Fixed Point, Phase Portrait for Âµ = Figure 4.1 to 4.3 show that in the case of Stable Fixed Point, the figures for X and Y versus t suggest that the oscillation damps away and the stochasticity takes over. The phase portrait is now a fuzzy circular trajectory: the trajectory fluctuate while reaching the stable fixed point. Right after the system reach the stable fixed point, the stochasticity/noise kick the system out of it and that cause the system to restart damping. PAGE 49 41 Figure IV.4: Stochastic Rayleigh van der Pol model, Limit Cycle, X versus t for Âµ = 0.01. Figure IV.5: Stochastic Rayleigh van der Pol model, Limit Cycle, Y versus t for Âµ = 0.01. PAGE 50 42 Figure IV.6: Stochastic Rayleigh van der Pol model, Limit Cycle, Phase Portrait for Figure 4.4 to 4.6 showed that in the case of limited cycle, the figures for X and Y versus t suggest that the model have a fuzzy limit cycle instead, the fuzzy trajectory and widened cycling region of the phase portrait agree with it. The black circle in the phase portrait indicate the constant amplitude solution where r = Âµ . PAGE 51 43 CHAPTER V BEHAVIOR NEAR THE BIFURCATION V.I Numerical Simulation In this chapter the numerical Euler Maruyama method is used to solved the stochastic Rayleigh van der Pol equation. V.I.I P oinca r Â´ e Section Study and statistical analysis of a stochastic system provides useful tools for un derstanding the behavior and propertie s of the system. In order to obtain and measure the statistical behavior of the nonlinear stochastic system, we chose to measure the statistical pro p erties of the P oinca r Â´ e section. Numerical si m ulation of P oinca r Â´ e section for deterministic system requires Jacob nian which need take the (partial) derivatives of the system. That bring the difficulty since stochastic differential equation is nowhere differentiable. In o r der to obtain the P oinca r Â´ e section for o ur system, w e set up a (partial) h y p er Every time the trajectory in phase plane cross the hyper plane, we mark it. This can give a good We only need to collect data for any one of the axis(i.e. the positive x axis) that the trajectory cross thanks to the circular phase portrait of Rayleigh van der Pol. PAGE 52 44 Figure V.1: Poin c a r Â´ e S e ction for Phase Port r ait of St o chastic R ayleigh van der Pol model under different noise level, Âµ = 0 . 1 = 0 . 01 , t > 1300 , for t [0 , 2000] . As figure 5.1 shows, under the influence of relatively weak noise ( = 0 . 01), the system starting became slightly unstable. The stochastic term causes the fluctuation of the limit cycle, which showed as the fuzzy path and the widened circular region of the phase portrait. As the strength (the level) of noise gets larger the fluctuation also increases. PAGE 53 45 Figure V.2: Poin c a r Â´ e S e ction for Phase Port r ait of St o chastic R ayleigh van der Pol model under different noise level, Âµ = 0 . 1 1300 , for t [0 , 2000] . Under the influence of relatively stronger noise ( = 0 . 05), figure 5.2 suggests that the system became way more unstable. We can also tell by the increasing of the width of P oinca r Â´ e section. PAGE 54 46 Figure V.3: Hist o g r am of Poin c a r Â´ e S e ction for Phase Port r ait of St o chastic R aylei gh van der Pol model under different noise level, Âµ = 0 . 1 = 0.01. W e obtained the histogram of P oinca r Â´ e section to obser v e the statistical pro p er t y of the system. It is a equi v ale n t to vi e w the cross section of th e P oinca r Â´ e section. Figure 5.3 showed that under the influence of relatively weak noise ( = 0 . 01), we can see an almost normal distribution. The red line indicate the constant amplitude solution where r = Âµ . PAGE 55 47 Figure V.4: Hist o g r am of Poin c a r Â´ e S e ction for Phase Port r ait of St o chastic R aylei gh van der Pol model under different noise level, Âµ = 0 . 1 = 0.05. Figure 5.4 showed that in the case of relatively stronger noise ( = 0 . 05), the distribution had increased the bell width and shift to the left while be less normal distributed. V.I.II Behavior of Mean and Variance with Bifurcation Parameter for Fixed Noise Level We found that the variance does not influence by dt for our system with mu / = 0 even though for a stochastic process [ V ar ] = c 2 Â· dt where c 2 = for our case. After the data collecting, w e computed the mean and v a r iance of the P oinca r Â´ e section versus Âµ to see if we can get any information for the bifurcation, as shown from figure II.12 to figure II.17. As we can tell from the figures, for small noise level (i.e. = 0 . 01), we can easily tell where the bifurcation happens (around Âµ = 0). From the sudden increase of mean in mean versus Âµ and the peak in variance versus Âµ , at PAGE 56 48 this noise level we can also observe the phase portrait for indications of bifurcation by changing Âµ to see the difference in phase portraits for limit cycle and fixed point. But for raised noise level (i.e. = 0 . 05), it is difficult to find out where the bifurcation happens by the phase portrait, but the change of mean in mean versus Âµ and the peak in variance versus Âµ . For intensified noise level ( 0 . 5), the mean versus Âµ cannot give us much indication about where the bifurcation happens, however the variance versus Âµ can still give us a region where the bifurcation happens. The theory of Hopf bifurcation defines where bifurcation happens as a critical point. However, after observing the experimental data, it seems like for the case of bifurcation under the influence o a higher amplitude while shifting to the right instead of stay at Âµ = 0 as the stochastic term increases. W e h a v e o btained the mean and v ariance of the P oinca r Â´ e section v ersus Âµ . It is worth noting that while we ran the simulation, we seeded the pseudo random number generator to ensure we were studying the same sequence of random variables. PAGE 57 49 (a) Mean v.s. Âµ for = 0.01. (b) Variance v.s. Âµ for = 0.01. Figure V.5: Stochastic Rayleigh van der Pol model, behavior of mean and variance, 1300 . PAGE 58 50 For = 0 . 01, we can tell from figure 5.5 where the bifurcation happened by the sudden increase in mean versus Âµ as well as the peak in variance versus Âµ . The peak of variance is around Âµ = 0 where the bifurcation happens for the deterministic model. PAGE 59 51 (a) Mean v.s. Âµ for = 0.02. (b) Variance v.s. Âµ for = 0.02. Figure V.6: Stochastic Rayleigh van der Pol model, behavior of mean and variance, 1300 . PAGE 60 52 For = 0 . 02, we are still be able to tell from figure 5.6 where the bifurcation happened by the sudden increase in mean versus Âµ as well as the peak in variance versus Âµ Âµ is a little more smooth than in the care of = 0 . 01. The variance also had a higher value at the peak and wider peak region. The peak of variance slightly shifts to the right side of Âµ = 0. PAGE 61 53 (a) Mean v.s. Âµ for = 0.05. (b) Variance v.s. Âµ for = 0.05. Figure V.7: Stochastic Rayleigh van der Pol model, behavior of mean and variance, 1300 . PAGE 62 54 In figure 5.7, for = 0 . 05, we are still be able to tell where the bifurcation happened by the sudden increase in mean versus Âµ as well as the peak in variance versus Âµ . However, Âµ is more smooth than in the care of = 0 . 2. The variance also had a higher value at the peak and wider peak region compare with the case for = 0 . 02. The peak of variance shifts to the right side of Âµ = 0 even more. PAGE 63 55 (a) Mean v.s. Âµ for = 0.1. (b) Variance v.s. Âµ for = 0.1. Figure V.8: Stochastic Rayleigh van der Pol model, behavior of mean and variance, 1300 . For = 0 . 1, as showed in figure 5.8, we had difficulty to tell where the bifurca PAGE 64 56 tion happened by mean versus Âµ , but we can tell the region where the bifurcation happened. We can still be able to tell where the bifurcation happened by the peak in variance versus Âµ . The variance also had a higher value at the peak and wider peak region compare with the previous cases. The peak of variance shifts to the right side of Âµ = 0 more significantly. PAGE 65 57 (a) Mean v.s. Âµ for = 0.7. (b) Variance v.s. Âµ for = 0.7. Figure V.9: Stochastic Rayleigh van der Pol model, behavior of mean and variance, 1300 . PAGE 66 58 Figure 5.9 shows mean and variance versus Âµ for significantly bigger noise/s tochastic term ( = 0 . 7). We cannot tell where bifurcation region is just by mean versus Âµ . The peak region of the variance, which indicates the bifurcation region, con tinued to widen with increased peak value compare with the previous cases. There seems to be a peak in the bifurcation region at Âµ = 0 . 7 as well. The surface plots in appendix A shows a better visualization for coupled system. PAGE 67 59 (a) Mean v.s. Âµ for = 0.8. (b) Variance v.s. Âµ for = 0.8. Figure V.10: Stochastic Rayleigh van der Pol model, behavior of mean and variance, 1300 . PAGE 68 60 For intensified noise/stochastic term ( = 0 . 8) showed in figure 5.10. We are still unable to tell where bifurcation region is just by mean versus Âµ . The bifurcation region can still be found by the peak region of the variance, which is continued to widen with increased peak value like in previous case. It is more clear to see the Âµ = 1 here. PAGE 69 61 CHAPTER VI CONCLUSION AND FUTURE WORK For the system of (modified) Rayleigh Van der Pol with additive stochastic term, the the right instead of staying at Âµ = 0 as the stochastic term increases. VI.I Future Work One of the applications of Rayleigh Van der Pol Oscillator is in the simulation of the Wein Bridge circuit[30]. VI.I.I Coupled System One of the applications of the coupled nonlinear system is in the simulation of the animal gaits[3][4][27]. VI.I.I.I Two Synaptic Coupled Stochastic Rayleigh van der Pol Model In the content of neuroscience, the neuron network can be represented in the form of coupled oscillators. There are different kind of coupling method such as synaptic coupling and diffusive coupling. We choose synaptic coupling since it is a more realisti c way to represents how different neurons connect each other. PAGE 70 62 2 2 1 1 2 2 2 2 System of Two Synaptic Coupled Stochastic Rayleigh van der Pol Equations: dX 1 = Y 1 Â· dt + 1 Â· X 2 Â· dt + dW 1 ( t ) , (VI.1a) dY 1 = ( Âµ Y X ) Â· Y 1 Â· dt X 1 Â· dt + 2 Â· Y 2 Â· dt + dW 2 ( t ) , (VI.1b) dX 2 = Y 2 Â· dt + 3 Â· X 1 Â· dt + dW 3 ( t ) , (VI.1c) dY 2 = ( Âµ Y X ) Â· Y 2 Â· dt X 2 Â· dt + 4 Â· Y 1 Â· dt + dW 4 ( t ) , (VI.1d) Where dW 1 ( t ) = dW 3 ( t ) = 0, i = i = 1 , 2 , 3 , 4. See appendix A for analysis and simulation for two coupled Rayleigh van der Pol system . PAGE 71 63 REFERENCES [1] Anishchenko Vadim S., Astakhov Vladimir, Neiman Alexander, Vadivasova Tatjana., Schimansky Geier Lutz, Nonlinear Dynamics of Chaotic and Stochastic Systems: Tutorial and Modern Developments , Springer, 2007. [2] Patrick Billingsley, Probability and Measure , Wiley Series in Probability and Statistics, 2012. [3] Pietro Luciano Buono, Martin Golubitsky , Models of central pattern generators for quadruped locomotion I. Primary gaits , J. Math. Biol, 2001. [4] Pietro Luciano Buono, Models of central pattern generators for quadruped locomotion II. Secondary gaits , J. Math. Biol, 2001. [5] Kevin Burrage, Ian Lenane, Grant Lythe, Numerical Methods for Second Order Stochastic Differential Equations , SIAM J. SCI. COMPUT Vol. 29, No. 1, pp. 245264, 2007. [6] T. K. Caughey, Response of Van , Journal of Applied Mechanics, 1959. [7] T. K. Caughey, On the Respones of a Class of Self excited Oscillators to Stochastic Excitation , Int. J. Non Linear Mechanics, 1967. [8] T. K. Caughey, Nonlinear Theory of Random Vibration , Advances in Applied Mechanics Volume 11, 1971. [9] Andrey Dmitrieva, Victor Dmitrieva, Oleg Sagaydakb, Olga Tsukanovaa, The Application of Stochastic Bifurcation Theory to the Early Detection of Economic Bubbles , Procedia Computer Science, 2017. [10] Thai Son Doan, Maximilian Engel, Jeroen S W Lamb, Martin Ras mussen, Hopf bifurcation with additive noise , Nonlinearity 31 4567, 2018. [11] Armin Fuchs, Nonlinear Dynamics in Complex Systems, Theory and Applica tions for the Life , Neuro and Natural Sciences , Springer, 2013. [12] M.A. Garstens, Noise in Non Linear Oscillators , J. Appl. Phys., 1957. [13] Kiyosi It o , St o chastic Int e g r al , Pr o c . Imp. Acad., 1 944. [14] Eugene M. Izhikevich, Dynamical Systems in Neuroscience: The Geometry of Excitability and Bursting , The MIT Press, 2007. PAGE 72 63 [15] Kurt Jacobs, Stochastic Processes for Physicists: Understanding Noisy Sys tems , Cambridge University Press, 2010 [16] Samuel Karlin, Howard M. Taylor, A First Course in Stochastic Pro cesses , Academic Press, 1975. [17] Tatsuyuki Kawakubo, Shigeharu Kabashima, Stochastic Process in Self Excited Oscillation , Journal of Physical Society Japan, 1974. [18] Peter Eris Kloeden, Eckhard Platen, Henri Schurz, Numerical Solu tion of SDE Through Computer Experiments , Springer, 2003. [19] Christian Kuehn, A Mathematical Framework for Critical Transitions: Bi furcations, Fast Slow Systems and Stochastic Dynamics , Physica D: Nonlinear Phenomena, Vol. 240, No. 12, pp. 1020 1035, 2011. [20] H.K. Leung, Stochastic Transient of a noisy Van der Pol Oscillator , Physica A, 1995. [21] H.K. Leung, Stochastic Hopf Bifurcation in a biased Van der Pol Model , Phys ica A, 1998. [22] Carlo Liang, Gabriel J. Lord, Stochastic Methods in Neuroscience , Oxford University Press, 2010. [23] Marc Mendler, Johannes Falk, Barbara D rosse , Analysis of stochastic bifurcations with phase portraits , arXiv, 2018. [24] John F. Moxnes, Kjell Hausken, Introducing Randomness into First Order and Second Order Deterministic Differential Equations , Advances in Mathemat ical Physics Volume 2010, 2010. [25] Bernt Oksendal, Stochastic Differential Equations (6th Ed.): An Introduction with Applications , Springer, 2010. [26] Sheldon M. Ross, Stochastic Processes , John Wiley and Sons, Inc., 1996. [27] Mozhdeh Saffari Parizi, Coupled Biological Oscillators and Animals Gaits , MIS Thesis, University of Colorado Denver, 2004. [28] P T. D. Spanos, Stochastic Analysis of Oscillators with Non Linear Damping , Int. J. Non Linear Mechanics, 1978. [29] Steven H. Strogatz, Nonlinear Dynamics and Chaos (2nd Ed.): With Appli cations to Physics, Biology, Chemistry, and Engineering , Westview Press, 2014. [30] Randall Tagg, UC Denver Physics Department, Private Communication. [31] Ruoshi Yuan, Xinan Wang, Yian Ma, Bo Yuan, Ping Ao , Explor ing a noisy van der Pol type oscillator with a stochastic approach , PHYSICAL REVIEW E 87, 062109, 2013. PAGE 73 64 APPENDIX A ANALYSIS AND SIMULATION FOR COUPLED SYSTEM This appendix contains the analysis and simulation for two coupled Rayleigh van der Pol system. For the coupled case, we are studying coupling strength v.s. Âµ for different noise level , find where the bifurcation happens; mean and v a riance of the P oinca r Â´ e section, fi nd where the bifurcation hap p ens; surface plot of mean and variance of each oscillator v.s. and Âµ for fixed noise level (i.e. 2 , 4 = 0 . 05). The results for the first oscillator in the system shows in the figures below. PAGE 74 65 1 Figure A.1: Coupled System of Stochastic Rayleigh 2 , 4 = 0 . For the deterministic twp coupled Rayleigh van der Pol equations, we can see from the Figure A.1 that the bifurcation happened at = 2 Âµ which is one of the eigenvalues for the system. PAGE 75 66 2 2 Figure A.2: Coupled System of Stochastic Rayleigh van 2 , 4 = 0 . 01 , the gray points indicate the bifurcation region. The intriguing part happened once we introduce stochastic term into the coupled system as showed in figure A.2. Even for = 0 . 01, the bifurcation transit from a line The bifurcation region is a combination of positive part of = 1 Âµ and negative part of = 1 Âµ both shift to the right of axis. Moreover, the bifurcation region widened more for positive . PAGE 76 67 1 1 2 2 Figure A.3: Coupled System of Stochastic Rayleigh 2 , 4 = 0 . 05 . For = 0 . 05 which showed in figure A.3, the bifurcation region is still a combina tion of positive part of = Âµ and negative part of = Âµ . The bifurcation region keep shifting to the right of axis while widened more for positive than negative . PAGE 77 68 1 1 2 2 Figure A.4: Coupled System of Stochastic Rayleigh 2 , 4 = 0 . 1 . For = 0 . 05 which showed in figure A.4, the bifurcation region is still a combina tion of positive part of = Âµ and negative part of = Âµ . The bifurcation region keep shifting to the right of axis while widened more for positive than negative . PAGE 78 69 Figure A.5: Coupled System of Stochastic Rayleigh van der Pol, mean of the first oscillator v.s. Âµ for different 2 , 4 = 0 . 05 . Now we fixed the stochastic term(i.e. = 0 . 05), we want to see how the statistical properties behave for different coupling strength . Figure A.5 suggests that only the magnitude of will influence the mean of the system. Increasing the magnitude of , has the opposite effect as increasing . The increase of the magnitude of value of mean bigger while shifting the bifurcation point/region to the left of the axis. PAGE 79 70 Figure A.6: Coupled System of Sto chastic Rayleigh van der Pol, variance of the first oscillator v.s. Âµ for different 2 , 4 = 0 . 05 . Meanwhile, Figure A.6 suggests similar tendency, that only the magnitude of will influence the mean of the system. Increasing the magnitude of , has the opposite effect as increasing . The increase of the magnitude of make the decrease the amplitude of the peak region of bifurcation region while narrowing and shifting the bifu rcation region to the left of the axis. PAGE 80 71 Figure A.7: Coupled System of Stochastic Rayleigh van der Pol, surface plot of the 2 , 4 = 0 . 05 . Figure A.8: Coupled System of Stochastic Rayleigh van der Pol, surface plot of the 2 , 4 = 0 . 05 . Figure A.7 and A.8 shows the surface plots for mean, Âµ and of the first oscillator at different angles. The surface plots give a better visualization about mean, Âµ , and bifurcation at fixed stochastic terms ( 2 , 4 = 0 . 05): mean increases as the magnitude of ion region is. PAGE 81 72 Figure A.9: Coupled System of Stochastic Rayleigh van der Pol, surface plot of the 2 , 4 = 0 . 05 . Figure A.10: Coupled System of Stochastic Rayleigh van der Pol, surface plot of the 2 , 4 = 0 . 05 . Figure A.9 and A.10 shows the surface plots for variance, Âµ and of the first oscil lator at different angles. The surface plots give a better visualization about variance, and bifurcation at fixed stochastic terms ( 2 , 4 = 0 . 05): variance decreases as the PAGE 82 73 magnitude of increase, the peaks indicate the bifurcation region. There seems to Figure A.11: Coupled System of Stochastic Rayleigh van der Pol, surface plot of the second oscillator for mean, Âµ and 2 , 4 = 0 . 05 . Figure A.12: Coupled System of Stochastic Rayleigh van der Pol, surface plot of the second oscillator for mean, 2 , 4 = 0 . 05 . Figure A.11 and A.12 shows the surface plots for mean, Âµ and of the second oscillator at different angles. The surface plots give a better visualization about PAGE 83 74 mean, and bifurcation at fixed stochastic terms ( 2 , 4 = 0 . 05): mean increases as the magnitude of bifurcation region is. Figure A.13: Coupled System of Stochastic Rayleigh van der Pol, surface plot of the second oscillator for variance, Âµ and 2 , 4 = 0 . 05 . Figure A.14: Coupled System of Stochastic Rayleigh van der Pol, surface plot of the second 2 , 4 = 0 . 05 . Figure A.13 and A.14 shows the surface plots for variance, Âµ and of the second PAGE 84 75 oscillator at different angles. The surface plots give a better visualization about variance and bifurcation at fixed stochastic terms ( 2 , 4 = 0 . 05): variance decreases as the magnitude of increase, the peaks indicate the bifurcation region. There seems both bifurcation regions are a little different compare to the first osc illator. The coupled system behaves a little different than we anticipated, the nature of the coupled stochastic nonlinear system and its behavior is not well understand yet. As we can tell from the plots, the bifurcation region has an excepted shift once the randomness has been introduced to the system, then the bifurcation region widens and shift to the right as the stochasticity increases. For the fixed , both positive and negative coupling shift the bifurcation region to the left, the bifurcation narrowed as the absolute value of increases. |