1. Introduction
Herein, we numerically investigate the fundamental properties of a canard explosion generated by a Bonhoeffer-van der Pol (BVP) oscillator using the parameter values obtained experimentally. Furthermore, we investigate two identical coupled canard-generating oscillators. However, the experimental results exhibit complex phenomena, which cannot be explained considering two identical autonomous coupled BVP oscillators. Thus, we construct a numerical model that is capable of accurately reproducing the experimentally observed results; this model includes two weak perturbations having the same frequency.
Canards are among the most important discoveries of the late 20th century, and they have been studied extensively for over four decades [1]-[11]. Here, we consider a van der Pol oscillator of the following form:
| \[\begin{equation*} \left\{ \begin{array}{l} \varepsilon\dot{x} = y+x-x^3, \\ \dot{y}= -x+B_0, \end{array} \right. \tag{1} \end{equation*}\] |
where \(\varepsilon\) is a small parameter and \(B_0\) represents the amplitude of the DC voltage source considered here. It has been demonstrated through a nonstandard analysis that there exists a case where the amplitude of the oscillator changes by a factor on the order of 1 when the value of \(B_0\) changes by a factor on the order of \(\exp(-1/\varepsilon)\) [2], [6]. The magnitude of \(\exp(-1/\varepsilon)\) is \(\sim\)\(4.5\times10^{-5}\) for \(\varepsilon=0.1\) and \(\sim\)\(10^{-44}\) for \(\varepsilon=0.01\). The observed oscillation in this system is referred to as a canard because the form of the oscillation resembles a duck. The oscillatory behavior observed here is notable because canard explosions are not phenomena in which the characteristic multipliers cross the threshold.
Herein, we first numerically analyze a slow-fast BVP oscillator that generates a canard explosion; we undertake this work utilizing the parameter values obtained experimentally. As the inductor in the system includes an internal resistance, we assert that BVP oscillators represent more natural and realistic systems than van der Pol oscillators. In our numerical study, we find the magnitude of the canard oscillation is changed by a factor on the order of 1 when parameter \(B_0\) is changed by 0.0001; this finding is in good qualitative agreement with the results of a previous nonstandard analysis.
We then discuss the synchronization of identical and nearly identical canard generating BVP oscillators. Recently, coupled canard- and mixed-mode oscillation (MMO)-generating oscillators have been the subject of intense research [12]-[17]. In such systems, the synchronization of canards and MMOs can be observed [13]-[15].
In our previous work, the complete and in-phase synchronizations of canards in identical and nearly identical coupled canard-generating BVP oscillators have been investigated [18]. Understanding how canard-synchronization occurs is of particular interest because canards are extremely sensitive to parameter values. Numerical simulations on coupled canard-generating oscillators have been undertaken previously, and both the complete and in-phase synchronizations have been observed. It has also been established why a strong coupling is necessary to achieve complete synchronization in experimental work [18]. Herein, we assume that capacitance mismatch of 10\(\%\) in the two oscillators; this mismatch is taken as it reflects the probable inaccuracies in the experimentally realized circuit. The 10\(\%\) mismatch is considered to be reasonable because it is probable that the nonlinear conductors include a 10\(\%\) parameter mismatch due to the difficulty of realizing elements with an identical nonlinear conductance using diode arrays. In a previous study, it has been confirmed that an error of \(\sim\)10\(\%\) can generate almost complete synchronization of the oscillators [18].
Herein, we consider a system that accurately represents an experimentally realized ciruit and demonstrates numerically that almost-complete synchronization is no longer stable for a coupling parameter of \(\alpha=0.4,\) and that the in-phase synchronization emerges instead. Here, \(\alpha\) is the coupling parameter that represents the coupling conductance, \(g\). Conversely, we observe the almost-complete synchronization for \(\alpha\simeq 1.2\). However, in the experimental work, we find that the oscillator demonstrates a stable in-phase synchronization for values of \(\alpha\) below \(1.2\); in this case, the phenomena observed experimentally cannot be replicated numerically.
In previous work, the complete synchronization of canards was observed across a wide range of coupling parameter values in numerical simulations in identical coupled canard-generating oscillators [18]; this synchronization can occur for relatively small coupling parameter values, such as \(\alpha=0.03\). However, in the circuit experiments undertanken in this work, we observe that the complete synchronization collapses as a result of even a slight decrease in the coupling parameter value at which complete synchronization can be seen. In this work, we have found a discrepancy between the numerical simulations and actual circuit experiments; indeed, we have observed a complex synchronization that emerges in experiments that had not previously been reported. Furthermore, the destruction of the complete synchronization of canards observed in experiments cannot be explained by considereing coupled canard-generating four-variable autonomous ordinary differential equations (ODEs) that include a 10% parameter mismatch, despite these ODEs would be a good approximation of the dynamics.
To explain such a discrepancy, in this paper, we investigate this phenomenon by constructing a set of four-variable nonautonomous coupled BVP ODEs with a 10% parameter mismatch and two weak periodic perturbations. By considering these non-autonomous ODEs, we were able to observe oscillatory behavior that was consistent with that present in the circuit experiments. These weak periodic perturbations could be considered to represent the noise that is present in an actual dynamical circuit. We refer to this complex synchronization behavior as “butterfly synchronization” of canards, a name that reflects the shape of the attractors. The butterfly synchronization of canards occurs when the phase difference of the canards is sufficiently large. We numerically confirm that such a phase difference can be achieved by considereing two weak periodic perturbations. We note that the butterfly synchronization of canards occurs not because of the change in shape of the attractors, but because the phases of the two canards are shifted.
2. Canard Explosion Generated by a BVP Oscillator
In this section, we accurately describe a canard explosion generated by a canard-generating oscillator using parameter values that reflects the properties of a circuit considered experimentally. The circuit considered here is shown in Fig. 1; it is referred to as a BVP oscillator. In Fig. 1, \(C\), \(L\), \(R\), \(E_0\), and N.C. represent the capacitance of the capacitor, the inductance of the inductor, a parasitic resistance of the inductor, a DC voltage source, and a nonlinear negative conductance, respectively.
 |
Fig. 1 The circuit diagram of the BVP oscillator considered here. |
The voltage-current characteristics of the nonlinear conductance (labeled N.C.), is assumed to be represented by the following third-order polynomial function:
| \[\begin{equation*} G(v) = -g_1v+g_3v^3. \\ \tag{2} \end{equation*}\] |
The circuit diagram for the element N.C. is given in Fig. 2(a). Herein, \(r_1\), \(r_2\), \(r_3\) and \(r_4\) in Fig. 2(a) are set to 2.0 k\(\Omega\), 3.0 k\(\Omega\), 3.2 k\(\Omega\), and 1.4 k\(\Omega\) respectively. Figure 2(b) shows the characteristic voltage-current curve obtained from the circuit shown in Fig. 2(a); in this figure, the experimental results are shown in blue and the polynomial function given in Eq. (2) is shown by a red curve.
 |
Fig. 2 (a) Realization of a nonlinear conductance via diode arrays (the values of the resistance used herein are \(r_1 = 2.0\) k\(\Omega\), \(r_2 = 3.0\) k\(\Omega\), \(r_3 = 3.2\) k\(\Omega\), and \(r_4 = 1.4\) k\(\Omega\)), (b) Voltage-current characteristics of the circuit shown in (a). |
The governing equation of the circuit shown in Fig. 1 can be written as
| \[\begin{equation*} \left\{ \begin{array}{l} C\displaystyle\frac{dv}{dt} = -i+g_1v-g_3v^3, \\ L\displaystyle\frac{di}{dt} = -v-iR+E_0. \\ \end{array} \right. \tag{3} \end{equation*}\] |
The parameter values measured experimentally are as follows.
| \[\begin{equation*} \begin{array}{l} C\simeq10~{\rm nF}, L\simeq100~{\rm mH}, R\simeq540~\Omega,\\ g_1\simeq{370~\mu {\rm S}},~{\rm and}~g_3\simeq2.3~{\mu {\rm S}}. \end{array} \tag{4} \end{equation*}\] |
Considering the following rescaling:
| \[\begin{equation*} \begin{array}{l} v=\sqrt{\displaystyle\frac{g_1}{g_3}}x,~i=-g_1\sqrt{\displaystyle\frac{g_1}{g_3}}y,~{t=Lg_1\tau}, \\ E_0=\sqrt{\displaystyle\frac{g_1}{g_3}}B_0,~k = g_1R, \varepsilon=\displaystyle\frac{C}{Lg_1^2}, \end{array} \tag{5} \end{equation*}\] |
we can normalize the equations of the system; the system can then be described by the following system of two-variable autonomous ordinary differential equations (ODEs):
| \[\begin{equation*} \left\{ \begin{array}{l} \varepsilon\dot{x} = y+x-x^3, \\ \dot{y}= -x-ky+B_0. \end{array} \right. \tag{6} \end{equation*}\] |
In Eq. (6), \(B_0\) describes the DC bias voltage \(E_0\), and \(\varepsilon\) is a parameter corresponding to the small capacitance \(C\). Herein, \(\varepsilon\) is assumed to be small. Considering the transformations described in Eq. (5) and the parameter values obtained from the experimental work, the normalized parameters are set as \(\varepsilon=0.1\) and \(k=0.1998\); the parameter \(B_0\) is used as the bifurcation parameter in the system. In the following work, we describe the solutions of this system for a range of values of \(B_0\).
Figure 3 shows the behavior of the solutions of Eq. (6) in the \(x-y\) plane. In Fig. 3, the lines in blue and green represent \(x\)- and \(y\)-nullclines, respectively. When the DC bias voltage is relatively large (for example, when \(B_0 = 0.5\)), the intersection of these nullclines becomes a stable equilibrium; this point is marked with the black dot in Fig. 3(a). As \(B_0\) decreases to 0.4898, a limit cycle with small amplitude emerges via a supercritical Hopf bifurcation as shown in Fig. 3(b). This attractor is a canard without a head. For \(B_0 {=} 0.4897\), which is 0.0001 smaller than the value of \(B_0\) for which a canard without a head can be observed, a canard with a head can be seen to emerge [19]. Once a canard with a head has been generated, the change in the magnitude of its amplitude becomes gradual with respect to the change in \(B_0\), as shown in Fig. 3(d).
 |
Fig. 3 Phase planes for various value of \(B_0\) (\(\varepsilon=0.1\) and \(k=0.1998\)): (a) Stable equilibrium (\(B_0=0.5\)), (b) canard without a head (\(B_0=0.4898\)), (c) canard with a head (\(B_0=0.4897\)), and (d) large-amplitude solution (\(B_0=0.4\)). |
Figure 4 shows a one-parameter bifurcation diagram where the horizontal axis represents the bifurcation parameter \(B_0\) and the vertical axis represents the largest and smallest values of \(x\). The equilibrium generated in the system described by Eq. (6) is stable for \(B_0 >0.4946604\) and unstable for \(0<B_0 < 0.4946604\). Stable and unstable equilibria are marked with solid and dashed lines, respectively. The transition of the equilibrium from stable to unstable occurs at \(B_0 = 0.4946604\) via a supercritical Hopf bifurcation. This value of \(B_0\) can be manually obtained. Immediately after this bifurcation, a small amplitude limit cycle is generated. Figure 4 shows that extremely rapid changes in the amplitude of the solution occur at around \(B_0=0.49\). This abrupt change in the magnitude of the oscillation is called a canard explosion [19]. A canard explosion is a notable phenomenon because it represents a bifurcation in the more general sense as an important change occurs; however, the canard explosion is not a phenomenon in which characteristic multipliers cross a threshold value.
 |
Fig. 4 One-parameter bifurcation diagram of \(x\) plotted as a function of \(B_0\) (\(\varepsilon=0.1\) and \(k=0.1998\)). |
3. Numerical Analysis for Two-Coupled BVP Oscillators with a Mismatch in Capacitances
In this section, we discuss the synchronization of canards that are generated by identical coupled BVP oscillators. Figure 5 shows a circuit diagram of two-coupled BVP oscillators with a conductance \(g\); in this circuit, the parameters describing and each oscillator are identical. However, notably, a slight parameter mismatch exists between the two oscillators because complete synchronization (that should not exist) may occur as a result of the finite precision of the computations undertaken. Based on the difficulty of realizing the qualitatively similar elements with a nonlinear conductance constructed using diode arrays, we intentionally introduce a 10\(\%\) parameter mismatch in the two capacitances \(C_1\) and \(C_2\).
 |
Fig. 5 Resistively coupled BVP oscillators. |
The governing equation of the coupled circuits can be expressed by the following system of four autonomous ODEs:
| \[\begin{equation*} \left\{ \begin{array}{l} C_1\displaystyle\frac{dv_1}{dt} = i_1- {G_1}(v_1)-g(v_1-v_2), \\ L\displaystyle\frac{di_1}{dt} = E_0-v_1- i_1R, \\ C_2\displaystyle\frac{dv_2}{dt} = i_2-G_2(v_2)-g(v_2-v_1), \\ L\displaystyle\frac{di_2}{dt} = E_0-v_2-i_2R. \end{array} \right. \tag{7} \end{equation*}\] |
The variable transformations,
| \[\begin{equation*} \begin{array}{@{}l@{}} v_1=\sqrt{\displaystyle\frac{g_1}{g_3}}x_1,~v_2=\sqrt{\displaystyle\frac{g_1}{g_3}}x_2,~i_1=-g_1\sqrt{\displaystyle\frac{g_1}{g_3}}y_1, \\ i_2=-g_1\sqrt{\displaystyle\frac{g_1}{g_3}}y_2,~t=Lg_1\tau,~E_0=\sqrt{\displaystyle\frac{g_1}{g_3}}B_0,\\[2mm] k = g_1R,~\varepsilon_1=\displaystyle\frac{C_1}{Lg_1^2},~ \varepsilon_2=\displaystyle\frac{C_2}{Lg_1^2},~\alpha=\displaystyle\frac{g}{g_1}, \end{array} \tag{8} \end{equation*}\] |
can be used to write the above expressions in a normalized form:
| \[\begin{equation*} \left\{ \begin{array}{l} \varepsilon_1\dot{x}_1 = y_1+x_1-x_1^3-\alpha(x_1-x_2), \\ \dot{y}_1 = -x_1-ky_1+B_0, \\ \varepsilon_2\dot{x}_2 = y_2+x_2-x_2^3+\alpha(x_1-x_2), \\ \dot{y}_2 = -x_2-ky_2+B_0. \end{array} \right. \tag{9} \end{equation*}\] |
Here, \(B_0\), which is the parameter that contains information related to the DC bias voltage \(E_0\), was set to 0.4897. This parameter value yields a canard in each oscillator. \(\varepsilon_1\) and \(\varepsilon_2\) are parameters that describe the capacitors of capacitances \(C_1\) and \(C_2\), respectively. As mentioned above, a 10\(\%\) parameter mismatch is introduced in the values of the capacitances, i.e., we set \(\varepsilon_1=0.1\) and \(\varepsilon_2=0.09\). We use these values of \(B_0\), \(\varepsilon_1\), and \(\varepsilon_2\) throughout this study. Furthermore, we note that \(\alpha\) is a parameter that corresponds to the conductance \(g\); the parameter \(g\) is used as the bifurcation parameter in this section.
We investigate the behavior of the two oscillators and investigate the extent of their synchronization for various values of the parameter \(\alpha\). Figure 6 shows the one-parameter bifurcation diagram for \(\alpha\) decreasing from 1.2 to 0. In Fig. 6, the intersection points between the orbits and Poincaré section (\(x_1 = 1/\sqrt3\), \(\dot{x}_1 \leq 0\)) are plotted as a function of \(\alpha\). As this figure demonstrates, after in-phase synchronization occurs, the in-phase synchronization is maintained over a wide range of \(\alpha\) values.
 |
Fig. 6 One-parameter bifurcation diagram for decreasing values of \(\alpha\). Here, the parameters used are \(\varepsilon_1=0.1\), \(\varepsilon_2 = 0.09\), \(B_0 = 0.4897\), and \(k = 0.1998\). |
The attractor of the in-phase synchronization for \(\alpha\) = 0.4 is shown in Fig. 7(a.1). As \(\alpha\) increases, the in-phase synchronization that more closely resembles complete synchronization emerges, as shown in Fig. 7(b.1). The attractors in both cases take the typical canard form, as shown in Figs. 7(a.2) and (b.2).
 |
Fig. 7 (a) In-phase and (b) almost complete synchronizations of canards in coupled BVP oscillators; in this case the two oscillators include a parameter mismatch of 10% (\(\varepsilon_1=0.1\) and \(\varepsilon_2=0.09\)): Attractors in the (a.1) \(x_1\)-\(x_2\) and (a.2) \(x_1\)-\(y_1\) planes for \(\alpha=0.4\) and attractors in the (b.1) \(x_1\)-\(x_2\) and (b.2) \(x_1\)-\(y_1\) planes for \(\alpha=1.2\). |
In the next section, we undertake circuit experiments for the two identical coupled oscillators to verify the results of the numerical analysis; we verify that in-phase synchronization can be observed across a wide range of values of the bifurcation parameter \(\alpha\). However, notably qualitative differences between the numerical results obtained via Eq. (9) and the experimental measurements presented here exist.
4. Experimental Study of Two Identical Canard-Generating Coupled Oscillators
In the circuit experiment, the conductance \(g\) is used as the bifurcation parameter. First, we investigate the system experimentally for large values of the coupling conductance, \(g\). As expected, for large values of \(g\), when \(g\simeq\) 0.455 mS (\(\alpha \simeq\) 1.2), the oscillators reach an almost complete (pseudo complete) synchronization state, as shown in Fig. 8(a.1). Figure 8(a.2) shows the attractor of the canard shape projected onto the \(v\)-\(i\) phase plane; this finding qualitatively agrees with the results of the numerical simulations.
 |
Fig. 8 Experimental measurement of synchronized canards. Trajectories of the almost-complete (pseudo complete) synchronization in the (a.1) \(v_1\)-\(v_2\) and (a.2) \(v_1\)-\(i_1\) planes in a system with coupling, \(g=445\,\mu\)S. Trajectories in the (b.1) \(v_1\)-\(v_2\) and (b.2) \(v_1\)-\(i_1\) planes in a system with coupling slightly higher than \(g=445\,\mu\)S. (The grid meshes represent 5 V/div in both the horizontal and vertical directions in both (a.1) and (b.1); the grid meshes represent 1 V/div and 5 V/div in the vertical and horizontal directions, respectively, in (a.2) and (b.2).) |
When the coupling parameter \(\alpha\) is decreased slightly, the phase difference between \(v_1\) and \(v_2\) can be clearly observed. At slightly smaller values of \(\alpha\), the in-phase synchronization of the two oscillators collapses, as shown in Fig. 8(b.1). Due to the shape of the Lissajous diagram of the synchronized attractor projected onto the \(v_1\)-\(v_2\) plane, we refer to this complex synchronization behavior as “butterfly synchronization” of canards, a name that reflects the shape of the attractors. The Lissajous diagrams observed here are explained by considering the corresponding numerical results in the next Section.
5. Modeling of a Dynamic Oscillator and the Collapse Mechanism of the in-Phase Synchronization
In the numerical simulation related to Eq. (9), the two nearly identical dynamic oscillators exhibit in-phase synchronization across a wide range of values of the parameter \(\alpha\). However, the experimental measurements indicate that the in-phase synchronization breaks down even for small decreases in the coupling conductance.
To explain the generation of complex synchronization seen in circuit experiments, we investigate the behavior of coupled canard-generating oscillators subject to two weak periodic perturbations, as shown in Fig. 9. Here, we assume that the two weak perturbations have the same angular frequency; this frequency is given by \(\omega\). Their phase angles are given by \(\phi_1\) and \(\phi_2\).
 |
Fig. 9 Resistively coupled BVP oscillators subject to weak periodic perturbations. |
The butterfly synchronization of canards was not observed numerically when considering two sinusoidal forcing waves with the same phase difference. We hypothesize that the dynamics with two sinusoidal perturbations with a moderate phase difference are sufficient to explain the butterfly synchronization of canards. To explain the butterfly synchronization of canards, we consider the following four nonautonomous ODEs:
| \[\begin{equation*} \left\{ \begin{array}{l} C_1\displaystyle\frac{dv_1}{dt} = i_1- G_1(v_1)-g(v_1-v_2), \\ L\displaystyle\frac{di_1}{dt} = E_0 +{E_1\sin(\omega t+\phi_1)}-v_1- i_1R, \\ C_2\displaystyle\frac{dv_2}{dt} = i_2-G_2(v_2)-g(v_2-v_1), \\ L\displaystyle\frac{di_2}{dt} = E_0 + {E_2\sin(\omega t+\phi_2)}-v_2-i_2R, \end{array} \right. \tag{10} \end{equation*}\] |
where the constant DC voltage source \(E_0\) in the second and fourth equations of Eq. (9) is replaced by \(E_0+E_1\sin(\omega t+\phi_1)\) and \(E_0+E_2\sin(\omega t+\phi_2)\), respectively. Here, \(E_1\) and \(E_2\) are assumed to be small, which permits the construction of a dynamic model that explains the phenomena observed experimentally.
In addition to Eq. (8), the following rescaling is used.
| \[\begin{equation*} \begin{array}{@{}l@{}} E_1=\sqrt{\displaystyle\frac{g_1}{g_3}}B_1, {E_2=\sqrt{\displaystyle\frac{g_1}{g_3}}B_2}, \mbox{ and } \omega = {\displaystyle\frac{\nu}{Lg_1}}. \end{array} \tag{11} \end{equation*}\] |
This yields the following normalized equations:
| \[\begin{equation*} \left\{ \begin{array}{l} \varepsilon_1\dot{x}_1 = y_1+x_1-x_1^3-\alpha(x_1-x_2), \\ \dot{y}_1 = -x_1-ky_1+B_0 + {B_1\sin(\nu \tau+\phi_1)}, \\ \varepsilon_2\dot{x}_2 = y_2+x_2-x_2^3+\alpha(x_1-x_2), \\ \dot{y}_2 = -x_2-ky_2+B_0 + {B_2\sin(\nu \tau+\phi_2)}, \end{array} \right. \tag{12} \end{equation*}\] |
where \(B_1\) and \(B_2\) represent the amplitudes, and \(\nu\) represents the angular frequency of the AC forcing terms, respectively. Here, we consider the case where \(B_1 = B_2\), \(\phi_1 = 0\), and \(\phi_2 = \pi/2\) hold; we found that phase difference of \(\pi/2\) between the two perturbations is sufficient to induce the butterfly synchronization of canards in numerical simulations that approximately resembles that observed in the circuit experiments.
The numerical experiments presented here are performed considering the variation of the amplitude and the angular frequency of the forcing terms. For larger values of \(\alpha\) (for example, \(\alpha = 1.2\)), the resistively coupled BVP oscillators with weak periodic perturbations exhibit almost-complete synchronization, as shown in the \(x_1-x_2\) plane of Fig. 10(a). The trajectory in the \(x_i\)-\(y_i\) plane (\(i = 1,~2\)) shows a canard for \(B_1 = 0.15\) and \(\nu_1 = 1.3\). However, for a smaller value of \(\alpha\) (for example, \(\alpha = 0.03\)), the phase difference between the two oscillators increases. As a result, the in-phase synchronization breaks down, as shown in Fig. 10(b). We have succeeded in tuning the set of parameters such that the numerically obtained attractor more closely resembles that obtained in the circuit experiments; this agreement is shown in Fig. 10(b.1).
 |
Fig. 10 (a.1)-(a.3) Numerically obtained almost complete synchronized canards for \(\alpha=1.2, B_1 = 0.15\) and \(\nu_1 = 1.3\). Trajectories projected onto the (a.1) \(x_1\)-\(x_2\), (a.2) \(x_1\)-\(y_1\) (blue) and \(x_2\)-\(y_2\) (green) planes, and (a.3) time-series waveform for \(x_1\) (blue) and \(x_2\) (green). (b.1)-(b.3) Numerically obtained butterfly synchronized canards for \(\alpha=0.03, B_1=0.15\) and \(\nu=1.3\). Trajectories projected onto the (b.1) \(x_1\)-\(x_2\), (b.2) \(x_1\)-\(y_1\) (blue) and \(x_2\)-\(y_2\) (green) planes, and (b.3) time-series waveforms for \(x_1\)(blue) and \(x_2\)(green) of the perturbed system. |
Our model, which is described by the dynamic expressions given in Eq. (12), suggests that weak periodic perturbations decrease the propensity of the system to achieve in-phase synchronizaiton at a lower value of \(\alpha\). Based on the comparison of the experimental and numerical results, we hypothesize that the modelling of coupled BVP dynamics subject to weak perturbations could be suitable to explain the phenomena observed in the circuit experiments. We see that qualitatively similar attractors to those observed in the experiments can be obtained in the numerical results corresponding to coupled oscillators subject to weak periodic perturbations. From Fig. 10(b.3), the butterfly synchronization of canards occurs when the phase difference of the canards is sufficiently large.
6. Conclusion
Herein, we have discussed the synchronization of identical canard-generating BVP oscillators with a parameter mismatch of 10\(\%\) in the capacitances of the coupled oscillators; this mismatch was selected to reflect the experimental errors that are probable to occur in the realized circuits.
We numerically demonstrated that the almost-complete synchronization, which occurs at a coupling parameter \(\alpha\simeq 1.2\), becomes unstable at \(\alpha = 0.4\) and the in-phase synchronization emerges at this parameter value. Converesly, it can be experimentally found that the almost-complete synchronization breaks down for values of the coupling parameter slightly smaller than \(\alpha=1.2\). From the shape of these attractors, we call them butterfly synchronization of canards. This finding could not be replicated in the numerical study of nearly identical coupled canard-generating autonomous ODEs.
To construct a model that explains the experimentally observed phenomena, we introduced nearly identical coupled BVP oscillators subject to two weak periodic perturbations of the same frequency. Considering the resultant nonautonomous dynamics, we showed that the weak perturbations weaken the propensity of the oscillators to reach in-phase synchronization.
In future works, we intend to numerically investigate alternative sets of parameter values with the aim of replicating the phenomena observed in experimental works.
References
[1] M. Diener, “The canard unchained or how fast/slow dynamical systems bifurcate,” The Mathematical Intelligencer, vol.6, pp.38-49, 1984.
CrossRef
[2] A.K. Zvonkin and M.A. Shubin, “Non-standard analysis and singular perturbations of ordinary differential equations,” Russ. Math. Surv., vol.39, no.2, pp.69-131, 1984.
CrossRef
[3] S.M. Baer and T. Erneux, “Singular Hopf bifurcation to relaxation oscillations,” SIAM J. Appl. Math., vol.46, no.5, pp.721-739, 1986.
CrossRef
[4] S.M. Baer and T. Erneux, “Singular Hopf bifurcation to relaxation oscillations II,” SIAM J. Appl. Math., vol.52, no.6, pp.1651-1664, 1992.
CrossRef
[5] B. Braaksma and J. Grasman, “Critical dynamics of the Bonhoeffer-van der Pol equation and its chaotic response to periodic stimulation,” Physica D, vol.68, no.2, pp.265-280, 1993.
CrossRef
[6] V.I. Arnol’d, ed., Encyclopedia of Mathematical Sciences 5, Springer-Verlag, 1994.
[7] J. Guckenheimer, K. Hoffman, and W. Weckesser, “Numerical computation of canards,” Int. J. Bifurc. Chaos, vol.10, no.23, pp.2669-2687, 2000.
CrossRef
[8] M. Desroches and V. Kirk, “Spike-adding in a canonical three time scale model: Superslow explosion & folded-saddle canards,” SIAM J. Appl. Dyn. Syst., vol.17, no.3, pp.1989-2017, 2018.
CrossRef
[9] J.U. Albizuri, M. Desroches, M. Krupa, and S. Rodrigues, “Inflection, canards and folded singularities in excitable systems ― Application to a 3D FitzHugh-Nagumo model,” J. Nonlinear Science, vol.30, pp.3265-3291, 2020.
CrossRef
[10] D. Avitabile, M. Desroches, R. Veltz, and M. Wechselberger, “Local theory for spatio-temporal canards and delayed bifurcations,” SIAM J. Math. Anal., vol.52, no.6, pp.5703-5747, 2020.
CrossRef
[11] M. Desroches, P. Kowalczyk, and S. Rodrigues, “Spike-adding and reset-induced canard cycles in adaptive integrate and fire models,” Nonlinear Dyn., vol.104, pp.2451-2470, 2021.
CrossRef
[12] W. Teka, J. Tabak, and R. Bertram, “The relationship between two fast/slow analysis techniques for bursting oscillations,” Chaos, vol.22, pp.043117-1-11, 2012.
CrossRef
[13] E.K. Ersöz, M. Desroches, M. Krupa, and F. Clément, “Canard-mediated (de)synchronization in coupled phantom bursters,” SIAM J. Appl. Dyn. Syst., vol.15, no.1, pp.580-608, 2016.
CrossRef
[14] E.K. Ersös, M. Desroches, and M. Krupa, “Synchronization of weakly coupled canard oscillators,” Physica D, vol.349, pp.46-61, 2017.
CrossRef
[15] N. Inaba, H. Ito, K. Shimizu, and H. Hikawa, “Complete mixed-mode oscillation synchronization in weakly coupled nonautonomous Bonhoeffer-van der Pol oscillators,” Prog. Theor. Phys., vol.2018, no.6, pp.063A01-1-15, 2018.
CrossRef
[16] S. Fernández-García and A. Vidal, “Symmetric coupling of multiple timescale systems with mixed-mode oscillations and synchronization,” Physica D, vol.401, pp.132129-1-22, 2020.
CrossRef
[17] N.M. Awal and I.R. Epstein, “Period-doubling route to mixed-mode chaos,” Phys. Rev. E, vol.104, pp.024211-1-25, 2021.
CrossRef
[18] K.L. Das, M. Sekikawa, T. Tsubone, N. Inaba, and H. Okazaki, “Experimental and numerical study of the synchronization of canards in identical coupled canard-generating Bonhoeffer-van der Pol oscillators,” Phys. Lett. A, vol.465, pp.128709-1-5, 2023.
CrossRef
[19] M. Krupa and P. Szmolyan, “Relaxation oscillation and canard explosion,” J. Differ. Equ., vol.174, no.2, pp.312-368, 2001.
CrossRef
Authors
Utsunomiya University
received a B.E. degree from the Institute of Engineering, Tribhuwan University, Kathmandu, Nepal in 2009 and an M.E. degree in Mechanical and Intelligent Engineering from Utsunomiya University, Utsunomiya, Japan in 2022. He is a Ph.D. candidate at Utsunomiya University and his research interests are in nonlinear dynamical systems.
 |
Utsunomiya University
received B.E. and M.E. degrees in information science and the Ph.D. degree in production and information science from Utsunomiya University, Utsunomiya, Japan, in 2000, 2002, and 2005, respectively. In 2013, he joined the Department of Mechanical and Intelligent Engineering at Utsunomiya University, Utsunomiya, Japan, where he is currently an Associate Professor. His research interests are in nonlinear dynamics, synchronization phenomena, and chaos.
 |
Nagaoka University of Technology
received B.E., M.E., and Ph.D. degrees in electrical engineering from Hosei University, Tokyo, Japan, in 1996, 1998, and 2001, respectively. He is currently Professor with the Department of Electrical, Electronics and Information Engineering, Nagaoka University of Technology, Niigata, Japan. His research interests include neural networks, chaos, and bifurcation.
 |
Shonan Institute of Technology
received B.E., M.E., and Ph.D. degrees in electrical engineering from Keio University in 1984, 1986, and 1989, respectively. He is currently a guest researcher at Shonan Institute of Technology. His research interests include bifurcations and chaos in nonlinear dynamical circuits.
 |
Shonan Institute of Technology
received his B.E., M.E., and Ph.D. degrees in mechanical engineering from Waseda University in 1984, 1986, and 1996, respectively. He was a research associate at Waseda University from 1989 to 1992, an assistant professor of the Department of Electronic Control Engineering, Gifu National College of Technology, from 1992 to 1996, and an associate professor of the Department of Electronic Control Engineering, Gifu National College of Technology, from 1996 to 2000. He was engaged in applied mathematics as a visiting scientist at Cornell University, USA from 1999 to 2001. He was made a professor in the Department of System and Communication Engineering, Shonan Institute of Technology (SIT), from 2001 to 2005, and a professor in the Department of Applied Computer Sciences at SIT, from 2001 to 2023. From 2023 to the present, he has been a professor in the Department of Informatics, Computer Science Course at SIT, and from 2013 to the present, he has been a professor in the Graduate School of Engineering, Electrical and Information Engineering at SIT. From 2019 to the present, He has been Director of Internet of Things-Inspired Interdisciplinary Research Center (IoT-IIRC) at SIT, and he also has been Director of Research Center for Quantum Computing (RCQC) at SIT, from 2022 to the present. His research interests include nonlinear circuit and system analysis and control and methods of rigorous computer-assisted proof for nonlinear circuits and systems.
 |
