Generation of stochastic mixed-mode oscillations in a pair of VDP oscillators with direct-indirect coupling

1.   Introduction

As we all know, any real-world dynamic system is influenced by random noise. In fact, noise can give rise to a wide range of complex phenomena, including coherence resonance ^[1,2,3], stochastic resonance ^[4,5], noise-induced bursting ^[6,7,8], noise-induced stochastic mixed-mode oscillations ^[9,10,11], noise-induce suppression of firing ^[12], noise-induced chaos and order ^[13,14], stochastic oscillation bistability in the zone of canard limit cycles ^[15] and noise-induced transitions between tonic spiking and bursting regimes ^[16].

For small random perturbations, Fredidlin and Wentzell proposed a large deviation theory to explain the long-term effects of them ^[17]. Specifically, it is based on the fact that when almost impossible events occur, they occurred along a specific path. This path is called the optimal path which produces the minimum value of some action function. It describes how difficult it is for a random process to pass through the neighborhood of a particular path. The quasi-potential of a particular point is then generated by the action functional along the optimal path of its connection to the attractor. The stochastic sensitivity function and confidence ellipse method was proposed by using the second-order approximation of the potential. It has been successfully used in many systems, such as Hodgkin-Huxley model ^[12], hair bundle model ^[9], prey-predator-plankton system ^[18], Morris-Lecar system ^[19], Higgins model ^[20], Hindmarsh-Rose model ^[21], etc.

The Van der Pol (VDP) oscillator is typical self-excited or self-sustained oscillations, which first originated in vacuum tube circuits. It has been studied in detail in the last century, and the main results can be found in the famous monograph by Nayfeh and Mook ^[22]. Koshcheev and Vladimir showed the white-noise-induced transition between the limiting cycle and state of rest of a VDP oscillator that has a threshold, and the threshold value is directly proportional to the product of the characteristic energy of self-oscillations and the friction coefficient ^[23]. For the coupled VDP oscillator, some coupling schemes are used, including weak coupling with delay ^[24], indirect coupling ^[25] and dynamic coupling ^[26]. In the direct-indirect coupling VDP oscillator coupling system, we have found a wealth of kinetic phenomena ^[26]. If we introduce noise into it, we find more dynamics in it. In this paper, we show the noise-induced escape on a direct-indirect coupled pair of VDP systems and generate mixed-mode stochastic oscillations.

This paper is organized as follows. We first propose a mathematical model of the pair of VDP oscillator systems with direct-indirect coupling and give its equilibrium distribution with changing parameters in Section 2. In Section 3, we show the oscillations of the corresponding random system under different noise intensities. In Section 4, we give the reason of this phenomenon by using random sensitivity matrix and confidence ellipse method. Conclusions and discussions are given in Section 5.

2.   Deterministic model

In this paper, we investigate the VDP coupled system with two types of connections, i.e., a direct and an indirect coupling. The direct coupling is diffusive mutually coupling, and the indirect coupling is carried by using external environment. The system dynamics is presented ^[27] as

where x_i1, x_i2 are firing activities of the i-th oscillator, a is a nonlinear damping ratio of the individual oscillator, d is a coupling strength of direct coupling, the state variable y denotes indirect coupling of the external environment, which is modeled by a one-dimensional over-damped oscillator with damping parameter k and the parameter e is a coupling strength between systems and environment.

We carry out the equilibrium analysis of system (1) and obtain the following equations:

By solving the above equation, we get four types of equlibrium points, including the origin O (0, 0, 0, 0, 0), inhomogeneous steady state (IHSS) ${A^ \pm }\left({ \pm x_{11}^A, \pm x_{12}^A, \mp x_{11}^A, \mp x_{12}^A, 0} \right)$, homogeneous steady state (HSS) ${B^ \pm }\left({ \pm x_{11}^B, \pm x_{12}^B, \pm x_{11}^B, \pm x_{12}^B, \mp {y^B}} \right)$ and no-pattern steady state (NPSS) ${C^ \pm }\left({ \pm x_{11}^C, \pm x_{12}^C, \mp x_{21}^C, \mp x_{22}^C, \pm {y^C}} \right)$, ${D^ \pm }\left({ \pm x_{11}^D, \pm x_{12}^D, \mp x_{21}^D, \mp x_{22}^D, \pm {y^D}} \right)$, where $x_{11}^A = \sqrt {2ad - 1} /\sqrt {2ad}$, $x_{12}^A = \pm \sqrt {2d} \sqrt {2ad - 1} /\sqrt a$, $x_{11}^B = \sqrt {a{e^2} - k} /e\sqrt a$, $x_{12}^B = e\sqrt {a{e^2} - k} /k\sqrt a$ and ${y^B} = \sqrt {a{e^2} - k} /k\sqrt a$. Because ${C^ \pm }$ and ${D^ \pm }$ are more complex, they are not listed here.

In order to show the distribution position and the properties of the equilibrium points, we fix $a = 1, d = 1, k = 1.3$ in this paper, and consider the dynamics of system (1) under different values of parameter $e$. When $e > 1.2465$, the system enters the equilibrium state. There are three types of equilibria: A pair of IHSS equilibria ${A^ \pm }$, a pair of HSS equilibria ${B^ \pm }$ and trivial equilibrium O. ${A^ \pm }$ are stable equilibrium points. ${B^ \pm }$ and O are unstable equilibrium points. A pitchfork bifurcation occurs, adding four unstable equilibrium points: NPSS equilibria ${C^ \pm }$and ${D^ \pm }$. Next, when $e = 1.6125$, the equilibria ${B^ \pm }$gains stability. Finally, an anti-pitchfork bifurcation appears, the points ${C^ \pm }$ and ${D^ \pm }$ disappear, and the points ${A^ \pm }$ lose their stability, leaving only a pair of stable equilibria ${B^ \pm }$.

Figure 1 shows three examples of different types of stable equilibria in the deterministic system characterized by the nullclines. As can be seen from Figure 1(a), when $e = 1.28$, there are five equilibrium points in the system: Unstable points O and ${B^ \pm }$ and stable fixed points ${A^ \pm }$. When $e = 1.9$, as shown in Figure 1(b), there are two pairs of stable fixed points of the system: ${A^ \pm }$ and ${B^ \pm }$, respectively. The unstable point is the point O, and there are two pairs of saddle points: ${C^ \pm }$ and ${D^ \pm }$. When $e = 3$, as shown in Figure 1(c), the system re-obtains the stable equilibrium points ${B^ \pm }$, where ${A^ \pm }$ are saddle points, and O is an unstable point. At the same time, the corresponding phase diagrams for different values of parameter $e$ with different initial values are given, as shown in Figure 2. It can be observed that for the fixed parameters of the system, different initial values can make the system fall into different equilibrium points.

Figure 1.  The intersection points of the system's nullclines and vector field show number of equilibria with the fixed parameter (a) e = 1.28, (b) e = 1.9, (c) e = 3, where O is the trivial equilibrium, ${A^ \pm }$ means the IHSS equilibrium, ${B^ \pm }$ means the HSS equilibrium, and ${C^ \pm }$ and ${D^ \pm }$ are the NPSS equilibria.

DownLoad: Full-Size Img PowerPoint

Figure 2.  Phase diagrams (a) e = 1.28 with initial values (−1, 0, 0, 0, −1), (−1, 0, 0, 0, 1), (1, 0, 0, 0, −1) and (1, 0, 0, 0, 1); (b) e = 1.9 with (−2, 0, 0, 0, −2), (−2, 0, 0, 0, 1), (2, 0, 0, 0, 1) and (2, 0, 0, 0, 2); (c) e = 3 with (−1, 0, 0, 0, −1), (−1, 0, 0, 0, 1), (1, 0, 0, 0, −1) and (1, 0, 0, 0, 1) for the fixed parameters a = 1, d = 1, k = 1.3.

DownLoad: Full-Size Img PowerPoint

3.   Stochastic model

Let us examine how random noise changes the system dynamics of model (1). Here we will study the following stochastic model:

where $\xi (t)$ is a standard Gaussian white noise with $E\xi (t) = 0, E\xi (t)\xi (\tau) = \delta (t - \tau)$ and $\varepsilon$ is a noise intensity parameter.

In order to illustrate the effect of Gaussian white noise, we focus on the stochastic phenomena for the different values of parameter $e$ with small and large noise intensity. For the numerical simulations of random trajectories, we use the Euler-Maruyama scheme with the time step 0.01 in accordance with practice.

Figure 3 illustrates the random phase trajectories and the corresponding history series with $e = 1.28, {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} 1.9, {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} 3$ for the low noise intensity $\varepsilon = 0.02$. We can obviously find that the random trajectories vibrate in the neighborhood of the stable fixed points and form dispersion. Further, the small noise intensity does not change the overall property of the model. However, the system enters a switch of equilibrium points under large noise intensity, as shown in Figure 4, where the history series are to the parameters $e = 1.28, {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} \varepsilon = 0.2$, $e = 1.9, {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} \varepsilon = 0.5$ and $e = 3, {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} \varepsilon = 1$, respectively. As can be seen from the figures, we can observe large amplitude stochastic oscillations around the equilibrium point and random switching between the two equilibrium points. The noise-induced transition phenomenon is observed.

Figure 3.  (Up) The history series of the stochastic system and corresponding deterministic system and (Down) phase diagram with (a) e = 1.28, (b) e = 1.9, (c) e = 3, where the system parameter values are a = 1, d = 1, k = 1.3 for the noise intensity 0.02.

DownLoad: Full-Size Img PowerPoint

Figure 4.  Time history diagrams of the stochastic system (2) with (a)$e = 1.28, {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} \varepsilon = 0.2$, (b)$e = 1.9, {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} \varepsilon = 0.5$, and (c) $e = 3, {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} \varepsilon = 1$.

DownLoad: Full-Size Img PowerPoint

Figure 5 shows the relationship between noise intensity and the coordinate of the stochastic system. The figure depicted here has its horizontal axis representing the intensity of noise, while the vertical axis corresponds to ${x_{11}}$. The system begins integrating at a noise intensity of zero, and subsequently carries out the same operation at intervals of ${\log _{10}}\left({0.01} \right)$. This process illustrates the range of coordinate ${x_{11}}$ of the system. As can be seen, when the noise intensity is very small, the oscillation amplitude caused by the random noise is also very small. The greater the noise intensity, the greater the amplitude of the oscillation. After the noise intensity of different threshold values, the random system will enter the noise-induced large amplitude oscillation.

Figure 5.  Bifurcation diagram of the noise intensity for (a) e = 1.28, (b) e = 1.9, and (c) e = 3.

DownLoad: Full-Size Img PowerPoint

4.   Stochastic sensitivity analysis

In this section, we first give the attraction basins of all equilibria for $e = 1.28, {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} 1.9, {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} 3$ in the deterministic system (1), as shown in Figure 6. As shown in the figure, areas of different colors represent the attraction basins of different equilibrium points. We take an interval of 0.01, with the value range of ${x_{11}}$ and ${x_{21}}$ being −2 to 2, and set the initial values to satisfy equation (2). Numerical integration is performed using the Runge-Kutta algorithm. If the system falls into the same equilibrium point, the corresponding point will be painted with the same color. As can be seen, the position of the stable fixed point, which the system finally stays in, depends on the position of the initial point of the deterministic system. Meanwhile, there are some boundary surfaces between different basins. If the trajectory of random system (2) crosses the basin of attraction, the system will inevitably fall into another equilibrium point. In stochastic systems, random terms will cause perturbations to phase trajectories. The larger the random coefficient, the larger the range of random system disturbance. If the perturbation of a random system breaks through the boundary of the attractive disk, the system will inevitably fall into another equilibrium point of the attractive disk range. The system vibrates up and down around another equilibrium point.

Figure 6.  The attraction basins of all equilibria for (a) $e = 1.28$, (b) $e = 1.9$, (c) $e = 3$.

DownLoad: Full-Size Img PowerPoint

In order to analyze the noise-induced random oscillation, we use the confidence ellipse method based on the stochastic sensitivity matrix technique ^[27]. The technique of the stochastic sensitivity matrix that we employ is fundamentally developed on the basis of the Hamilton-Jacobi equation. Before analyzing, we give the Jacobian matrix of the deterministic system (1).

The random sensitivity matrix W of the stable equilibrium point $\left({{x_{11}}, {x_{12}}, {x_{21}}, {x_{22}}, y} \right)$ for the deterministic system (1) can be obtained by the following matrix equation:

The random sensitivity matrix of the different stable fixed points can be obtained by solving the corresponding matrix equation. Using the matrix W, one can write an asymptotic of the time-varying probability density function

The eigenvalues and eigenvectors can form the confidence ellipse of the corresponding points. The eigenvalue determines the semi-axis of the confidence ellipse, and the eigenvector determines the direction axis of the corresponding confidence ellipse. Let ${\mu _1}, {\mu _2}, {\mu _3}, {\mu _4}, {\mu _5}$ be eigenvalues and ${u_1}, {u_2}, {u_3}, {u_4}, {u_5}$ be the normalized eigenvectors of the stochastic sensitivity matrix. For coordinates ${z_1}, {z_2}, {z_3}, {z_4}, {z_5}$ of the confidence ellipse in the basis ${u_1}, {u_2}, {u_3}, {u_4}, {u_5}$ with the equilibrium as the origin, we have the following equation:

where $\varepsilon$ is the noise intensity and $K{\text{ (}}P)$ is the confidence probability function.

4.1.   Random sensitivity analysis of the equailibrium points ${A^ \pm }$

Substituting the equilibrium point ${A^ \pm }$ into the matrix equation, we can solve the corresponding random sensitivity matrix:

where ${\Delta _3} = 2({e^2} - 6k)(- 6 + k + 4{e^2}k - 2{k^2})$. It is easy to know that when $e = \sqrt {6k}$ or $e = \sqrt {(6 - k + 2{k^2})/4k}$, i.e., $e = 1.2465$ or $e = 2.7928$, the random system will immediately enter the random oscillation state. Here, when we set $e = 1.2465$, points ${A^ \pm }$ turn into the stable equilibrium points and when we set $e = 2.7928$, points ${A^ \pm }$ lose their stability. Figure 7 shows the eigenvalue of the stochastic sensitivity matrix of the equilibrium point ${A^ \pm }$ described by the function with parameter $e$ as the independent variable. It is obvious that the random sensitivity of the system will change with the system parameter $e$. However, when $e$ approaches $1.2465$ or $2.7928$, the random sensitivity of the system increases indefinitely.

Figure 7.  Eigenvalue function of the stochastic sensitivity matrix for the equilibrium point ${A^ \pm }$.

DownLoad: Full-Size Img PowerPoint

Figure 8 also shows the confidence ellipses when we set $e = 1.28, {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} 1.9$ and $\varepsilon = 0.1, {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} 0.2$. As can be seen from the figure, when we set $\varepsilon = 0.1$, the confidence ellipses are still in the basin of attraction for the corresponding equilibrium, but when we set $\varepsilon = 0.2$, the confidence ellipses become larger and cross the dividing line into the attraction basin of another equilibrium state. As can be seen from Figure 8(a), the random trajectory will enter into another IHSS equilibrium. The random trajectory oscillates back and forth between the equilibrium ${A^ \pm }$. In Figure 8(b), when we set $\varepsilon = 0.2$, the random trajectory will be equally possible to enter into the HSS equilibria ${B^ \pm }$. The motion of the stochastic system at the equilibria ${B^ \pm }$ will be given in the calculation results in the following subsection. Figure 9 shows the time history diagram and stochastic bifurcation diagram of the corresponding stochastic system.

Figure 8.  The attraction basin of equilibrium ${A^ + }$ (by color) and the corresponding confidence ellipse for $\varepsilon = 0.1$ (small ellipse) and $\varepsilon = 0.2$ (large ellipse) in the deterministic system with the fixed parameter (a) $e = 1.28$ and (b) $e = 1.9$.

DownLoad: Full-Size Img PowerPoint

Figure 9.  (a) Time history diagrams for $e = 1.9$ and $\varepsilon = 0.2$. (b) Bifurcation diagram between the noise intensity and the coordinate $y$.

DownLoad: Full-Size Img PowerPoint

4.2.   Random sensitivity analysis of the equailibrium points ${B^ \pm }$

Considering the random sensitivity matrix equation corresponding to the equilibrium point ${B^ \pm }$, we get the following random senstivity matrix:

where ${\Delta _4} = 4({e^2} - k)(- 2{e^4} + {e^6}k + {k^3} + {e^2}(k - {k^3}))$. However, when HSS equilibria ${B^ \pm }$ enters the stable state, ${\Delta _4} > 0$. Figure 10 shows the eigenvalue of the stochastic sensitivity matrix of the equilibrium points ${B^ \pm }$ described by the function with parameter $e$ as the independent variable. Similar to Figure 7, the random sensitivity of the system will change as $e$ changes.

Figure 10.  Eigenvalue function of the random sensitivity matrix for the equilibrium point ${B^ \pm }$.

DownLoad: Full-Size Img PowerPoint

Figure 11 shows the confidence ellipses of the large and small noise intensity when the point ${B^ \pm }$ are stable fixed points. The confidence ellipse then degenerates into a line segment. As can be seen from the figure, when we set $e = 1.9$ and $\varepsilon = 0.3$, the confidence ellipse is still in the basin of attraction of the corresponding equilibrium, but when $\varepsilon = 0.5$, the confidence ellipse becomes larger and crosses the boundary entering into the attraction basin of the another equilibrium. The random trajectory of system (2) jumps back and forth between ${B^ + }$ and ${B^ - }$. However, due to the property of the confidence ellipse, the system cannot jump into the equilibrium point ${A^ \pm }$. When we set $e = 3$ and $\varepsilon = 0.5$, the random will not cause the transition between two equilibria. Further, when $\varepsilon = 1$, the phenomenon of transition of the equilibrium point can be observed. The corresponding time history is shown in Figure 4 (c).

Figure 11.  The confidence ellipse described by the attraction basin of the equilibrium point ${B^ \pm }$ for the deterministic system with (a) $e = 1.9$, $\varepsilon = 0.3$ (b) $e = 1.9$, $\varepsilon = 0.5$, (c) $e = 3$, $\varepsilon = 0.5$, and (d) $e = 3$, $\varepsilon = 1$.

DownLoad: Full-Size Img PowerPoint

5.   Conclusions and discussions

The phenomenon of a pair of VDP oscillator affected by Gaussian white noise was studied. In this paper, we showed that the random noise can produce mixed-mode oscillations in the parameter region for the deterministic model having only the stable equilibrium point as the attractor. The phenomenon of spontaneous mixed-mode oscillation induced by noise was studied using the semi-analytic method based on the random sensitivity function and confidence ellipse.

The results showed that the mechanism of this phenomenon can be explained by the specificity of the attraction basin and the random sensitivity of the equilibrium point. When a random noise term was introduced into the system, the system vibrated around the equilibrium point. The greater the random noise intensity, the larger the confidence ellipse of system vibration. When the confidence ellipse crossed the boundary of the attraction basin at that point, the system would no longer vibrate around this equilibrium point. It jumped into a different equilibrium point. It is worth noting that when the system parameters allowed the system to have two equilibrium points of ${A^ \pm }$ and ${B^ \pm }$, the deterministic system at ${A^ \pm }$ would vibrate around the equilibrium point of ${B^ \pm }$ and would not return to ${A^ \pm }$ after the random strength became more larger. This article can give us a better understanding of the behavior of nonlinear dynamical systems. This study opened up new avenues for further research. For example, it would be very interesting to study how different types of noise (e.g., Gaussian noise, Poisson noise, etc.) affect the behavior of VDP oscillators and other similar systems.

Use of AI tools declaration

The authors declare they have not used Artificial Intelligence (AI) tools in the creation of this article.

Acknowledgments

This research is supported by the National Natural Science Foundation of China (No. 12172212) and the Fundamental Research Funds for the Central Universities (No. 22120220588).

Conflict of interest

The authors declare there is no conflict of interest.