Math. Model. Nat. Phenom. 14 (2019) 409 Mathematical Modelling of Natural Phenomena https://doi.org/10.1051/mmnp/2019012 www.mmnp-journal.org CANARDS EXISTENCE IN THE HINDMARSH–ROSE MODEL Jean-Marc Ginoux1, Jaume Llibre2, * and Kiyoyuki Tchizawa3 Abstract. In two previous papers we have proposed a new method for proving the existence of “canard solutions” on one hand for three and four-dimensional singularly perturbed systems with only one fast variable and, on the other hand for four-dimensional singularly perturbed systems with two fast variables [J.M. Ginoux and J. Llibre, Qual. Theory Dyn. Syst. 15 (2016) 381–431; J.M. Ginoux and J. Llibre, Qual. Theory Dyn. Syst. 15 (2015) 342010]. The aim of this work is to extend this method which improves the classical ones used till now to the case of three-dimensional singularly perturbed systems with two fast variables. This method enables to state a unique generic condition for the existence of “canard solutions” for such three-dimensional singularly perturbed systems which is based on the stability of folded singularities (pseudo singular points in this case) of the normalized slow dynamics deduced from a well-known property of linear algebra. Applications of this method to a famous neuronal bursting model enables to show the existence of “canard solutions” in the Hindmarsh-Rose model. Mathematics Subject Classification. 34C23, 34C25, 37G10 Received October 23, 2018. Accepted February 13, 2019. 1. Introduction The concept of “canard solutions” for three-dimensional singularly perturbed systems with two slow vari- ables and one fast has been introduced in the beginning of the eighties by Benoı̂t and Lobry [5], Benoı̂t [4]. Their existence has been proved by Benoı̂t ([4], p. 170) in the framework of “Non-Standard Analysis” according to a theorem which states that canard solutions exist in such systems provided that the pseudo singular point of the slow dynamics, i.e., of the reduced vector field is of saddle type. Nearly twenty years later, while using the so-called “blow-up” technique they introduced, Dumortier and Roussarie [6] and then, Szmolyan and Wechselberger [31] provided a “standard version” of Benoı̂t’s theorem [4]. Recently, Wechsel- berger [37] generalized this theorem for n-dimensional singularly perturbed systems with k slow variables and m fast (where n = k + m). The method they used require to implement a “desingularization proce- dure” which can be summarized as follows: first, they compute the normal form of such singularly perturbed systems which is expressed according to some coefficients (a and b for dimension three and ã, b̃ and c̃1 for Keywords and phrases: Hindmarsh–Rose model, singularly perturbed dynamical systems, canard solutions. 1 Laboratoire d’Informatique et des Systèmes, UMR, CNRS 7020, Université de Toulon, BP 20132, 83957 La Garde cedex, France. 2 Departament de Matemàtiques, Universitat Autònoma de Barcelona, 08193 Bellaterra, Barcelona, Catalonia, Spain. 3 Institute of Administration Engineering, Ltd., Tokyo 101-0021, Japan. * Corresponding author: jllibre@mat.uab.cat Article published by EDP Sciences c EDP Sciences, 2019 2 J.-M. GINOUX ET AL. dimension four) depending on the functions defining the original vector field and their partial derivatives with respect to the variables. Secondly, they project the “desingularized vector field” (originally called “normal- ized slow dynamics” by Eric Benoı̂t [4], p. 166) of such a normal form on the tangent bundle of the critical manifold. Finally, they evaluate the Jacobian of the projection of this “desingularized vector field” at the folded singularity (originally called pseudo singular points by José Argémi [1], p. 336). This lead Szmolyan and Wechselberger ([31], p. 427) and Wechselberger ([37], p. 3298) to a “classification of folded singularities (pseudo singular points)”. Thus, they showed that for three-dimensional (resp. four-dimensional) singularly perturbed systems such folded singularity is of saddle type if the following condition is satisfied: a < 0 (resp. ã < 0). In a first paper entitled: “Canards Existence in Memristor’s Circuits” (see Ginoux and Llibre [16]) we pre- sented a method enabling to state a unique “generic” condition for the existence of “canard solutions” for three and four-dimensional singularly perturbed systems with only one fast variable which is based on the stability of folded singularities of the normalized slow dynamics deduced from a well-known property of linear algebra. We proved that this unique condition is completely identical to that provided by Benoı̂t [4], Szmolyan and Wechselberger [31] and Wechselberger [37]. In a second paper entitled: “Canards Existence in FitzHugh-Nagumo and Hodgkin-Huxley Neuronal Mod- els” (see Ginoux and Llibre [15]) we extended this method to the case of four-dimensional singularly perturbed systems with k = 2 slow and m = 2 fast variables. Then, we stated that the provided condition for the existence of canards is “generic” since it is exactly the same for singularly perturbed systems of dimension three and four with one or two fast variables. The method we used led us to the following proposition: If the normal- ized slow dynamics has a pseudo singular point of saddle type, i.e. if the sum σ2 of all second-order diagonal minors of the Jacobian matrix of the normalized slow dynamics evaluated at the pseudo singular point is neg- ative, i.e. if σ2 < 0 then, the three-dimensional (resp. four-dimensional) singularly perturbed system exhibits a canard solution which evolves from the attractive part of the slow manifold towards its repelling part. Then, we proved on one hand for three-dimensional singularly perturbed systems with only one fast variable that the condition for which the pseudo singular point is of saddle type, i.e. σ2 < 0 is identical to that proposed by Benoı̂t ([4], p. 171) in his theorem, i.e. D < 0 and also to that provided by Szmolyan and Wechselberger [31], i.e. a < 0. On the other hand, we proved for four-dimensional singularly perturbed systems with one or two fast variables that the condition for which the folded singularity (resp. the pseudo singular point) is of saddle type, i.e. σ2 < 0 is identical to that proposed by Wechselberger ([37], p. 3298) in his theorem, i.e. ã < 0. Notice that there is no proof of the approximation. It is not established that the time-scaled reduced system holds on the approximation for the original system in the case of k slow variables (k ≥ 3), m fast variables (m ≥ 2). It was proved in the case k = 2 and m = 1 by Benoı̂t; constructing a local model and obtaining its solution, and in the case k = 2 and m = 2 was also proved extensively by Tchizawa [33, 34]. For the case k = 1 and m = 2 (Hindmarsch–Rose model), we shall construct a local model again and we shall obtain their solutions, providing a constructive proof for the approximation (see also the seminal works of Belikov and Samborskii [3], Gol’dshtein et al. [17], Sobolev and Shchepakina [30] and Kolesov and Rozov [25]). Being the pseudo-singular point a saddle, or a node it does not ensure the existence of canards, because it may not satisfy the approximation. The aim of this work is to extend this method to the case of three-dimensional singularly perturbed systems with one slow and two fast variables and to show that the provided condition for the existence of canards, i.e. σ2 < 0 still holds and is consequently “generic”. The outline of this paper is as follows. In Section 2, definitions of singularly perturbed system, critical manifold, reduced system, “constrained system”, canard cycles, folded singularities and pseudo singu- lar points are recalled. The method proposed in this article is presented in Section 3 for the case of three-dimensional singularly perturbed systems with two fast variables. In Section 4, applications of this method to the famous Hindmarsh-Rose model enables to show the existence of “canard solutions” in such system. CANARDS EXISTENCE IN THE HINDMARSH–ROSE MODEL 3 2. Definitions 2.1. Singularly perturbed systems According to Tikhonov [35], Jones [22] and Kaper [23] singularly perturbed systems (or slow-fast) are defined as: ~x0 = εf~ (~x, ~y , ε) , (2.1) ~y 0 = ~g (~x, ~y , ε) , where ~x ∈ Rk , ~y ∈ Rm , ε ∈ R+ , and the prime denotes differentiation with respect to the independent variable t0 . The functions f~ and ~g are assumed to be C ∞ functions1 of ~x, ~y and ε in U × I, where U is an open subset of Rk × Rm and I is an open interval containing ε = 0. In the case when 0 < ε  1, i.e. ε is a small positive number, the variable ~x is called slow variable, and ~y is called fast variable. Using Landau’s notation2 : O (εp ) represents a function f of u and ε such that f (u, ε)/εp is bounded for positive ε going to zero, uniformly for u in the given domain. In general we consider that ~x evolves at an O (ε) rate; while ~y evolves at an O (1) slow rate. Reformulating system (2.1) in terms of the rescaled variable t = εt0 , we obtain ~x˙ = f~ (~x, ~y , ε) , (2.2) ε~y˙ = ~g (~x, ~y , ε) . The dot represents the derivative with respect to the new independent variable t. The independent variables t0 and t are referred to the fast and slow times, respectively, and (2.1) and (2.2) are called the fast and slow systems, respectively. These systems are equivalent whenever ε 6= 0, and they are labeled singular perturbation problems when 0 < ε  1. The term “singular perturbation” was introduced by Friedrichs and Wasow [12] and the label “singular” stems in part from the discontinuous limiting behavior in system (2.1) as ε → 0. 2.2. Reduced slow system When ε → 0 system (2.2) leads to a system of differential-algebraic equations (D.A.E.) called reduced slow system whose dimension decreases from k + m = n to m. Then, the slow variable ~x ∈ Rk partially evolves in the submanifold n o M0 := (~x, ~y ) : ~g (~x, ~y , 0) = ~0 . (2.3) called the critical manifold3 . The reduced slow system is ~x˙ = f~ (~x, ~y , ε) , (2.4) ~0 = ~g (~x, ~y , ε) . 2.3. Slow invariant manifold Such a normally hyperbolic invariant manifold (2.3) of the reduced slow system (2.4) persists as a locally invariant slow manifold of the full system (2.1) for ε sufficiently small. Following to Kaper [23]: “The definition of normal hyperbolicity for general systems involves verifying that the exponential growth rates of orbits in the 1 In certain applications these functions will be supposed to be C r , r > 1. 2 Edmund Georg Hermann Landau (1877–1938) (not to be confused with Lev Landau (1908–1968), Nobel Prize in Physics in 1962) specifies having adopted in 1909 the symbol O introduced by Paul Bachmann (1837–1920) in 1894. 3 It represents the approximation of the slow invariant manifold, with an error of O(ε). 4 J.-M. GINOUX ET AL. directions normal to M0 are stronger than the growth rates in the tangential directions, as may be measured, for example, using Lyapunov type numbers.” (For more details see Zadiraka [39, 40], Baris [2], Fenichel [7], Hirsch et al. [20], Knobloch [24] and Wiggins [38]). This locally slow invariant manifold is O(ε) close to the critical manifold. When D~x f~ is invertible, thanks to the Implicit Function Theorem, M0 is given by the graph of a C ∞ function ~x = G~ 0 (~y ) for ~y ∈ D, where D ⊆ Rk is a compact, simply connected domain and the boundary of D is a (k − 1)–dimensional C ∞ submanifold4 . ~ (~y , ε) According to Fenichel [7, 10] theory if 0 < ε  1 is sufficiently small, then there exists a function G defined on D such that the manifold n o ~ (~y , ε) , Mε := (~x, ~y ) : ~x = G (2.5) is locally invariant under the flow of system (2.1). Moreover, there exist perturbed local stable (or attracting) M a and unstable (or repelling) M r branches of the slow invariant manifold Mε . Thus, normal hyperbolicity of Mε is lost via a saddle-node bifurcation of the reduced slow system (2.4). Then, it gives rise to solutions of “canard” type. 2.4. Canards, singular canards and maximal canards A canard is a solution of a singularly perturbed dynamical system (2.1) following the attracting branch M a of the slow invariant manifold, passing near a bifurcation point located on the fold of this slow invariant manifold, and then following the repelling branch M r of the slow invariant manifold. A singular canard is a solution of a reduced slow system (2.4) following the attracting branch M0a of the critical manifold, passing near a bifurcation point located on the fold of this critical manifold, and then following the repelling branch M0r of the critical manifold. A maximal canard corresponds to the intersection of the attracting and repelling branches Mεa ∩ Mεr of the slow manifold in the vicinity of a non-hyperbolic point. According to Wechselberger ([37], p. 3302): “Such a maximal canard defines a family of canards nearby which are exponentially close to the maximal canard, i.e. a family of solutions of (2.1) that follow an attracting branch Mεa of the slow manifold and then follow, rather surprisingly, a repelling/saddle branch Mεr of the slow manifold for a considerable amount of slow time. The existence of this family of canards is a consequence of the non-uniqueness of Mεa and Mεr . However, in the singular limit ε → 0, such a family of canards is represented by a unique singular canard.” Canards are a special class of solution of singularly perturbed dynamical systems for which normal hyper- bolicity is lost. Canards in singularly perturbed systems with two or more slow variables (~x ∈ Rk , k > 2) and one fast variable (~y ∈ Rm , m = 1) are robust, since maximal canards generically persist under small parameter changes5 . 2.5. Constrained system In order to characterize the “slow dynamics”, i.e. the slow trajectory of the reduced slow system (2.4) (obtained by setting ε = 0 in (2.2)), Floris Takens [32] introduced the “constrained system” defined as follows: ~x˙ = f~ (~x, ~y , 0) , D~y~g .~y˙ = −(D~x~g .f~) (~x, ~y , 0) , (2.6) ~0 = ~g (~x, ~y , 0) . 4 The set D is overflowing invariant with respect to (2.2) when ε = 0. See Kaper [23] and Jones [22]. 5 See Benoı̂t [4], Szmolyan and Wechselberger [31] and Wechselberger [36, 37]. CANARDS EXISTENCE IN THE HINDMARSH–ROSE MODEL 5 Since, according to Fenichel [7, 10], the critical manifold ~g (~x, ~y , 0) may be considered as locally invariant under the flow of system (2.1), we have: d~g (~x, ~y , 0) = 0 ⇐⇒ D~x~g .~x˙ + D~y~g .~y˙ = ~0. dt By replacing ~x˙ by f~ (~x, ~y , 0) leads to: D~x~g .f~ (~x, ~y , 0) + D~y~g .~y˙ = ~0. This justifies the introduction of the constrained system. Now, let adj(D~y~g ) denote the adjoint of the matrix D~y~g which is the transpose of the co-factor matrix D~y~g , then while multiplying the left hand side of (2.6) by the inverse matrix (D~y~g )−1 obtained by the adjoint method we have: ~x˙ = f~ (~x, ~y , 0) , det(D~y~g )~y˙ = −(adj(D~y~g ).D~x~g .f~) (~x, ~y , 0) , (2.7) ~0 = ~g (~x, ~y , 0) . 2.6. Normalized slow dynamics By rescaling the time by setting t = −det(D~y~g )τ in system (2.7) we obtain the following system which has been called by Eric Benoı̂t ([4], p. 166) “normalized slow dynamics”: ~x˙ = −det(D~y~g )f~ (~x, ~y , 0) , ~y˙ = (adj(D~y~g ).D~x~g .f~) (~x, ~y , 0) , (2.8) ~0 = ~g (~x, ~y , 0) . where the overdot now denotes the time derivation with respect to τ . We notice that Argémi [1] proposed to rescale the time by setting t = −det(D~y~g )sgn(det(D~y~g ))τ in order to keep the same flow direction in (2.8) as in (2.7). 2.7. Desingularized vector field By application of the Implicit Function Theorem, we suppose that we can explicitly express from equa- tion (2.3), say without loss of generality, x1 as a function φ1 of the other variables. This implies that M0 is locally the graph of a function φ1 : Rk → Rm over the base U = (~ χ, ~y ) where χ ~ = (x2 , x3 , ..., xk ). Thus, we can span the “normalized slow dynamics” on the tangent bundle at the critical manifold M0 at the pseudo singular point. This leads to the so-called desingularized vector field: ~˙ = −det(D~y~g )f~ (~ χ χ, ~y , 0) , (2.9) ~y˙ = (adj(D~y~g ).D~x~g .f~) (~ χ, ~y , 0) . 2.8. Pseudo singular points and folded singularities As recalled by Guckenheimer and Haiduc ([18], p. 91), pseudo-singular points have been introduced by the late José Argémi [1] for low-dimensional singularly perturbed systems and are defined as singular points of the “nor- malized slow dynamics” (2.8). Twenty-three years later, Szmolyan and Wechselberger ([31], p. 428) called such pseudo singular points, folded singularities. In a recent publication entitled “A propos de canards” Wechselberger 6 J.-M. GINOUX ET AL. ([37], p. 3295) proposed to define such singularities for n-dimensional singularly perturbed systems with k slow variables and m fast as the solutions of the following system: det(D~y~g ) = 0, (adj(D~y~g ).D~x~g .f~) (~x, ~y , 0) = ~0, (2.10) ~g (~x, ~y , 0) = ~0. Thus, for dimensions higher than three, his concept encompasses that of Argémi. Moreover, for k > 2), Wechselberger ([37], p. 3296) proved that folded singularities form a (k − 2)-dimensional manifold. Thus, for k = 2 the folded singularities are nothing else than the pseudo singular points defined by Argémi [1]. Neverthe- less, for the degenerate case k = 1, folded singularities still form a zero-dimensional manifold. So, we will see in Section 3 that the stability analysis of the pseudo singular points will give rise to a condition for the existence of canard solutions in the original system (2.1). 3. Three-dimensional singularly perturbed systems with two fast variables A three-dimensional singularly perturbed dynamical system (2.2) with k = 1 slow variables and m = 2 fast may be written as: ẋ1 = f1 (x1 , y1 , y2 ) , (3.1a) εẏ1 = g1 (x1 , y1 , y2 ) , (3.1b) εẏ2 = g2 (x1 , y1 , y2 ) , (3.1c) where x1 ∈ R, ~y = (y1 , y2 )t ∈ R2 , 0 < ε  1 and the functions fi and gi are assumed to be C 2 functions of (x1 , y1 , y2 ). 3.1. Critical manifold The critical manifold equation of system (3.1) is defined by setting ε = 0 in equations (3.1b) and (3.1c). Thus, we obtain: g1 (x1 , y1 , y2 ) = 0, (3.2a) g2 (x1 , y1 , y2 ) = 0. (3.2b) By application of the Implicit Function Theorem, we suppose that we can explicitly express from equa- tions (3.2a) and (3.2b), say without loss of generality, x1 and y1 as functions of the others variables: y1 = φ1 (x1 , y2 ) , (3.3a) x1 = φ2 (y1 , y2 ) . (3.3b) 3.2. Constrained system The constrained system is obtained by equating to zero the time derivative of g1,2 (x1 , y1 , y2 ): dg1 ∂g1 ∂g1 ∂g1 = ẋ1 + ẏ1 + ẏ2 = 0, (3.4a) dt ∂x1 ∂y1 ∂y2 CANARDS EXISTENCE IN THE HINDMARSH–ROSE MODEL 7 dg2 ∂g2 ∂g2 ∂g2 = ẋ1 + ẏ1 + ẏ2 = 0. (3.4b) dt ∂x1 ∂y1 ∂y2 Equations (3.4a) and (3.4b) may be written as: ∂g1 ∂g1 ∂g1 ẏ1 + ẏ2 = − ẋ1 , (3.5a) ∂y1 ∂y2 ∂x1 ∂g2 ∂g2 ∂g2 ẏ1 + ẏ2 = − ẋ1 . (3.5b) ∂y1 ∂y2 ∂x1 By solving the system of two equations (3.5a) and (3.5b) with two unknowns (ẏ1 , ẏ2 ) we find:   ∂g1 ∂g2 ∂g1 ∂g2 − ẋ1 ∂x1 ∂y2 ∂y2 ∂x1 ẏ1 = −   , (3.6a) det J(y1 ,y2 )   ∂g1 ∂g2 ∂g1 ∂g2 − ẋ1 ∂y1 ∂x1 ∂x1 ∂y1 ẏ2 = −   . (3.6b) det J(y1 ,y2 ) So, we have the following constrained system: ẋ1 = f1 (x1 , y1 , y2 ) ,   ∂g1 ∂g2 ∂g1 ∂g2 − ẋ1 ∂x1 ∂y2 ∂y2 ∂x1 ẏ1 = −   , det J(y1 ,y2 )   ∂g1 ∂g2 ∂g1 ∂g2 (3.7) − ẋ1 ∂y1 ∂x1 ∂x1 ∂y1 ẏ2 = −   , det J(y1 ,y2 ) 0 = g1 (x1 , y1 , y2 ) , 0 = g2 (x1 , y1 , y2 ) . 3.3. Normalized slow dynamics   By rescaling the time by setting t = −det J(y1 ,y2 ) τ we obtain the “normalized slow dynamics”:   ẋ1 = −f1 (x1 , y1 , y2 ) det J(y1 ,y2 ) = F1 (x1 , y1 , y2 ) ,   ∂g1 ∂g2 ∂g1 ∂g2 ẏ1 = f1 (x1 , y1 , y2 ) − = G1 (x1 , y1 , y2 ) , ∂x1 ∂y2 ∂y2 ∂x1   ∂g1 ∂g2 ∂g1 ∂g2 (3.8) ẏ2 = f1 (x1 , y1 , y2 ) − = G2 (x1 , y1 , y2 ) , ∂y1 ∂x1 ∂x1 ∂y1 0 = g1 (x1 , y1 , y2 ) , 0 = g2 (x1 , y1 , y2 ) , where the overdot now denotes the time derivation with respect to τ . 8 J.-M. GINOUX ET AL. 3.4. Desingularized system on the critical manifold Since we have supposed that y1 and y2 may be explicitly expressed as functions of the others variables (3.3a) and (3.3b), they can be used to project the normalized slow dynamics (3.8) on the tangent bundle of the critical manifold. So, we have:   ẋ1 = −f1 (x1 , y1 , y2 ) det J(y1 ,y2 ) = F1 (x1 , y2 ) , (3.9)   ∂g1 ∂g2 ∂g1 ∂g2 ẏ2 = f1 (x1 , y1 , y2 ) − = G2 (x1 , y2 ) . ∂y1 ∂x1 ∂x1 ∂y1 3.5. Pseudo singular points Pseudo-singular points are defined as singular points of the “normalized slow dynamics”, i.e. as the set of points for which we have:   det J(y1 ,y2 ) = 0, (3.10a)   ∂g1 ∂g2 ∂g1 ∂g2 − = 0, (3.10b) ∂x1 ∂y2 ∂y2 ∂x1   ∂g1 ∂g2 ∂g1 ∂g2 − = 0, (3.10c) ∂y1 ∂x1 ∂x1 ∂y1 g1 (x1 , y1 , y2 ) = 0, (3.10d) g2 (x1 , y1 , y2 ) = 0. (3.10e) Remark. We notice on the one hand that equations (3.10b) and (3.10c) are linearly dependent and on the other hand that contrary to previous works we do not use the “desingularized vector field” (3.9) but the “normalized slow dynamics” (3.8). The Jacobian matrix of system (3.8) reads: ∂F1 ∂F1 ∂F1    ∂x1 ∂y1 ∂y2     ∂G ∂G1 ∂G1   1  J(F1 ,G1 ,G2 ) =  ∂x1 ∂y1 ∂y2  .  (3.11)    ∂G2 ∂G2 ∂G2    ∂x1 ∂y1 ∂y2 3.6. Extension of Benoı̂t’s generic hypothesis Without loss of generality, it seems reasonable to extend Benoı̂t’s generic hypotheses introduced for three- dimensional singularly perturbed systems with only one fast variable to the case of two variables. So, first, we suppose that by a “standard translation” the pseudo singular point can be shifted at the origin O(0, 0, 0) and that by a “standard rotation” of y1 -axis that the slow manifold is tangent to (x1 , y1 )-plane, so we have g1 (0, 0, 0) = g2 (0, 0, 0) = 0, (3.12a) ∂g1 ∂g2 = = 0, (3.12b) ∂y1 (0,0,0) ∂y1 (0,0,0) ∂g1 ∂g2 = = 0. (3.12c) ∂y2 (0,0,0) ∂y2 (0,0,0) CANARDS EXISTENCE IN THE HINDMARSH–ROSE MODEL 9 Then, we make the following assumptions for the non-degeneracy of the pseudo singular point: f1 (0, 0, 0) 6= 0. (3.13) Thus, we have the following Cayley-Hamilton eigenpolynomial associated with such a Jacobian matrix (3.11) evaluated at the pseudo singular point, i.e. at the origin: λ3 − σ1 λ2 + σ2 λ − σ3 = 0. (3.14) First, it can be proved that the sum of all first-order diagonal minors of J(F1 ,G1 ,G2 ) , i.e. σ1 = T r(J(F1 ,G1 ,G2 ) ) = 0 due to circular permutations of the partial derivatives. Secondly, it can also be proved that σ3 = |J(F1 ,G1 ,G2 ) | = 0 vanishes at a pseudo singular point provided that generic condition (3.12b) is satisfied. So, the eigenpolynomial (3.14) is reduced to λ λ2 + σ2 = 0  (3.15) Let λi be the eigenvalues of the eigenpolynomial (3.15) and we denote by λ3 = 0 the obvious root of this polynomial. We have: 3 X ii σ2 = J(F 1 ,G1 ,G2 ) = λ1 λ2 . (3.16) i=1 P3 ii where σ2 = i=1 J(F 1 ,G1 ,G2 ) = q represents the sum of all second-order diagonal minors of J(F1 ,G1 ,G2 ) . Obviously, the pseudo singular point is of saddle-type if and only if σ2 < 0. This leads to the the following condition: C1 : q < 0. (3.17) 3.7. Canard existence in R1+2 In an article entitled “Systèmes lents-rapides dans R3 et leurs canards”, Benoı̂t ([4], p. 171) has stated in the framework of “non-standard analysis” a theorem that can be written as follows: Benoı̂t’s theorem [1983]. If the desingularized vector field (3.7) has a pseudo singular point of saddle type, then system (3.1) exhibits a canard solution which evolves from the attractive part of the slow manifold towards its repelling part. A few years later, Szmolyan and Wechselberger [31] gave a “standard version” of Benoı̂t’s theorem [4] (see Benoı̂t’s theorem above) for three-dimensional singularly perturbed systems with k = 2 slow variables and m = 1 fast. While using “standard analysis” and blow-up technique they introduced, Dumortier and Roussarie [6], and then, Szmolyan and Wechselberger ([31], p. 427) stated in their Lemma 2.1, while using “a smooth change of coordinates” (see Ginoux et al. [15, 16]), that the original system can be transformed into a “normal form” (4.1) from which they deduced that the condition for the pseudo singular point to be of saddle type is a < 0. Then, they proved the existence of canard solutions for the original system according to their Theorem 4.1(a). In our previous papers (see Ginoux et al. [15, 16]), we have established the following Proposition 3.4 for three and four-dimensional singularly perturbed systems with k = 1 fast variable and for four-dimensional singularly perturbed systems with k = 2 fast variables. In this work we show that this proposition still holds for three- dimensional singularly perturbed systems with k = 2 fast variables and m = 1 slow variable. (See also the seminal works of Belikov and Samborskii [3], Gol’dshtein et al. [17], Sobolev and Shchepakina [30] and Kolesov and Rozov [25].) 10 J.-M. GINOUX ET AL. As we have pointed out the approximation regarding the reduced system is done constructing a local model in the case k = 1 and m = 3. Let the origin O be a saddle or a node. Blowing up the coordinates x = (x1 , x2 , x3 ) ∈ R3 and y ∈ R as follows: x1 = α2 u1 , x2 = α2 u2 , x3 = αu3 and y = α2 v with α ≈ 0, the system is reduced to du1 1 ε = h1 (x, y, ε), dt α2 du2 1 ε = 2 h2 (x, y, ε), dt α (3.18) du3 1 ε = h3 (x, y, ε), dt α dv 1 = 2 f (x, y, ε). dt α Assume that the rank of the matrix (hx ) is 3, i.e. there exists a function x = ξ(y) which is invertible. When ε = 0 the time-scaled reduced system is dx dy = − det(hx )(hx )−1 (hy ) , dt dt dy = f (x, y, 0). dt Doing the scaling t = α2 τ and ε/α2 ≈ 0 we obtain the local model 1 δ u̇1 = h1 (x, y, ε) α2 h1 (0) ∂h1 ∂h1 ∂h1 ∂h1 1 ∂ 2 h1 = + (0)u1 + (0)u2 + (0)u3 + (0)v + (0)u23 + · · · , α2 ∂x1 ∂x2 ∂x3 ∂y 2 ∂x23 1 δ u̇2 = h2 (x, y, ε) α2 h2 (0) ∂h2 ∂h2 ∂h2 ∂h2 1 ∂ 2 h2 = + (0)u1 + (0)u2 + (0)u3 + (0)v + (0)u23 + · · · , α2 ∂x1 ∂x2 ∂x3 ∂y 2 ∂x23 1 δ u̇3 = h3 (x, y, ε) α (3.19) h3 (0) ∂h3 ∂h3 ∂h3 ∂h3 = +α (0)u1 + α (0)u2 + (0)u3 + α (0)v + · · · , α ∂x1 ∂x2 ∂x3 ∂y v̇ = f (x, y, ε) ∂f ∂f ∂f ∂f = f (0) + α2 (0)u1 + α2 (0)u2 + α (0)u3 + α2 (0)v + · · · , ∂x1 ∂x2 ∂x3 ∂y where here the dot denotes derivative with respect to τ , and δ = ε/α2 ≈ 0. Now we denote ∂hi ∂ 2 hi aij = (0) and bij = (0). ∂xj ∂x2j CANARDS EXISTENCE IN THE HINDMARSH–ROSE MODEL 11 Then, under the assumptions a13 = a23 = 0, a33 < 0 f (0) 6= 0, and h(0) = 0, and taking δ = 0 we have the system a13 1 a11 u1 + a12 u2 + u3 + a14 v + b13 u23 = 0, α 2 a23 1 a21 u1 + a22 u2 + u3 + a24 v + b23 u23 = 0, α 2 a33 u3 = 0, v = f (0)t. Solving this system we get 1 u1 = u01 = (a12 a24 − a14 a22 )f (0)t, a11 a22 − a12 a21 1 u2 = u02 = (a14 a21 − a11 a24 )f (0)t, a11 a22 − a12 a21 u3 = u03 = 0. Now the solution with δ ≈ 0 is of the form u1 = u01 + L(δ), u2 = u02 + L(δ), u3 = u03 + L(δ). Lemma 3.1. The time-scaled-reduced system in R1+3 is an approximation of system (3.18). If the singular point of the reduced system is stable, a33 < 0 and the trace[∂h(0, 0)/∂x] < 0, then there exists a canard in the system R1+3 . For the system in R1+2 taking x2 = 0 and h2 (x, y, ε) = 0, then x = (x1 , x3 ) ∈ R2 and y ∈ R and system (3.18) becomes du1 1 ε = 2 h1 (x, y, ε), dt α du3 1 ε = h3 (x, y, ε), (3.20) dt α dv 1 = 2 f (x, y, ε), dt α and consequently the local model given by system (3.19) now writes 1 a13 1 δ u̇1 = h1 (x, y, ε) = a11 u1 + u3 + a14 v + b13 u23 , α2 α 2 1 (3.21) δ u̇3 = h3 (x, y, ε) = a33 u3 , α v̇ = f (x, y, ε) = f (0). 12 J.-M. GINOUX ET AL. Therefore, under the assumptions a13 = a23 = 0, a33 < 0 f (0) 6= 0, and h(0) = 0, taking δ = 0 we obtain the system a11 u1 + a14 v = 0, a33 u3 = 0, v = f (0)t. Solving this system we get the solution u1 = u01 = −a14 f (0)t, u3 = u03 = 0, v = f (0)t. Hence the solution with δ ≈ 0 is of the form u1 = u01 + L(δ), u3 = u03 + L(δ). Lemma 3.2. The time-scaled-reduced system in R3+1 is an approximation for the system (3.20). If the pseudo singular point is a saddle (all the eigenvalues are negative in the slow vector field), or a node (all the eigenvalues are positive) and a11 < 0 (or a11 > 0) and f (0) 6= 0, then there exists a canard in the system R3+1 . Now we use Theorems 4.1 and 4.2 from [33] for the following model in R2+2 : εẋ = h(x, y, ε), ẏ = f (x, y, ε), where x = (x1 , x2 ), y = (y1 , y2 ), h = (h1 , h2 ) and f = (f1 , f2 ). Assume that the origin of this system is a singular point, and that the slow variable y2 ≡ 0 and the corresponding function f2 ≡ 0. Then we obtain the following system in R1+2 , i.e. with k = 1 and m = 2. Blowing up the variables as x1 = α2 u1 , x2 = αu2 and y1 = α2 v1 , we write the system as 1 εu̇1 = h1 (x, y, ε), α2 1 εu̇2 = h2 (x, y, ε), (3.22) α 1 v̇1 = 2 f1 (x, y, ε). α CANARDS EXISTENCE IN THE HINDMARSH–ROSE MODEL 13 Changing the time scale as t = α2 τ this system becomes 1 δ u̇1 = h1 (x, y, ε), α2 1 δ u̇2 = h2 (x, y, ε), α 1 v̇1 = 2 f1 (x, y, ε), α where δ = ε/α2 is very small and now the dot denotes derivative with respect to τ . Under the assumptions (1) a12 = 0 and a11 6= 0, (2) f1 (0, 0, 0) 6= 0, and (3) there exists ϕ(x) = y1 such that the trace of the Jacobian matrix hx (0, ϕ(0), 0) is negative. Then the last system can be written into the form 1 δ u̇1 = a11 u1 + a12 u2 + a13 v1 + b12 u22 + L(δ), 2 δ u̇2 = a21 u2 + L(δ), 1 v̇1 = f1 (0, 0, 0), α2 Taking δ = 0 we get that a13 u1 = − f1 (0, 0, 0)t, a11 u2 = 0, v1 = f1 (0, 0, 0)t. Therefore, for δ 6= 0 but sufficiently small we have a13 u1 = − f1 (0, 0, 0)t + L(δ), a11 u2 = +L(δ), v1 = f1 (0, 0, 0)t + L(δ). Lemma 3.3. Using the case R2+2 , and taking y2 = 0 and f2 = 0, the time-scaled reduced system has the approximation of system (3.22). Proposition 3.4. If the normalized slow dynamics (3.8) has a pseudo singular point of saddle type, i.e. if the sum σ2 of all second-order diagonal minors of the Jacobian matrix of the normalized slow dynamics (3.8) evaluated at the pseudo singular point is negative, i.e. if σ2 < 0 then, system (3.1) exhibits a canard solution which evolves from the attractive part of the slow manifold towards its repelling part. Proof. The method used by Benoı̂t [4], Szmolyan and Wechselberger [31] and Wechselberger [37] requires to implement a “desingularization procedure” which implies to project the “normalized slow dynamics” (3.8) on the tangent bundle of the critical manifold and, to evaluate the Jacobian of the projection of this “desingularized system” (3.9) at the pseudo singular points. As previously recalled, the method presented in this paper does not use the “desingularized system” (3.8) but the “normalized slow dynamics” (3.9). So, to prove the Proposition 3.4, 14 J.-M. GINOUX ET AL. we have just to show that the determinant of the Jacobian of the “desingularized system” (3.8) is identical to the sum σ2 of all second-order diagonal minors of the Jacobian matrix of the “normalized slow dynamics” (3.9). According to equation (3.16): 3 X ii 11 22 33 σ2 = J(F 1 ,G1 ,G2 ) = J(F 1 ,G1 ,G2 ) + J(F1 ,G1 ,G2 ) + J(F 1 ,G1 ,G2 ) . i=1 While using the generic hypotheses equations (3.12), it is easy to prove that: 11 33 J(F 1 ,G1 ,G2 ) = J(F1 ,G1 ,G2 ) =0 22 The remaining determinant J(F 1 ,G1 ,G2 ) is exactly that of the “desingularized system” (3.9). So, Proposition 3.4 can be also used to state the existence of canard solution for such systems. 4. Canards existence in the Hindmarsh–Rose model The Hindmarsh-Rose model [19] describes the basic properties of individual neurons and appears as a reduc- tion of the conductance based in the Hodgkin-Huxley model for neural spiking, see for more details [21]. Thus, the three-dimensional Hindmarsh-Rose polynomial ordinary differential system was originally written as: dx = y − ax3 + bx2 − z + I, dt dy = c − dx2 − y, (4.1) dt dz = r [s (x − α) − z] , dt where x is a transmembrane neuron potential, y and z are the characteristics of ionic currents dynamic, I is ambient current. The other parameters (a, b, c, d, I, s, α and r) reflect the physical features of the neurons and the dot indicates derivative with respect to the time t. We notice that the parameter r  1. Existence of canard solutions in such system (4.1) has been originally suspected by Shilnikov et al. ([29], p. 2149) and highlighted by Shchepakina [28]. Thus, by posing x → y2 , y → y1 , z → x1 and t0 → εt with ε = r and according to the previous definitions, the Hindmarsh-Rose model may be written as a three-dimensional singularly perturbed system with k = 1 slow variable and m = 2 fast variables: ẋ1 = f1 (x1 , y1 , y2 ) = s (y2 − α) − x1 , (4.2a) εẏ1 = g1 (x1 , y1 , y2 ) = c − dy22 − y1 , (4.2b) εẏ2 = g2 (x1 , y1 , y2 ) = y1 − ay23 + by22 − x1 + I, (4.2c) where x1 ∈ R, ~y = (y1 , y2 )t ∈ R2 , 0 < ε  1 and the functions fi and gi are assumed to be C 2 functions of (x1 , y1 , y2 ) and the dot now indicates derivative with respect to the time t0 . 4.1. Critical manifold The critical manifold equation of system (4.2) is defined by setting ε = 0 in equations (4.2b) and (4.2c). Thus, we have: g1 (x1 , y1 , y2 ) = c − dy22 − y1 = 0, (4.3a) CANARDS EXISTENCE IN THE HINDMARSH–ROSE MODEL 15 g2 (x1 , y1 , y2 ) = y1 − ay23 + by22 − x1 + I = 0. (4.3b) This leads to the following critical manifold equation: y1 = φ (y2 ) = c + I − ay23 + (b − d) y22 . (4.4) 4.2. Constrained system According to equations (3.7), we have the following constrained system: ẋ1 = f1 (x1 , y1 , y2 ) = s (y2 − α) − x1 , −2dy2 f1 (x1 , y1 , y2 ) ẏ1 = − , 3ay22 − 2 (b − d) y2 f1 (x1 , y1 , y2 ) (4.5) ẏ2 = − , 3ay22 − 2 (b − d) y2 0 = c − dy22 − y1 , 0 = y1 − ay23 + by22 − x1 + I. 4.3. Normalized slow dynamics By rescaling the time by setting t0 = −(3ay22 −2 (b − d) y2 )τ we obtain the “normalized slow dynamics”: ẋ1 = − 3ay22 − 2 (b − d) y2 [s (y2 − α) − x1 ] = F1 (x1 , y1 , y2 ) ,   ẏ1 = −2dy2 [s (y2 − α) − x1 ] = G1 (x1 , y1 , y2 ) , ẏ2 = [s (y2 − α) − x1 ] = G2 (x1 , y1 , y2 ) , (4.6) 0=c− dy22 − y1 , 0 = y1 − ay23 + by22 − x1 + I. 4.4. Desingularized system on the critical manifold The projection of the normalized slow dynamics (4.5) on the tangent bundle of the critical manifold reads: ẋ1 = − 3ay22 − 2 (b − d) y2 [s (y2 − α) − x1 ] = F1 (x1 , y2 ) ,   (4.7) ẏ2 = [s (y2 − α) − x1 ] = G2 (x1 , y2 ) . 4.5. Pseudo singular points According to equations (3.10) the pseudo-singular points of system (4.2) are: (x̃1 , ỹ1 , ỹ2 ) = (c + I, c, 0) , (4.8a)   3 2 2 2 (x̃1 , ỹ1 , ỹ2 ) = ỹ1 − aỹ2 + bỹ2 + I, c − dỹ2 , (b − d) . (4.8b) 3a 16 J.-M. GINOUX ET AL. 4.6. Canard existence The Jacobian matrix of system (4.6) evaluated at the pseudo singular points (4.8a) reads:   0 0 −2 (b − d) (c + I + sα)   J(F1 ,G1 ,G2 ) =  0 0 2d (c + I + sα) . (4.9)     −1 0 s According to equations (3.16) we find that: p = σ1 = T r [J] = s, q = σ2 = −2 (b − d) (c + I + sα) . Thus, according to Proposition 3.4, the pseudo singular point is of saddle-type if: − 2 (b − d) (c + I + sα) < 0. (4.10) In her work Shchepakina [28] used the following parameter set: a = 1, b = 3, c = 1, d = 0.275255, I = 2.7 and α = −1.2. She found a canard without head (see Fig. 1) for the “duck parameter” value s = 3.0810445478558141214. According to equation (4.10) and with such a parameter set, i.e. b − d > 0, the pseudo singular point is of saddle-type if and only if: c+I s< . α With c = 1, I = 2.7 and α = −1.2, we find that: s < 3.0833. Thus, Shchepakina highlighted a canard without head in the Hindmarsh-Rose model (see Fig. 1) for the “duck parameter” value s = 3.0810445478558141214 < 3.0833. In the inset of Figure 1, the zoom in highlights a large distance between the canard solution and that of the critical manifold (4.4). This is due to the fact that this latter corresponds to zero-order approximation in ε of the slow invariant manifold. Nevertheless, while using the so-called Flow Curvature Method Ginoux and Rossetto [13] have already provided a second-order approximation in ε of the slow invariant manifold of the Hindmarsh-Rose model (4.1). The result is presented in Figure 2. Doing the change of variables (x1 , y1 , y2 ) → (y, x1 , x3 ) system (4.2) writes ẏ = f (x1 , x3 , y) = s(x3 − α) − y, εẋ1 = h1 (x1 , x3 , y) = c − dx23 − x1 , (4.11) εẋ3 = h3 (x1 , x3 , y) = x1 − αx33 + bx23 − y + I. Using the previous notation we get ∂h1 a11 = (x̃1 , x̃3 , ỹ) = −1, ∂x1 ∂h1 a14 = (x̃1 , x̃3 , ỹ) = 0, ∂y ∂h3 a33 = (x̃1 , x̃3 , ỹ) = −3αx̃23 + 2bx̃3 = 0, ∂x3 CANARDS EXISTENCE IN THE HINDMARSH–ROSE MODEL 17 Figure 1. Singular canard solution of the Hindmarsh-Rose (4.1) model in the (x3 , y) plane phase with the following parameter set: a = 1, b = 3, c = 1, d = 0.275255, I = 2.7, α = −1.2 and for the “duck parameter” value s = 3.0810445478558141214. where x̃1 = c, x̃3 = 0 and ỹ = c + I. The solution of the corresponding equations (3.21) for system (4.11) with δ = 0 are u1 = u01 = 0, u2 = u02 = 0, v = −(αs + c + I)t. So the solutions for the corresponding equations (3.21) for system (4.11) with δ ≈ 0 are u1 = u01 = L(δ), u2 = u02 = L(δ), v = −(αs + c + I)t + L(δ). Remark. We notice on the one hand that the Flow Curvature Method provides a better of the approximation of the slow invariant manifold of such system and, on the other hand, that it could also be used to highlight the bifurcation leading to a canard solution as emphasized by Ginoux and Llibre [14] in the case of the Van der Pol model. The Jacobian matrix of system (4.6) evaluated at the pseudo singular points (4.8b) reads:   0 0 −2 (b − d) [s (ỹ2 − α) − x̃1 ]   J(F1 ,G1 ,G2 ) = −2dỹ2 0 2d [s (α − 2ỹ2 ) + x̃1 ]  . (4.12)    −1 0 s 18 J.-M. GINOUX ET AL. Figure 2. Singular canard solution of the Hindmarsh-Rose (4.1) model in the (x3 , y) plane phase, its critical manifold (in green) and the second-order approximation in ε of the slow invariant manifold (in blue) with the following parameter set: a = 1, b = 3, c = 1, d = 0.275255, I = 2.7, α = −1.2 and for the “duck parameter” value s = 3.0810445478558141214. According to equations (3.16) we find that: p = σ1 = T r [J] = s, q = σ2 = −2 (b − d) [s (ỹ2 − α) − x̃1 ] . Thus, according to Proposition 3.4, the pseudo singular point is of saddle-type if and only if: − 2 (b − d) [s (z̃ − α) − x̃] < 0. (4.13) In her work Shchepakina [28] used the following parameter set: a = 1, b = 3, c = 1, d = 0.275255, I = 2.7 and α = −1.2. According to equation (4.13) and with such a parameter set, i.e. b − d > 0, the pseudo singular point is of saddle-type if and only if: 3 4 (b − d) + 27a2 (c + I) s< 9a [2 (b − d) − 3aα] With c = 1, I = 2.7 and α = −1.2, we find that: s < 2.2200954. Thus, we have highlighted a canard with head in the Hindmarsh-Rose model (see Fig. 3) for the “duck parameter” value s = 2.220095 < 2.2200954. For this parameters set the second-order approximation in ε of the slow invariant manifold of the Hindmarsh-Rose model (4.1) can be provided while using the Flow Curvature Method introduced by Ginoux and Rossetto [13]. The result is presented in Figure 3. CANARDS EXISTENCE IN THE HINDMARSH–ROSE MODEL 19 Figure 3. Maximal canard solution of the Hindmarsh-Rose (4.1) model in the (x3 , y) plane phase, its critical manifold (in green) and the second-order approximation in ε of the slow invariant manifold (in blue) with the following parameter set: a = 1, b = 3, c = 1, d = 0.275255, I = 2.7, α = −1.2 and for the “duck parameter” value s = 2.220095. 5. Discussion In a previous paper entitled: “Canards Existence in Memristor’s Circuits” (see Ginoux & Llibre [16]) we have proposed a new method for proving the existence of “canard solutions” for three and four-dimensional singularly perturbed systems with only one fast variable which is based on the stability of folded singularities of the normalized slow dynamics (also called “desingularized vector field” in Llibre et al. [26]) deduced from a well-known property of linear algebra. Thus, we proved that this unique condition is completely identical to that provided by Benoı̂t [4], Szmolyan and Wechselberger [31] and Wechselberger [37]. In a second paper entitled: “Canards Existence in FitzHugh-Nagumo and Hodgkin-Huxley Neuronal Models” (see Ginoux and Llibre [15]) we extended this method to the case of four-dimensional singularly perturbed systems with k = 2 slow and m = 2 fast variables. In this work we have extended this new method to the case of three-dimensional singularly perturbed systems with one slow and two fast variables and we have stated that the condition for the existence of “canard solutions” in such systems is exactly identical to those proposed in our previous paper. This result confirms the genericity of the condition (σ2 < 0) that we have highlighted and provides a simple and efficient tool for testing the occurrence of “canard solutions” in any three or four-dimensional singularly perturbed systems with one or two fast variables. Applications of this method to the famous coupled Hindmarsh-Rose model has enabled to confirm the existence of “canard solutions” in such systems as already stated by Shchepakina [28]. However, in this paper, only the case of pseudo singular points or folded singularities of saddle-type has been analyzed. Of course, the case of of pseudo singular points or folded singularities of node-type and focus-type could be also studied with the same method. Acknowledgements. Authors would like to thank the reviewers for their fruitful comments and suggestions which helped us to improve the presentation of this paper. This work is supported by the Ministerio de Economı́a, Industria y 20 J.-M. GINOUX ET AL. Competitividad, Agencia Estatal de Investigación grant MTM2016-77278-P (FEDER), the Agència de Gestió d’Ajuts Universitaris i de Recerca grant 2017SGR1617, and the H2020 European Research Council grant MSCA-RISE-2017- 777911. The first author would like to thank Pr. J. Llibre for is kind invitation at Universitat Autònoma de Barcelona. References [1] J. Argémi, Approche qualitative d’un problème de perturbations singulières dans R4 , in Equadiff 1978, edited by R. Conti, G. Sestini and G. Villari (1978) 330–340. [2] J.S. Baris, On an integral manifold of a nonregularly perturbed differential system. Usp. Mat. Nauk 20 (1968) 439–448. [3] S.A. Belikov and S.N. Samborskii, ‘Canard’ cycles of fast-slow fields with one-dimensional slow component. Mat. Zametki 49 (1991) 3–14; Math. Notes 49 (1991) 339–346. [4] E. Benoı̂t, Systèmes lents-rapides dans R3 et leurs canards. Société Mathématique de France, Astérisque 190–110 (1983) 159–191. [5] E. Benoı̂t and C. Lobry, Les canards de R3 . CR. Acad. Sci. Paris, Série I 294 (1982) 483–488. [6] F. Dumortier and R. Roussarie, Canard cycles and center manifolds. Mem. Amer. Math. Soc. 577 (1996) 1–100. [7] N. Fenichel, Persistence and smoothness of invariant manifolds for flows. Ind. Univ. Math. J. 21 (1971) 193–225. [8] N. Fenichel, Asymptotic stability with rate conditions. Ind. Univ. Math. J. 23 (1974) 1109–1137. [9] N. Fenichel, Asymptotic stability with rate conditions II. Ind. Univ. Math. J. 26 (1977) 81–93. [10] N. Fenichel, Geometric singular perturbation theory for ordinary differential equations. J. Diff. Eq. 31 (1979) 53–98. [11] R. FitzHugh, Impulses and physiological states in theoretical models of nerve membranes. Biophys. J. 1 (1961) 445–466. [12] K.O. Friedrichs and W.R. Wasow, Singular perturbations of non-linear oscillations. Duke Math. J. 13 (1946) 367–381. [13] J.M. Ginoux and B. Rossetto, Slow Manifold of a Neuronal Bursting Model, in Emergent Properties in Natural and Artificial Dynamical Systems, edited by M.A. Aziz-Alaoui and C. Bertelle, Understanding Complex Systems, Series. Springer-Verlag, Heidelberg (2006), pp. 119–128. [14] J.M. Ginoux and J. Llibre, Flow curvature method applied to canard explosion. J. Phys. A: Math. Theor. 44 (2011) 465203. [15] J.M. Ginoux and J. Llibre, Canards existence in FitzHugh-Nagumo and Hodgkin-Huxley neuronal models. Math. Probl. Eng. 15 (2015) 342010. [16] J.M. Ginoux and J. Llibre, Canards in Memristor’s Circuits. Qual. Theory Dyn. Syst. 15 (2016) 381–431. [17] V. Gol’dshtein, A. Zinoviev, V. Sobolev and E. Shchepakina, Criterion for thermal explosion with reactant consumption in a dusty gas. Proc. London Roy. Soc. Series A 452 (1996) 2103–2119. [18] J. Guckenheimer and R. Haiduc, Canards at folded nodes. Mosc. Math. J. 5 (2005) 91–103. [19] J.L. Hindmarsh and R.M. Rose, A model of neuronal bursting using three coupled first order differential equations. Proc. R. Soc. Lond. Ser. B: Biol. Sci. 221 (1984) 87–102; 246 (2009) 541–551. [20] M.W. Hirsch, C.C. Pugh and M. Shub, Invariant manifolds. Springer-Verlag, New-York (1977). [21] A.L. Hodgkin and A.F. Huxley, A quantitative description of membrane current and its application to conduction and excitation in nerve. J. Physiol. Lond. 117 (1952) 500–544. [22] C.K.R.T. Jones, Geometric Singular Perturbation Theory in Dynamical Systems, Montecatini Terme, edited by L. Arnold. Vol. 1609 of Lecture Notes in Mathematics. Springer-Verlag (1994), pp. 44–118. [23] T. Kaper, An introduction to geometric methods and dynamical systems theory for singular perturbation problems, in Ana- lyzing multiscale phenomena using singular perturbation methods. Amer. Math. Soc., Baltimore, MD/Providence, RI (1998) 85–131. [24] H.W. Knobloch and B. Aulback, Singular perturbations and integral manifolds. J. Math. Phys. Sci. 18 (1984) 415–424. [25] A.Yu. Kolesov and N.Kh. Rozov, Duck cycles of three-dimensional relaxation systems with one fast and two slow variables. Differ. Uravn. 32 (1996) 180–184; Differ. Equ. 32 (1996) 181–185. [26] J. Llibre, J. Sotomayor and M. Zhitomirskii, Impasse bifurcations of constrained systems, Differential Equations and Dynamical Systems, Fields Institute Communications. Amer. Mat. Soc. A 31 (2002) 235–255. [27] J.S. Nagumo, S. Arimoto and S. Oshizawa, An active pulse transmission line simulating nerve axon. Proc. Inst. Radio Eng. 50 (1962) 2061–2070. [28] E.A. Shchepakina, Three scenarios for changing of stability in the dynamic model of nerve conduction. Math. Model. Inf. Technol. Nanotechnol. (ITNT-2016) 1638 (2016) 664–673. [29] A. Shilnikov and M. Kolomiets, Methods of the qualitiative theory for the Hindmarsh-Rose model: a case study. A tutorial. Int. J. Bifurcat. Chaos 18 (2008) 2141–2168. [30] V.A. Sobolev and E.A. Shchepakina, Duck trajectories in a problem of combustion theory. Differ. Uravn. 32 (1996) 1175–1184; Differ. Equ. 32 (1996) 1177–1186. [31] P. Szmolyan and M. Wechselberger, Canards in R3 . J. Differ. Eqs. 177 (2001) 419–453. [32] F. Takens, Constrained equations, a study of implicit differential equations and their discontinuous solutions. Structural stability, the theory of catastrophes and applications in the sciences, in Vol. 525 of Springer Lecture Notes in Math. (1976) 143–234. [33] K. Tchizawa, On relative stability in 4-dimensional duck solutions. J. Math. Syst. Sci. 2 (2012) 558–563. [34] K. Tchizawa, On the two methods for finding 4-dimensional duck solutions. Appl. Math. Sci. Res. Publ. 5 (2014) 16–24. CANARDS EXISTENCE IN THE HINDMARSH–ROSE MODEL 21 [35] A.N. Tikhonov, On the dependence of solutions of differential equations on a small parameter. Mat. Sbornik N.S. 31 (1948) 575–586. [36] M. Wechselberger, Existence and Bifurcation of Canards in R3 in the case of a Folded Node. SIAM J. Appl. Dyn. Syst. 4 (2005) 101–139. [37] M. Wechselberger, À propos de canards. Trans. Amer. Math. Soc. 364 (2012) 3289–3309. [38] S. Wiggins, Normally hyperbolic invariant manifolds in dynamical systems. Springer-Verlag, New York (1994). [39] K.V. Zadiraka, On the Integral Manifold of a System of Differential Equations containing a small parameters. Dokl. Akad. Nauk. SSSR 115 (1957) 646–649. [40] K.V. Zadiraka, On a Nonlocal Integral Manifold of a Nonregularly Perturbed Differential System. Ukrain. Math. Z. 17 (1695) 47–63; translated in AMS Transl. Ser. 2 89 (1970) 29–49.