# Canards Existence in FitzHugh-Nagumo and Hodgkin-Huxley Neuronal Models

Authors: Jean-Marc Ginoux; Jaume Llibre
Journal: Mathematical Problems in Engineering
 (2015)

---

## Abstract

In a previous paper we have proposed a new method for proving the existence of “canard solutions” for three- and four-dimensional singularly perturbed systems with only onefast variable which improves the methods used until now. The aim of this work is to extend this method to the case of four-dimensional singularly perturbed systems with twoslow and twofast variables. This method enables stating a unique generic condition for the existence of “canard solutions” for such four-dimensional singularly perturbed systems which is based on the stability offolded singularities (pseudo singular points in this case) of thenormalized slow dynamics deduced from a well-known property of linear algebra. This unique generic condition is identical to that provided in previous works. Application of this method to the famous coupled FitzHugh-Nagumo equations and to the Hodgkin-Huxley model enables showing the existence of “canard solutions” in such systems.


---

## Body

## 1. Introduction

In the beginning of the eighties, Benoît and Lobry [1], Benoît [2], and then Benoît [3] in his PhD- thesis studied canard solutions in R3. In the article entitled “Systèmes Lents-Rapides dans R3 et Leurs Canards,” Benoît [2, p. 170] proved the existence of canards solution for three-dimensional singularly perturbed systems with twoslow variables and onefast variable while using “nonstandard analysis” according to a theorem which stated that canard solutions exist in such systems provided that thepseudo singular point (this concept has been originally introduced by Argémi [4]; see Section 2.8) of theslow dynamics, that is, of thereduced vector field, is ofsaddle-type. Nearly twenty years later, Szmolyan and Wechselberger [5] extended “Geometric Singular Perturbation Theory (see Fenichel [6, 7], O’Malley [8], Jones [9], and Kaper [10])” to canards problems in R3 and provided a “standard version” of Benoît’s theorem [2]. Very recently, Wechselberger [11] generalized this theorem for n-dimensional singularly perturbed systems with kslow variables and mfast (1). The methods used by Szmolyan and Wechselberger [5] and Wechselberger [11] require implementing a “desingularization procedure” which can be summarized as follows: first, they compute thenormal form of such singularly perturbed systems which is expressed according to some coefficients (a and b for dimension three and a~, b~, and c~j for 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 “normalized slow dynamics” by Benoît [2, p. 166]) of such anormal form on the tangent bundle of the critical manifold. Finally, they evaluate the Jacobian of the projection of this “desingularized vector field” at thefolded singularity (originally calledpseudo singular points by Argémi [4, p. 336]). This leads Szmolyan and Wechselberger [5, p. 427] and Wechselberger [11, p. 3298] to a “classification offolded singularities (pseudo singular points).” Thus, they showed that for three-dimensional singularly perturbed systems suchfolded singularity is ofsaddl-type if the following condition is satisfied, a<0 while for four-dimensional singularly perturbed systems suchfolded singularity is ofsaddle-type if a~<0. Then, Szmolyan and Wechselberger [5, p. 439] and Wechselberger [11, p. 3304] established their Theorem 4.1. which states that “In the folded saddle and in the folded node case singular canards perturb to maximal canard for sufficiently small ε.” However, in their works neither Szmolyan and Wechselberger [5] nor Wechselberger [11] did not provide (to our knowledge) the expression of these constants (a and a~) which are necessary to state the existence of canard solutions in such systems.In a previous paper entitled “Canards Existence in Memristor’s Circuits” (see Ginoux and Llibre [12]) we first provided the expression of these constants and then showed that they can be directly determined starting from thenormalized slow dynamics and not from the projection of the “desingularized vector field” of thenormal form. This method enabled stating a unique “generic” condition for the existence of “canard solutions” for such three- and four-dimensional singularly perturbed systems which is based on the stability offolded singularities of thenormalized slow dynamics deduced from a well-known property of linear algebra. This unique condition which is completely identical to that provided by Benoît [2] and then by Szmolyan and Wechselberger [5] and finally by Wechselberger [11] is “generic” since it is exactly the same for singularly perturbed systems of dimensions three and four with only onefast variable.The aim of this work is to extend this method to the case of four-dimensional singularly perturbed systems withk=2slow and m=2fast variables. Since the dimension of the system is m=k+m, such problem is known as “canards existence in R2+2.” Moreover, in this particular case, where k=m=2, thefolded singularities of Wechselberger [11, p. 3298] are nothing else but thepseudo singular points of the late Argémi [4] as we will see below. Following the previous works, we show that for such four-dimensional singularly perturbed systemspseudo singular points are ofsaddle-type if a~<0. Then, according to Theorem 4.1. of Wechselberger [11, p. 3304] we provide the expression of this constant a~ which is necessary to establish the existence of canard solutions in such systems. So, we can state that the condition a~<0 for existence of canards in such R2+2 is “generic” since it is exactly the same for singularly perturbed systems of dimensions three and four with only onefast variable.The outline of this paper is as follows. In Section1, definitions of singularly perturbed system, critical manifold, reduced system, “constrained system,” canard cycles, folded singularities, and pseudo singular points are recalled. The method proposed in this paper is presented in Section 2 for the case of four-dimensional singularly perturbed systems with twofast variables. In Section 3, applications of this method to the famous coupled FitzHugh-Nagumo equations and to the Hodgkin-Huxley model enable showing the existence of “canard solutions” in such systems.


## 2. Definitions

### 2.1. Singularly Perturbed Systems

According to Tikhonov [13], Jones [9], and Kaper [10]singularly perturbed systems are defined as(1)x→′=εf→x→,y→,ε,y→′=g→x→,y→,ε,where x→∈Rk, y→∈Rm, ε∈R+, and the prime denotes differentiation with respect to the independent variable t′. The functions f→ and g→ are assumed to be C∞ functions (in certain applications these functions will be supposed to be Cr, r⩾1) 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 when0<ε≪1, that is, ε is a small positive number, the variable x→ is calledslow variable and y→ is calledfast variable. Using Landau’s notation, 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 thatx→ evolves at an Oε rate while y→ evolves at an O1slow rate. Reformulating system (1) in terms of the rescaled variable t=εt′, we obtain(2)x→˙=f→x→,y→,ε,εy→˙=g→x→,y→,ε.The dot represents the derivative with respect to the new independent variablet.The independent variablest′ and t are referred to thefast andslow times, respectively, and (1) and (2) are called thefast andslow systems, respectively. These systems are equivalent whenever ε≠0, and they are labeledsingular perturbation problems when 0<ε≪1. The label “singular” stems in part from the discontinuous limiting behavior in system (1) as ε→0.


### 2.2. Reduced Slow System

In such case system (2) leads to a differential-algebraic system (DAE) calledreduced slow system whose dimension decreases from k+m=n to m. Then, theslow variable x→∈Rk partially evolves in the submanifold M0 called thecritical manifold (it represents the approximation of the slow invariant manifold, with an error of O(ε)). Thereduced slow system is(3)x→˙=f→x→,y→,ε,0→=g→x→,y→,ε.


### 2.3. Slow Invariant Manifold

Thecritical manifold is defined by(4)M0≔x→,y→:g→x→,y→,0=0→.Such a normally hyperbolic invariant manifold (4) of thereduced slow system (3) persists as a locally invariantslow manifold of the full problem (1) for ε sufficiently small. This locallyslow invariant manifold is O(ε) close to thecritical manifold.WhenDx→f→ is invertible, thanks to the Implicit Function Theorem, M0 is given by the graph of a C∞ function x→=G→0y→ for y→∈D, where D⊆Rk is a compact, simply connected domain and the boundary of D is a (k-1)-dimensional C∞ submanifold (the set D is overflowing invariant with respect to (2) when ε=0; see Kaper [10] and Jones [9]).According to Fenichel [6, 7] theory if 0<ε≪1 is sufficiently small, then there exists a function G→y→,ε defined on D such that the manifold,(5)Mε≔x→,y→:x→=G→y→,ε,is locally invariant under the flow of system (1). Moreover, there exist perturbed local stable (or attracting) Ma and unstable (or repelling) Mr branches of theslow invariant manifoldMε. Thus, normal hyperbolicity of Mε is lost via a saddle-node bifurcation of thereduced slow system (3). Then, it gives rise to solutions of “canard” type.


### 2.4. Canards: Singular Canards and Maximal Canards

Acanard is a solution of singularly perturbed dynamical system (1) following theattracting branch Ma of theslow invariant manifold, passing near a bifurcation point located on the fold of thisslow invariant manifold, and then following therepelling branch Mr of theslow invariant manifold.Asingular canard is a solution ofreduced slow system (3) following theattracting branch Ma,0 of thecritical manifold, passing near a bifurcation point located on the fold of thiscritical manifold, and then following therepelling branch Mr,0 of thecritical manifold.Amaximal canard corresponds to the intersection of the attracting and repelling branches Ma,ε∩Mr,ε of the slow manifold in the vicinity of a nonhyperbolic point.According to Wechselberger [11, 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 (1) that follow an attracting branch Ma,ε of the slow manifold and then follow, rather surprisingly, a repelling/saddle branch Mr,ε 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 Ma,ε and Mr,ε. However, in the singular limit ε→0, such a family of canards is represented by a unique singular canard.”Canards are special class of solutions of singularly perturbed dynamical systems for which normal hyperbolicity 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 changes (see Benoît [2, 14], Szmolyan and Wechselberger [5], and Wechselberger [11, 15]).


### 2.5. Constrained System

In order to characterize the “slow dynamics,” that is, the slow trajectory of thereduced slow system (3) (obtained by setting ε=0 in (2)), Takens [16] introduced the “constrained system” defined as follows:(6)x→˙=f→x→,y→,0,Dy→g→·y→˙=-Dx→g→·f→x→,y→,0,0→=g→x→,y→,0.Since, according to Fenichel [6, 7], thecritical manifoldg→x→,y→,0 may be considered as locally invariant under the flow of system (1), we have(7)dg→dtx→,y→,0=0⟺Dx→g→·x→˙+Dy→g→·y→˙=0→.By replacingx→˙ by f→x→,y→,0 leads to(8)Dx→g→·f→x→,y→,0+Dy→g→·y→˙=0→.This justifies the introduction of theconstrained system.Now, letadj(Dy→g→) denote the adjoint of the matrix Dy→g→ which is the transpose of the cofactor matrix Dy→g→; then while multiplying the left-hand side of (6) by the inverse matrix (Dy→g→)-1 obtained by the adjoint method we have(9)x→˙=f→x→,y→,0,det⁡Dy→g→y→˙=-adjDy→g→·Dx→g→·f→x→,y→,0,0→=g→x→,y→,0.


### 2.6. Normalized Slow Dynamics

Then, by rescaling the time by settingt=-det⁡(Dy→g→)τ we obtain the following system which has been called by Benoît [2, p. 166] “normalized slow dynamics”:(10)x→˙=-det⁡Dy→g→f→x→,y→,0,y→˙=adjDy→g→·Dx→g→·f→x→,y→,0,0→=g→x→,y→,0,where the overdot now denotes the time derivation with respect to τ.Let us notice that Argémi [4] proposed to rescale time by setting t=-det⁡Dy→g→sgn⁡(det⁡(Dy→g→))τ in order to keep the same flow direction in (10) as in (9).


### 2.7. Desingularized Vector Field

By application of the Implicit Function Theorem, let suppose that we can explicitly express from (4), 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 thecritical manifoldM0 at thepseudo singular point. This leads to the so-calleddesingularized vector field:(11)χ→˙=-det⁡Dy→g→f→χ→,y→,0,y→˙=adjDy→g→·Dx→g→·f→χ→,y→,0.


### 2.8. Pseudo Singular Points and Folded Singularities

As recalled by Guckenheimer and Haiduc [17, p. 91],pseudo singular points have been introduced by the late Argémi [4] for low-dimensional singularly perturbed systems and are defined as singular points of the “normalized slow dynamics” (10). Twenty-three years later, Szmolyan and Wechselberger [5, p. 428] called suchpseudo singular points folded singularities. In a recent publication entitled “A Propos de Canards” Wechselberger [11, p. 3295] proposed to define such singularities for n-dimensional singularly perturbed systems with kslow variables and mfast as the solutions of the following system:(12)det⁡Dy→g→=0,adjDy→g→·Dx→g→·f→x→,y→,0=0→,g→x→,y→,0=0→.Thus, for dimensions higher than three, his concept encompasses that of Argémi. Moreover, Wechselberger [11, p. 3296] proved thatfolded singularities form a (k-2)-dimensional manifold. Thus, for k=2 thefolded singularities are nothing else than thepseudo singular points defined by Argémi [4], while for k⩾3 thefolded singularities are no more points but a (k-2)-dimensional manifold. Moreover, let us notice on one hand that the original system (1) includes n=k+m variables and on the other hand that system (12) comprises p=2m+1 equations. However, in the particular case k=m=2, two equations of system (12) are linearly dependent. So, such system only comprises p=2m=2k equations. So, all the variables (unknowns) of system (12) can be determined. The solutions of this system are calledpseudo singular points. We will see in Section 3 that the stability analysis of thesepseudo singular points will give rise to a condition for the existence of canard solutions in the original system (1).


## 2.1. Singularly Perturbed Systems

According to Tikhonov [13], Jones [9], and Kaper [10]singularly perturbed systems are defined as(1)x→′=εf→x→,y→,ε,y→′=g→x→,y→,ε,where x→∈Rk, y→∈Rm, ε∈R+, and the prime denotes differentiation with respect to the independent variable t′. The functions f→ and g→ are assumed to be C∞ functions (in certain applications these functions will be supposed to be Cr, r⩾1) 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 when0<ε≪1, that is, ε is a small positive number, the variable x→ is calledslow variable and y→ is calledfast variable. Using Landau’s notation, 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 thatx→ evolves at an Oε rate while y→ evolves at an O1slow rate. Reformulating system (1) in terms of the rescaled variable t=εt′, we obtain(2)x→˙=f→x→,y→,ε,εy→˙=g→x→,y→,ε.The dot represents the derivative with respect to the new independent variablet.The independent variablest′ and t are referred to thefast andslow times, respectively, and (1) and (2) are called thefast andslow systems, respectively. These systems are equivalent whenever ε≠0, and they are labeledsingular perturbation problems when 0<ε≪1. The label “singular” stems in part from the discontinuous limiting behavior in system (1) as ε→0.


## 2.2. Reduced Slow System

In such case system (2) leads to a differential-algebraic system (DAE) calledreduced slow system whose dimension decreases from k+m=n to m. Then, theslow variable x→∈Rk partially evolves in the submanifold M0 called thecritical manifold (it represents the approximation of the slow invariant manifold, with an error of O(ε)). Thereduced slow system is(3)x→˙=f→x→,y→,ε,0→=g→x→,y→,ε.


## 2.3. Slow Invariant Manifold

Thecritical manifold is defined by(4)M0≔x→,y→:g→x→,y→,0=0→.Such a normally hyperbolic invariant manifold (4) of thereduced slow system (3) persists as a locally invariantslow manifold of the full problem (1) for ε sufficiently small. This locallyslow invariant manifold is O(ε) close to thecritical manifold.WhenDx→f→ is invertible, thanks to the Implicit Function Theorem, M0 is given by the graph of a C∞ function x→=G→0y→ for y→∈D, where D⊆Rk is a compact, simply connected domain and the boundary of D is a (k-1)-dimensional C∞ submanifold (the set D is overflowing invariant with respect to (2) when ε=0; see Kaper [10] and Jones [9]).According to Fenichel [6, 7] theory if 0<ε≪1 is sufficiently small, then there exists a function G→y→,ε defined on D such that the manifold,(5)Mε≔x→,y→:x→=G→y→,ε,is locally invariant under the flow of system (1). Moreover, there exist perturbed local stable (or attracting) Ma and unstable (or repelling) Mr branches of theslow invariant manifoldMε. Thus, normal hyperbolicity of Mε is lost via a saddle-node bifurcation of thereduced slow system (3). Then, it gives rise to solutions of “canard” type.


## 2.4. Canards: Singular Canards and Maximal Canards

Acanard is a solution of singularly perturbed dynamical system (1) following theattracting branch Ma of theslow invariant manifold, passing near a bifurcation point located on the fold of thisslow invariant manifold, and then following therepelling branch Mr of theslow invariant manifold.Asingular canard is a solution ofreduced slow system (3) following theattracting branch Ma,0 of thecritical manifold, passing near a bifurcation point located on the fold of thiscritical manifold, and then following therepelling branch Mr,0 of thecritical manifold.Amaximal canard corresponds to the intersection of the attracting and repelling branches Ma,ε∩Mr,ε of the slow manifold in the vicinity of a nonhyperbolic point.According to Wechselberger [11, 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 (1) that follow an attracting branch Ma,ε of the slow manifold and then follow, rather surprisingly, a repelling/saddle branch Mr,ε 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 Ma,ε and Mr,ε. However, in the singular limit ε→0, such a family of canards is represented by a unique singular canard.”Canards are special class of solutions of singularly perturbed dynamical systems for which normal hyperbolicity 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 changes (see Benoît [2, 14], Szmolyan and Wechselberger [5], and Wechselberger [11, 15]).


## 2.5. Constrained System

In order to characterize the “slow dynamics,” that is, the slow trajectory of thereduced slow system (3) (obtained by setting ε=0 in (2)), Takens [16] introduced the “constrained system” defined as follows:(6)x→˙=f→x→,y→,0,Dy→g→·y→˙=-Dx→g→·f→x→,y→,0,0→=g→x→,y→,0.Since, according to Fenichel [6, 7], thecritical manifoldg→x→,y→,0 may be considered as locally invariant under the flow of system (1), we have(7)dg→dtx→,y→,0=0⟺Dx→g→·x→˙+Dy→g→·y→˙=0→.By replacingx→˙ by f→x→,y→,0 leads to(8)Dx→g→·f→x→,y→,0+Dy→g→·y→˙=0→.This justifies the introduction of theconstrained system.Now, letadj(Dy→g→) denote the adjoint of the matrix Dy→g→ which is the transpose of the cofactor matrix Dy→g→; then while multiplying the left-hand side of (6) by the inverse matrix (Dy→g→)-1 obtained by the adjoint method we have(9)x→˙=f→x→,y→,0,det⁡Dy→g→y→˙=-adjDy→g→·Dx→g→·f→x→,y→,0,0→=g→x→,y→,0.


## 2.6. Normalized Slow Dynamics

Then, by rescaling the time by settingt=-det⁡(Dy→g→)τ we obtain the following system which has been called by Benoît [2, p. 166] “normalized slow dynamics”:(10)x→˙=-det⁡Dy→g→f→x→,y→,0,y→˙=adjDy→g→·Dx→g→·f→x→,y→,0,0→=g→x→,y→,0,where the overdot now denotes the time derivation with respect to τ.Let us notice that Argémi [4] proposed to rescale time by setting t=-det⁡Dy→g→sgn⁡(det⁡(Dy→g→))τ in order to keep the same flow direction in (10) as in (9).


## 2.7. Desingularized Vector Field

By application of the Implicit Function Theorem, let suppose that we can explicitly express from (4), 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 thecritical manifoldM0 at thepseudo singular point. This leads to the so-calleddesingularized vector field:(11)χ→˙=-det⁡Dy→g→f→χ→,y→,0,y→˙=adjDy→g→·Dx→g→·f→χ→,y→,0.


## 2.8. Pseudo Singular Points and Folded Singularities

As recalled by Guckenheimer and Haiduc [17, p. 91],pseudo singular points have been introduced by the late Argémi [4] for low-dimensional singularly perturbed systems and are defined as singular points of the “normalized slow dynamics” (10). Twenty-three years later, Szmolyan and Wechselberger [5, p. 428] called suchpseudo singular points folded singularities. In a recent publication entitled “A Propos de Canards” Wechselberger [11, p. 3295] proposed to define such singularities for n-dimensional singularly perturbed systems with kslow variables and mfast as the solutions of the following system:(12)det⁡Dy→g→=0,adjDy→g→·Dx→g→·f→x→,y→,0=0→,g→x→,y→,0=0→.Thus, for dimensions higher than three, his concept encompasses that of Argémi. Moreover, Wechselberger [11, p. 3296] proved thatfolded singularities form a (k-2)-dimensional manifold. Thus, for k=2 thefolded singularities are nothing else than thepseudo singular points defined by Argémi [4], while for k⩾3 thefolded singularities are no more points but a (k-2)-dimensional manifold. Moreover, let us notice on one hand that the original system (1) includes n=k+m variables and on the other hand that system (12) comprises p=2m+1 equations. However, in the particular case k=m=2, two equations of system (12) are linearly dependent. So, such system only comprises p=2m=2k equations. So, all the variables (unknowns) of system (12) can be determined. The solutions of this system are calledpseudo singular points. We will see in Section 3 that the stability analysis of thesepseudo singular points will give rise to a condition for the existence of canard solutions in the original system (1).


## 3. Four-Dimensional Singularly Perturbed Systems with Two Fast Variables

Four-dimensionalsingularly perturbed dynamical system (2) with k=2slow variables and m=2fast may be written as(13a)x˙1=f1x1,x2,y1,y2,(13b)x˙2=f2x1,x2,y1,y2,(13c)εy˙1=g1x1,x2,y1,y2,(13d)εy˙2=g2x1,x2,y1,y2,where x→=(x1,x2)t∈R2, y→=(y1,y2)∈R2, 0<ε≪1, and the functions fi and gi are assumed to be C2 functions of (x1,x2,y1,y2).


### 3.1. Critical Manifold

The critical manifold equation of system (13a), (13b), (13c), and (13d) is defined by setting ε=0 in ((13c) and (13d)). Thus, we obtain(14a)g1x1,x2,y1,y2=0,(14b)g2x1,x2,y1,y2=0.By application of the Implicit Function Theorem, let us suppose that we can explicitly express from ((14a) and (14b)), say, without loss of generality, x1 and y1 as functions of the other variables:(15a)x1=ϕ1x2,y1,y2,(15b)y1=ϕ2x1,x2,y2.


### 3.2. Constrained System

Theconstrained system is obtained by equating with zero the time derivative of g1,2x1,x2,y1,y2:(16a)dg1dt=∂g1∂x1x˙1+∂g1∂x2x˙2+∂g1∂y1y˙1+∂g1∂y1y˙2=0,(16b)dg2dt=∂g2∂x1x˙1+∂g2∂x2x˙2+∂g2∂y1y˙1+∂g2∂y1y˙2=0.Equations (16a) and (16b) may be written as(17a)∂g1∂y1y˙1+∂g1∂y1y˙2=-∂g1∂x1x˙1+∂g1∂x2x˙2,(17b)∂g2∂y1y˙1+∂g2∂y1y˙2=-∂g2∂x1x˙1+∂g2∂x2x˙2.By solving the system of (17a) and (17b) with two unknowns (y˙1,y˙2) we find(18a)y˙1=-∂g1/∂x1x˙1+∂g1/∂x2x˙2∂g2/∂y2+∂g2/∂x1x˙1+∂g2/∂x2x˙2∂g1/∂y2det⁡Jy1,y2,(18b)y˙2=-∂g1/∂x1x˙1+∂g1/∂x2x˙2∂g2/∂y1+∂g2/∂x1x˙1+∂g2/∂x2x˙2∂g1/∂y1det⁡Jy1,y2.So, we have the following constrained system:(19)x˙1=f1x1,x2,y1,y2,x˙2=f2x1,x2,y1,y2,y˙1=-∂g1/∂x1x˙1+∂g1/∂x2x˙2∂g2/∂y2+∂g2/∂x1x˙1+∂g2/∂x2x˙2∂g1/∂y2det⁡Jy1,y2,y˙2=-∂g1/∂x1x˙1+∂g1/∂x2x˙2∂g2/∂y1+∂g2/∂x1x˙1+∂g2/∂x2x˙2∂g1/∂y1det⁡Jy1,y2,0=g1x1,x2,y1,y2,0=g2x1,x2,y1,y2.


### 3.3. Normalized Slow Dynamics

By rescaling the time by settingt=-det⁡J(y1,y2)τ we obtain the “normalized slow dynamics”:(20)x˙1=-f1x1,x2,y1,y2det⁡Jy1,y2=F1x1,x2,x3,y1,x˙2=-f2x1,x2,y1,y2det⁡Jy1,y2=F2x1,x2,x3,y1,y˙1=∂g1∂x1x˙1+∂g1∂x2x˙2∂g2∂y2-∂g2∂x1x˙1+∂g2∂x2x˙2∂g1∂y2=G1x1,x2,x3,y1,y˙2=∂g1∂x1x˙1+∂g1∂x2x˙2∂g2∂y1-∂g2∂x1x˙1+∂g2∂x2x˙2∂g1∂y1=G2x1,x2,x3,y1,0=g1x1,x2,y1,y2,0=g2x1,x2,y1,y2,where the overdot now denotes the time derivation with respect to τ.


### 3.4. Desingularized System on the Critical Manifold

Then, since we have supposed thatx1 and y1 may be explicitly expressed as functions of the other variables ((15a) and (15b)), they can be used to project the normalized slow dynamics (20) on the tangent bundle of the critical manifold. So, we have(21)x˙2=-f2x1,x2,y1,y2det⁡Jy1,y2=F2x2,y2,y˙2=∂g1∂x1x˙1+∂g1∂x2x˙2∂g2∂y1-∂g2∂x1x˙1+∂g2∂x2x˙2∂g1∂y1=G2x2,y2.


### 3.5. Pseudo Singular Points

Pseudo singular points are defined as singular points of the “normalized slow dynamics,” that is, as the set of points for which we have(22a)det⁡Jy1,y2=0,(22b)∂g1∂x1x˙1+∂g1∂x2x˙2∂g2∂y2-∂g2∂x1x˙1+∂g2∂x2x˙2∂g1∂y2=0,(22c)∂g1∂x1x˙1+∂g1∂x2x˙2∂g2∂y1-∂g2∂x1x˙1+∂g2∂x2x˙2∂g1∂y1=0,(22d)g1x1,x2,y1,y2=0,(22e)g2x1,x2,y1,y2=0.Remark 1.
Let us notice on one hand that (22b) and (22c) are linearly dependent and on the other hand that contrary to previous works we do not use the “desingularized vector field” (21) but the “normalized slow dynamics” (20).The Jacobian matrix of system (20) reads(23)JF1,F2,G1,G2=∂F1∂x1∂F1∂x2∂F1∂y1∂F1∂y2∂F2∂x1∂F2∂x2∂F2∂y1∂F2∂y2∂G1∂x1∂G1∂x2∂G1∂y1∂G1∂y2∂G2∂x1∂G2∂x2∂G2∂y1∂G2∂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 the three-dimensional case to the four-dimensional case. So, first, let us suppose that by a “standard translation” thepseudo singular point can be shifted at the origin O(0,0,0,0) and that by a “standard rotation” of y1-axis that theslow manifold is tangent to (x2,x3,y1)-hyperplane, so we have(24)f10,0,0,0=g10,0,0,0=0,∂g1∂x20,0,0,0=∂g1∂x30,0,0,0=∂g1∂y10,0,0,0=0.Then, let us make the following assumptions for the nondegeneracy of thefolded singularity:(25)f20,0,0,0≠0;∂g1∂x10,0,0,0≠0;∂2g1∂y120,0,0,0≠0.According to these generic hypotheses ((24)-(25)), the Jacobian matrix (23) reads(26)JF1,F2,G1,G2=0000-f2∂P∂x1-f2∂P∂x2-f2∂P∂y1-f2∂P∂y2a31a32a33a34a41a42a43a44,where(27)P=det⁡Jy1,y2,a3i=-f2∂g2∂x2∂2g1∂y2∂xi+∂g2∂y2f2∂2g1∂x2∂xi+∂g1∂x1∂f1∂xifori=1,2,a3i=-f2∂g2∂x2∂2g1∂y2∂yi+∂g2∂y2f2∂2g1∂x2∂yi+∂g1∂x1∂f1∂yifori=3,4,a4i=f2∂g2∂x2∂2g1∂y1∂xi-∂g2∂y1f2∂2g1∂x2∂xi+∂g1∂x1∂f1∂xifori=1,2,a4i=f2∂g2∂x2∂2g1∂y1∂yi-∂g2∂y1f2∂2g1∂x2∂yi+∂g1∂x1∂f1∂yifori=3,4.Thus, we have the following Cayley-Hamilton eigenpolynomial associated with such Jacobian matrix (26) evaluated at thepseudo singular point, that is, at the origin:(28)λ4-σ1λ3+σ2λ2-σ3λ+σ4=0,where σ1=Tr(J) is the sum of all first-order diagonal minors of J, that is, the trace of the Jacobian matrix J, σ2 represents the sum of all second-order diagonal minors of J, and σ3 represents the sum of all third-order diagonal minors of J. It appears that σ4=|J|=0 since one row of the Jacobian matrix (26) is null. So, the eigenpolynomial reduces to(29)λλ3-σ1λ2+σ2λ-σ3=0.But, according to Wechselberger [11], σ3 vanishes at apseudo singular point as it is easy to prove it. So, the eigenpolynomial (29) is reduced to(30)λ2λ2-σ1λ+σ2=0.Letλi be the eigenvalues of eigenpolynomial (30) and let us denote by λ3,4=0 the obvious double root of this polynomial. We have(31)σ1=TrJF1,F2,G1,G2=λ1+λ2=∂g2∂x1∂g1∂y1∂f1∂y2,σ2=∑i=13JF1,F2,G1,G2ii=λ1λ2=∂g1∂y12f22∂2g2∂x22∂2g2∂y22-∂2g2∂x2∂y22+f2∂g2∂x1∂2g2∂y22∂f1∂x2-∂2g2∂x2∂y2∂f1∂y2,where σ1=Tr(J(F1,F2,G1,G2))=p is is the sum of all first-order diagonal minors of J(F1,F2,G1,G2), that is, the trace of the Jacobian matrix J(F1,F2,G1,G2), and σ2=∑i=13J(F1,F2,G1,G2)ii=q represents the sum of all second-order diagonal minors of J(F1,F2,G1,G2). Thus, thepseudo singular point is of saddle-typeiff the following conditions C1 and C2 are verified:(32)C1:Δ=p2-4q>0,C2:q<0.ConditionC1 is systematically satisfied provided that condition C2 is verified. Thus, thepseudo singular point is of saddle-typeiffq<0.


### 3.7. Canard Existence inR2+2

Following the works of Wechselberger [11] it can be stated, while using a standard polynomial change of variables, that any n-dimensional singularly perturbed systems with kslow variables (k⩾2) and mfast (m⩾1) (1) can be transformed into the following “normal form”:(33)x1˙=a~x2+b~y2+Ox1,ϵ,x22,x2y2,y22,x2˙=1+Ox1,x2,y2,ϵ,ϵy1˙=c~y1+Oϵx1,ϵx2,ϵy2,x12,x22,y22,x2y2,ϵy2˙=-x1+y22+Oϵx1,ϵx2,ϵy2,ϵ2,x12y2,y23,x1x2y2.We establish in the Appendix for any four-dimensional singularly perturbed systems (13a), (13b), (13c), and (13d) with k=2slow and m=2fast variables that(34)a~=12f22∂2g2∂x22∂2g2∂y22-∂2g2∂x2∂y22+f2∂g2∂x1∂2g2∂y22∂f1∂x2-∂2g2∂x2∂y2∂f1∂y2,b~=-∂g2∂x1∂f1∂y2,c~=∂g1∂y1.Thus, in his paper Wechselberger [11, p. 3304] provided in the framework of “standard analysis” a generalization of Benoît’s theorem [2] for any n-dimensional singularly perturbed systems with kslow variables (k⩾2) and mfast (m⩾1). According to his theorem presented below as Theorem 2 he proved the existence of canard solutions for the original system (1).Theorem 2.
In the folded saddle case of system (33) singular canards perturb to maximal canards solutions for sufficiently small ε≪1.Proof.
See Wechselberger [11].Since our method does not use the “desingularized vector field” (21) but the “normalized slow dynamics” (20), we have the following proposition.Proposition 3.
If the normalized slow dynamics (20) has apseudo singular point of saddle-type, that is, if the sum σ2 of all second-order diagonal minors of the Jacobian matrix of the normalized slow dynamics (20) evaluated at thepseudo singular point is negative, that is, if σ2<0 then, according to Theorem 2, system (13a), (13b), (13c), and (13d) exhibits a canard solution which evolves from the attractive part of the slow manifold towards its repelling part.Proof.
By making some smooth changes of time and smooth changes of coordinates (see the Appendix) we brought system (13a), (13b), (13c), and (13d) to the following “normal form”:(35)x1˙=a~x2+b~y2+Ox1,ϵ,x22,x2y2,y22,x2˙=1+Ox1,x2,y2,ϵ,ϵy1˙=c~y1+Oϵx1,ϵx2,ϵy2,x12,x22,y22,x2y2,ϵy2˙=-x1+y22+Oϵx1,ϵx2,ϵy2,ϵ2,x12y2,y23,x1x2y2.
Then, we deduce that the condition for thepseudo singular point to be of saddle-type is a~<0. According to (32) it is easy to verify that(36)σ1=TrJF1,F2,G1,G2=λ1+λ2=-b~c~,σ2=∑i=13JF1,F2,G1,G2ii=λ1λ2=2a~c~2.
So, the condition for which thepseudo singular point is of saddle-type, that is, σ2<0, is identical to that proposed by Wechselberger [11, p. 3298] in his theorem; that is, a~<0.So, Proposition3 can be used to state the existence of canard solution for such systems. Application of Proposition 3 to the coupled FitzHugh-Nagumo equations, presented in Section 4, which is a four-dimensional singularly perturbed system with twoslow and twofast variables will enable proving, as many previous works such as those of Tchizawa and Campbell [18] and Tchizawa [19–24], the existence of “canard solutions” in such system. According to Tchizawa [25], it is very important to notice on one hand that the fast equation has 2 dimensions in the system R2+2 and on the other hand that the fast system can give attractive, repulsive, or attractive-repulsive part at eachpseudo singular point. Then, Tchizawa [25] has established that the jumping direction can be shown using the eigenvectors. In the same way we will find again the results of Rubin and Wechselberger [26] concerning the existence of “canard solutions” in the Hodgkin-Huxley model but with a set of more realistic parameters used in Chua et al. [27, 28].


## 3.1. Critical Manifold

The critical manifold equation of system (13a), (13b), (13c), and (13d) is defined by setting ε=0 in ((13c) and (13d)). Thus, we obtain(14a)g1x1,x2,y1,y2=0,(14b)g2x1,x2,y1,y2=0.By application of the Implicit Function Theorem, let us suppose that we can explicitly express from ((14a) and (14b)), say, without loss of generality, x1 and y1 as functions of the other variables:(15a)x1=ϕ1x2,y1,y2,(15b)y1=ϕ2x1,x2,y2.


## 3.2. Constrained System

Theconstrained system is obtained by equating with zero the time derivative of g1,2x1,x2,y1,y2:(16a)dg1dt=∂g1∂x1x˙1+∂g1∂x2x˙2+∂g1∂y1y˙1+∂g1∂y1y˙2=0,(16b)dg2dt=∂g2∂x1x˙1+∂g2∂x2x˙2+∂g2∂y1y˙1+∂g2∂y1y˙2=0.Equations (16a) and (16b) may be written as(17a)∂g1∂y1y˙1+∂g1∂y1y˙2=-∂g1∂x1x˙1+∂g1∂x2x˙2,(17b)∂g2∂y1y˙1+∂g2∂y1y˙2=-∂g2∂x1x˙1+∂g2∂x2x˙2.By solving the system of (17a) and (17b) with two unknowns (y˙1,y˙2) we find(18a)y˙1=-∂g1/∂x1x˙1+∂g1/∂x2x˙2∂g2/∂y2+∂g2/∂x1x˙1+∂g2/∂x2x˙2∂g1/∂y2det⁡Jy1,y2,(18b)y˙2=-∂g1/∂x1x˙1+∂g1/∂x2x˙2∂g2/∂y1+∂g2/∂x1x˙1+∂g2/∂x2x˙2∂g1/∂y1det⁡Jy1,y2.So, we have the following constrained system:(19)x˙1=f1x1,x2,y1,y2,x˙2=f2x1,x2,y1,y2,y˙1=-∂g1/∂x1x˙1+∂g1/∂x2x˙2∂g2/∂y2+∂g2/∂x1x˙1+∂g2/∂x2x˙2∂g1/∂y2det⁡Jy1,y2,y˙2=-∂g1/∂x1x˙1+∂g1/∂x2x˙2∂g2/∂y1+∂g2/∂x1x˙1+∂g2/∂x2x˙2∂g1/∂y1det⁡Jy1,y2,0=g1x1,x2,y1,y2,0=g2x1,x2,y1,y2.


## 3.3. Normalized Slow Dynamics

By rescaling the time by settingt=-det⁡J(y1,y2)τ we obtain the “normalized slow dynamics”:(20)x˙1=-f1x1,x2,y1,y2det⁡Jy1,y2=F1x1,x2,x3,y1,x˙2=-f2x1,x2,y1,y2det⁡Jy1,y2=F2x1,x2,x3,y1,y˙1=∂g1∂x1x˙1+∂g1∂x2x˙2∂g2∂y2-∂g2∂x1x˙1+∂g2∂x2x˙2∂g1∂y2=G1x1,x2,x3,y1,y˙2=∂g1∂x1x˙1+∂g1∂x2x˙2∂g2∂y1-∂g2∂x1x˙1+∂g2∂x2x˙2∂g1∂y1=G2x1,x2,x3,y1,0=g1x1,x2,y1,y2,0=g2x1,x2,y1,y2,where the overdot now denotes the time derivation with respect to τ.


## 3.4. Desingularized System on the Critical Manifold

Then, since we have supposed thatx1 and y1 may be explicitly expressed as functions of the other variables ((15a) and (15b)), they can be used to project the normalized slow dynamics (20) on the tangent bundle of the critical manifold. So, we have(21)x˙2=-f2x1,x2,y1,y2det⁡Jy1,y2=F2x2,y2,y˙2=∂g1∂x1x˙1+∂g1∂x2x˙2∂g2∂y1-∂g2∂x1x˙1+∂g2∂x2x˙2∂g1∂y1=G2x2,y2.


## 3.5. Pseudo Singular Points

Pseudo singular points are defined as singular points of the “normalized slow dynamics,” that is, as the set of points for which we have(22a)det⁡Jy1,y2=0,(22b)∂g1∂x1x˙1+∂g1∂x2x˙2∂g2∂y2-∂g2∂x1x˙1+∂g2∂x2x˙2∂g1∂y2=0,(22c)∂g1∂x1x˙1+∂g1∂x2x˙2∂g2∂y1-∂g2∂x1x˙1+∂g2∂x2x˙2∂g1∂y1=0,(22d)g1x1,x2,y1,y2=0,(22e)g2x1,x2,y1,y2=0.Remark 1.
Let us notice on one hand that (22b) and (22c) are linearly dependent and on the other hand that contrary to previous works we do not use the “desingularized vector field” (21) but the “normalized slow dynamics” (20).The Jacobian matrix of system (20) reads(23)JF1,F2,G1,G2=∂F1∂x1∂F1∂x2∂F1∂y1∂F1∂y2∂F2∂x1∂F2∂x2∂F2∂y1∂F2∂y2∂G1∂x1∂G1∂x2∂G1∂y1∂G1∂y2∂G2∂x1∂G2∂x2∂G2∂y1∂G2∂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 the three-dimensional case to the four-dimensional case. So, first, let us suppose that by a “standard translation” thepseudo singular point can be shifted at the origin O(0,0,0,0) and that by a “standard rotation” of y1-axis that theslow manifold is tangent to (x2,x3,y1)-hyperplane, so we have(24)f10,0,0,0=g10,0,0,0=0,∂g1∂x20,0,0,0=∂g1∂x30,0,0,0=∂g1∂y10,0,0,0=0.Then, let us make the following assumptions for the nondegeneracy of thefolded singularity:(25)f20,0,0,0≠0;∂g1∂x10,0,0,0≠0;∂2g1∂y120,0,0,0≠0.According to these generic hypotheses ((24)-(25)), the Jacobian matrix (23) reads(26)JF1,F2,G1,G2=0000-f2∂P∂x1-f2∂P∂x2-f2∂P∂y1-f2∂P∂y2a31a32a33a34a41a42a43a44,where(27)P=det⁡Jy1,y2,a3i=-f2∂g2∂x2∂2g1∂y2∂xi+∂g2∂y2f2∂2g1∂x2∂xi+∂g1∂x1∂f1∂xifori=1,2,a3i=-f2∂g2∂x2∂2g1∂y2∂yi+∂g2∂y2f2∂2g1∂x2∂yi+∂g1∂x1∂f1∂yifori=3,4,a4i=f2∂g2∂x2∂2g1∂y1∂xi-∂g2∂y1f2∂2g1∂x2∂xi+∂g1∂x1∂f1∂xifori=1,2,a4i=f2∂g2∂x2∂2g1∂y1∂yi-∂g2∂y1f2∂2g1∂x2∂yi+∂g1∂x1∂f1∂yifori=3,4.Thus, we have the following Cayley-Hamilton eigenpolynomial associated with such Jacobian matrix (26) evaluated at thepseudo singular point, that is, at the origin:(28)λ4-σ1λ3+σ2λ2-σ3λ+σ4=0,where σ1=Tr(J) is the sum of all first-order diagonal minors of J, that is, the trace of the Jacobian matrix J, σ2 represents the sum of all second-order diagonal minors of J, and σ3 represents the sum of all third-order diagonal minors of J. It appears that σ4=|J|=0 since one row of the Jacobian matrix (26) is null. So, the eigenpolynomial reduces to(29)λλ3-σ1λ2+σ2λ-σ3=0.But, according to Wechselberger [11], σ3 vanishes at apseudo singular point as it is easy to prove it. So, the eigenpolynomial (29) is reduced to(30)λ2λ2-σ1λ+σ2=0.Letλi be the eigenvalues of eigenpolynomial (30) and let us denote by λ3,4=0 the obvious double root of this polynomial. We have(31)σ1=TrJF1,F2,G1,G2=λ1+λ2=∂g2∂x1∂g1∂y1∂f1∂y2,σ2=∑i=13JF1,F2,G1,G2ii=λ1λ2=∂g1∂y12f22∂2g2∂x22∂2g2∂y22-∂2g2∂x2∂y22+f2∂g2∂x1∂2g2∂y22∂f1∂x2-∂2g2∂x2∂y2∂f1∂y2,where σ1=Tr(J(F1,F2,G1,G2))=p is is the sum of all first-order diagonal minors of J(F1,F2,G1,G2), that is, the trace of the Jacobian matrix J(F1,F2,G1,G2), and σ2=∑i=13J(F1,F2,G1,G2)ii=q represents the sum of all second-order diagonal minors of J(F1,F2,G1,G2). Thus, thepseudo singular point is of saddle-typeiff the following conditions C1 and C2 are verified:(32)C1:Δ=p2-4q>0,C2:q<0.ConditionC1 is systematically satisfied provided that condition C2 is verified. Thus, thepseudo singular point is of saddle-typeiffq<0.


## 3.7. Canard Existence inR2+2

Following the works of Wechselberger [11] it can be stated, while using a standard polynomial change of variables, that any n-dimensional singularly perturbed systems with kslow variables (k⩾2) and mfast (m⩾1) (1) can be transformed into the following “normal form”:(33)x1˙=a~x2+b~y2+Ox1,ϵ,x22,x2y2,y22,x2˙=1+Ox1,x2,y2,ϵ,ϵy1˙=c~y1+Oϵx1,ϵx2,ϵy2,x12,x22,y22,x2y2,ϵy2˙=-x1+y22+Oϵx1,ϵx2,ϵy2,ϵ2,x12y2,y23,x1x2y2.We establish in the Appendix for any four-dimensional singularly perturbed systems (13a), (13b), (13c), and (13d) with k=2slow and m=2fast variables that(34)a~=12f22∂2g2∂x22∂2g2∂y22-∂2g2∂x2∂y22+f2∂g2∂x1∂2g2∂y22∂f1∂x2-∂2g2∂x2∂y2∂f1∂y2,b~=-∂g2∂x1∂f1∂y2,c~=∂g1∂y1.Thus, in his paper Wechselberger [11, p. 3304] provided in the framework of “standard analysis” a generalization of Benoît’s theorem [2] for any n-dimensional singularly perturbed systems with kslow variables (k⩾2) and mfast (m⩾1). According to his theorem presented below as Theorem 2 he proved the existence of canard solutions for the original system (1).Theorem 2.
In the folded saddle case of system (33) singular canards perturb to maximal canards solutions for sufficiently small ε≪1.Proof.
See Wechselberger [11].Since our method does not use the “desingularized vector field” (21) but the “normalized slow dynamics” (20), we have the following proposition.Proposition 3.
If the normalized slow dynamics (20) has apseudo singular point of saddle-type, that is, if the sum σ2 of all second-order diagonal minors of the Jacobian matrix of the normalized slow dynamics (20) evaluated at thepseudo singular point is negative, that is, if σ2<0 then, according to Theorem 2, system (13a), (13b), (13c), and (13d) exhibits a canard solution which evolves from the attractive part of the slow manifold towards its repelling part.Proof.
By making some smooth changes of time and smooth changes of coordinates (see the Appendix) we brought system (13a), (13b), (13c), and (13d) to the following “normal form”:(35)x1˙=a~x2+b~y2+Ox1,ϵ,x22,x2y2,y22,x2˙=1+Ox1,x2,y2,ϵ,ϵy1˙=c~y1+Oϵx1,ϵx2,ϵy2,x12,x22,y22,x2y2,ϵy2˙=-x1+y22+Oϵx1,ϵx2,ϵy2,ϵ2,x12y2,y23,x1x2y2.
Then, we deduce that the condition for thepseudo singular point to be of saddle-type is a~<0. According to (32) it is easy to verify that(36)σ1=TrJF1,F2,G1,G2=λ1+λ2=-b~c~,σ2=∑i=13JF1,F2,G1,G2ii=λ1λ2=2a~c~2.
So, the condition for which thepseudo singular point is of saddle-type, that is, σ2<0, is identical to that proposed by Wechselberger [11, p. 3298] in his theorem; that is, a~<0.So, Proposition3 can be used to state the existence of canard solution for such systems. Application of Proposition 3 to the coupled FitzHugh-Nagumo equations, presented in Section 4, which is a four-dimensional singularly perturbed system with twoslow and twofast variables will enable proving, as many previous works such as those of Tchizawa and Campbell [18] and Tchizawa [19–24], the existence of “canard solutions” in such system. According to Tchizawa [25], it is very important to notice on one hand that the fast equation has 2 dimensions in the system R2+2 and on the other hand that the fast system can give attractive, repulsive, or attractive-repulsive part at eachpseudo singular point. Then, Tchizawa [25] has established that the jumping direction can be shown using the eigenvectors. In the same way we will find again the results of Rubin and Wechselberger [26] concerning the existence of “canard solutions” in the Hodgkin-Huxley model but with a set of more realistic parameters used in Chua et al. [27, 28].


## 4. Coupled FitzHugh-Nagumo Equations

The FitzHugh-Nagumo model [29, 30] is a simplified version of the Hodgkin-Huxley model [31] which models in a detailed manner activation and deactivation dynamics of a spiking neuron. By coupling two FitzHugh-Nagumo models Tchizawa and Campbell [18] and Tchizawa [19, 24] obtained the following four-dimensional singularly perturbed system with twoslow and twofast variables:(37a)dx1dt=1cy1+bx1,(37b)dx2dt=1cy2+bx2,(37c)εdy1dt=x1-y133+y2,(37d)εdy2dt=x2-y233+y1,where 0<ε≪1 and b is the “canard parameter” or “duck parameter” while c is a scale factor.


### 4.1. Slow Manifold and Constrained System

The slow manifold equation of system (37a), (37b), (37c), and (37d) is defined by setting ε=0 in (37c) and (37d). Thus, we obtain(38)dx1dt=1cy1+bx1,dx2dt=1cy2+bx2,dy1dt=-1/cy2+bx2+y22/cy1+bx1y12y22-1,dy2dt=-1/cy1+bx1+y12/cy2+bx2y12y22-1,0=x1-y133+y2,0=x2-y233+y1.


### 4.2. Normalized Slow Dynamics

Then, by rescaling the time by settingt=-det⁡J(y1,y2)τ=-(y12y22-1) we obtain the “normalized slow dynamics”:(39)dx1dt=-1cy1+bx1y12y22-1=F1x1,x2,y1,y2,dx2dt=-1cy2+bx2y12y22-1=F2x1,x2,y1,y2,dy1dt=1cy2+bx2+y22cy1+bx1=G1x1,x2,y1,y2,dy2dt=1cy1+bx1+y12cy2+bx2G2x1,x2,y1,y2,0=x1-y133+y2,0=x2-y233+y1.


### 4.3. Pseudo Singular Points

From (22a), (22b), (22c), (22d), and (22e), thepseudo singular points of system (37a), (37b), (37c), and (37d) are defined by(40a)det⁡Jy1,y2=y12y22-1=0,(40b)∂g1∂x1x˙1+∂g1∂x2x˙2∂g2∂y2-∂g2∂x1x˙1+∂g2∂x2x˙2∂g1∂y2=1cy2+bx2+y22cy1+bx1=0,(40c)∂g1∂x1x˙1+∂g1∂x2x˙2∂g2∂y1-∂g2∂x1x˙1+∂g2∂x2x˙2∂g1∂y1=1cy1+bx1+y12cy2+bx2=0,(40d)g1x1,x2,y1,y2=x1-y133+y2=0,(40e)g2x1,x2,y1,y2=x2-y233+y1=0.According to Tchizawa and Campbell [18] and Tchizawa [19, 20], there are sixpseudo singular points with the last four depending on the parameter b:(41a)x~1,x~2,y~1,y~2=±43,∓43,±1,∓1,(41b)x~1,x~2,y~1,y~2=±3-9-4b2/b3+29-4b232b,∓3-9-4b2/b9-8b2+39-4b262b2,∓3-9-4b22b,∓2b3-9-4b2,(41c)x~1,x~2,y~1,y~2=±3-29-4b23+9-4b2/b32b,∓3+9-4b2/b9-8b2-39-4b262b2,∓3-9-4b22b,∓2b3-9-4b2.


### 4.4. Canard Existence in Coupled FitzHugh-Nagumo Equations

The Jacobian matrix of system (39) evaluated at thepseudo singular points (41a) reads(42)JF1,F2,G1,G2=0023+4b3c-23+4b3c00-23+4b3c23+4b3cbcbc1c-3+8b3cbcbc-3+8b3c1c.Remark 4.
Although thepseudo singular points have not been shifted at the origin extension of Benoît’s generic hypotheses (24)-(25) are satisfied. In other words, we have σ4=σ3=0.According to (31) we find that(43)p=σ1=TrJ=+2c,q=σ2=-16b3+4b9c2.Thus, according to Proposition3, thepseudo singular points are of saddle-type if and only if:(44)-16b3+4b9c2<0.So, we have conditionsC1 and C2:(45)C1:Δ=43+8b29c2>0,C2:q=-16b3+4b9c2<0.Let us choose arbitrarilyb as the “canard parameter” or “duck parameter.” Obviously, it appears that the condition C1 is still satisfied. Finally, thepseudo singular points are of saddle-type if and only if we have(46)b>0orb<-34.Remark 5.
Let us notice that thepseudo singular points are of node-type if -3/4<b<0 as stated by Tchizawa and Campbell [18] and Tchizawa [19, 20].The Jacobian matrixJ(F1,F2,G1,G2) of system (39) evaluated at thepseudo singular points (41b) reads(47)00-49-4b23c-2-9+4b2+39-4b23bc002-9+4b2+39-4b23bc49-4b23c3+9-4b22cbc3+9-4b22bc3-49-4b23cbc3-9-4b22c3+49-4b23c3-9-4b22bc.Remark 6.
Although thepseudo singular points have not been shifted at the origin extension of Benoît’s generic hypotheses (24)-(25) are satisfied. In other words, we have σ4=σ3=0.According to (31) we find that(48)p=σ1=TrJ=+3bc,q=σ2=169-4b29c2.Thus, according to Proposition3, thepseudo singular points are of saddle-type if and only if(49)169-4b29c2<0Δ=p2-4q>0,q<0.So, we have conditionsC1 and C2:(50)C1:Δ=3bc2-649-4b29c2>0,C2:q=169-4b29c2<0.Let us choose arbitrarilyb as the “canard parameter” or “duck parameter.” Obviously, it appears that if the condition C2 is verified, then the condition C1 isde facto satisfied. Finally, thepseudo singular points are of saddle-type if and only if we have(51)b>32orb<-32.Remark 7.
Because of the symmetry of these coupled FitzHugh-Nagumo equations, the Jacobian matrix of system (39) evaluated at thepseudo singular points (41c) provides the same result as just above.


## 4.1. Slow Manifold and Constrained System

The slow manifold equation of system (37a), (37b), (37c), and (37d) is defined by setting ε=0 in (37c) and (37d). Thus, we obtain(38)dx1dt=1cy1+bx1,dx2dt=1cy2+bx2,dy1dt=-1/cy2+bx2+y22/cy1+bx1y12y22-1,dy2dt=-1/cy1+bx1+y12/cy2+bx2y12y22-1,0=x1-y133+y2,0=x2-y233+y1.


## 4.2. Normalized Slow Dynamics

Then, by rescaling the time by settingt=-det⁡J(y1,y2)τ=-(y12y22-1) we obtain the “normalized slow dynamics”:(39)dx1dt=-1cy1+bx1y12y22-1=F1x1,x2,y1,y2,dx2dt=-1cy2+bx2y12y22-1=F2x1,x2,y1,y2,dy1dt=1cy2+bx2+y22cy1+bx1=G1x1,x2,y1,y2,dy2dt=1cy1+bx1+y12cy2+bx2G2x1,x2,y1,y2,0=x1-y133+y2,0=x2-y233+y1.


## 4.3. Pseudo Singular Points

From (22a), (22b), (22c), (22d), and (22e), thepseudo singular points of system (37a), (37b), (37c), and (37d) are defined by(40a)det⁡Jy1,y2=y12y22-1=0,(40b)∂g1∂x1x˙1+∂g1∂x2x˙2∂g2∂y2-∂g2∂x1x˙1+∂g2∂x2x˙2∂g1∂y2=1cy2+bx2+y22cy1+bx1=0,(40c)∂g1∂x1x˙1+∂g1∂x2x˙2∂g2∂y1-∂g2∂x1x˙1+∂g2∂x2x˙2∂g1∂y1=1cy1+bx1+y12cy2+bx2=0,(40d)g1x1,x2,y1,y2=x1-y133+y2=0,(40e)g2x1,x2,y1,y2=x2-y233+y1=0.According to Tchizawa and Campbell [18] and Tchizawa [19, 20], there are sixpseudo singular points with the last four depending on the parameter b:(41a)x~1,x~2,y~1,y~2=±43,∓43,±1,∓1,(41b)x~1,x~2,y~1,y~2=±3-9-4b2/b3+29-4b232b,∓3-9-4b2/b9-8b2+39-4b262b2,∓3-9-4b22b,∓2b3-9-4b2,(41c)x~1,x~2,y~1,y~2=±3-29-4b23+9-4b2/b32b,∓3+9-4b2/b9-8b2-39-4b262b2,∓3-9-4b22b,∓2b3-9-4b2.


## 4.4. Canard Existence in Coupled FitzHugh-Nagumo Equations

The Jacobian matrix of system (39) evaluated at thepseudo singular points (41a) reads(42)JF1,F2,G1,G2=0023+4b3c-23+4b3c00-23+4b3c23+4b3cbcbc1c-3+8b3cbcbc-3+8b3c1c.Remark 4.
Although thepseudo singular points have not been shifted at the origin extension of Benoît’s generic hypotheses (24)-(25) are satisfied. In other words, we have σ4=σ3=0.According to (31) we find that(43)p=σ1=TrJ=+2c,q=σ2=-16b3+4b9c2.Thus, according to Proposition3, thepseudo singular points are of saddle-type if and only if:(44)-16b3+4b9c2<0.So, we have conditionsC1 and C2:(45)C1:Δ=43+8b29c2>0,C2:q=-16b3+4b9c2<0.Let us choose arbitrarilyb as the “canard parameter” or “duck parameter.” Obviously, it appears that the condition C1 is still satisfied. Finally, thepseudo singular points are of saddle-type if and only if we have(46)b>0orb<-34.Remark 5.
Let us notice that thepseudo singular points are of node-type if -3/4<b<0 as stated by Tchizawa and Campbell [18] and Tchizawa [19, 20].The Jacobian matrixJ(F1,F2,G1,G2) of system (39) evaluated at thepseudo singular points (41b) reads(47)00-49-4b23c-2-9+4b2+39-4b23bc002-9+4b2+39-4b23bc49-4b23c3+9-4b22cbc3+9-4b22bc3-49-4b23cbc3-9-4b22c3+49-4b23c3-9-4b22bc.Remark 6.
Although thepseudo singular points have not been shifted at the origin extension of Benoît’s generic hypotheses (24)-(25) are satisfied. In other words, we have σ4=σ3=0.According to (31) we find that(48)p=σ1=TrJ=+3bc,q=σ2=169-4b29c2.Thus, according to Proposition3, thepseudo singular points are of saddle-type if and only if(49)169-4b29c2<0Δ=p2-4q>0,q<0.So, we have conditionsC1 and C2:(50)C1:Δ=3bc2-649-4b29c2>0,C2:q=169-4b29c2<0.Let us choose arbitrarilyb as the “canard parameter” or “duck parameter.” Obviously, it appears that if the condition C2 is verified, then the condition C1 isde facto satisfied. Finally, thepseudo singular points are of saddle-type if and only if we have(51)b>32orb<-32.Remark 7.
Because of the symmetry of these coupled FitzHugh-Nagumo equations, the Jacobian matrix of system (39) evaluated at thepseudo singular points (41c) provides the same result as just above.


## 5. Hodgkin-Huxley Model

The original Hodgkin-Huxley model [31] is described by the following system of four nonlinear ordinary differential equations:(52a)dVdt=1CMI-g-Kn4V-VK-g-Nam3hV-VNa-g-LV-VL,(52b)dndt=αnV1-n-βnVn,(52c)dmdt=αmV1-m-βmVm,(52d)dhdt=αhV1-h-βhVh,where:(53a)αnV=0.01V+10exp⁡V+10/10⁡-1,(53b)βnV=0.125exp⁡V80,(53c)αmV=0.1V+25exp⁡V+25/10-1,(53d)βmV=4exp⁡V18,(53e)αhV=0.07exp⁡V20,(53f)βnV=1exp⁡V+30/10+1.The first equation (52a) results from the application of Kirchhoff’s law to the space clamped squid giant axon. Thus, the total membrane current CMdV/dt for which CM represents the specific membrane capacity and V the displacement of the membrane potential from its resting value is equal to the sum of the following intrinsic currents:(54)IK=g-Kn4V-VK,INa=g-Nam3hV-VNa,IL=g-LV-VL,where IK is a delayed rectifier potassium current, INa is fast sodium current, and IL is the “leakage current.” The parameter I is the total membrane current density, inward positive, that is, the total current injected into the space clamped squid giant axon, and VK, VNa, and VL are the equilibrium potentials of potassium, sodium, and “leakage current,” respectively. The maximal specific conductance degrees of the ionic currents are denoted by g-K, g-Na, and g-L, respectively. Functions αn,m,h and βn,m,h are gates’ opening and closing rates depending on V. Variable m denotes the activation of the sodium current, variable h the inactivation of the sodium current, and variable n the activation of the potassium current. These dimensionless gating variables vary in the range [0,1].Let us notice that the variables and symbols in ((52a), (52b), (52c), and (52d) and (53a), (53b), (53c), (53d), (53e), and (53f)) originally chosen by Hodgkin and Huxley and are different from those found in recent literatures, where the reference polarity of the voltage V and the reference direction of the current I are defined as the negative of the voltages and currents. We have opted to adopt the reference assumption in Hodgkin and Huxley [31] for ease in comparison of our results with those from Hodgkin and Huxley (for more details see Chua et al. [27, 28]). The parameter values are exactly those chosen in the original Hodgkin-Huxley [31] works:(55)CM=1.0μF/cm2,VNa=-115mV,VK=12mV,VL=-10.613mV,g-Na=120mS/cm2,g-K=36mS/cm2,g-L=0.3mS/cm2.According to Suckley and Biktashev [32] and Suckley [33], dimensionless functions n-, h-, and m- called gates’ instant equilibrium values, that is, steady-state relation for gating variables n, h, and m, respectively, as well as τn, τh, and τm called gates dynamics timescales in ms, that is, time constant for gating variables n, h, and m, respectively, may be defined as follows:(56a)n-V=αnVαnV+βnV,(56b)h-V=αhVαhV+βhV,(56c)m-V=αmVαmV+βmV,(56d)τnV=1αnV+βnV,(56e)τhV=1αhV+βhV,(56f)τmV=1αmV+βmV.By using (56a), (56b), (56c), (56d), (56e), and (56f), the original Hodgkin-Huxley model [31] reads(57a)dVdt=1CMI-g-Kn4V-VK-g-Nam3hV-VNa-g-LV-VL,(57b)dndt=n--nτn,(57c)dhdt=h--hτh,(57d)dmdt=m--mτm.Now, in order to apply thesingular perturbation method to the Hodgkin-Huxley model, two small multiplicative parameters ε≪1 are introduced. According to Suckley and Biktashev [32], Suckley [33], and Rubin and Wechselberger [26], the existence of two different timescales of evolution for couples of dynamic variables (n,h) and (m,V) enables justifying such an introduction. So, in order to differentiateslow variables fromfast variables, Suckley and Biktashev [32], Suckley [33], and Rubin and Wechselberger [26] have plotted the inverse of “time constant for gating variable i,” that is, τi-1 according to V with i=n,h,m. In Figure 1, they have been plotted for the original functions αi and βi (53a). However, let us notice that this plot is exactly the same as those presented by Rubin and Wechselberger [26] (Figure 1) for a nondimensionalized three-dimensional Hodgkin-Huxley singularly perturbed system obtained after the following variable changes: V→-V and I-→-I-; then V→V+65 and finally V→V/100.Figure 1
Graph of1/τi (ms−1) against V (mV).Figure1 shows a plot of the functions τi-1 according to V with i=n,h,m over the physiological range. We observe that τm-1 is of an order of magnitude bigger than τh-1 and τn-1, which are of comparable size. Indeed, we can deduce that the values of times scales are approximately τm-1≈10ms-1 while τn-1≈τh-1≈1ms-1. Then, it appears that m corresponds to thefast variable while n and h correspond toslow variables. Moreover, since the activation of the sodium channel m is directly related to the dynamics of the membrane (action) potential V, Rubin and Wechselberger [26] consider that m and V evolve on the samefast timescale. So, the Hodgkin-Huxley model may be transformed into asingularly perturbed system with two timescales in which theslow variables are n,h and thefast variables are m,V.So, according to Awiszus et al. [34], Suckley and Biktashev [32], Suckley [33], and Rubin and Wechselberger [26] small multiplicative parameters 0<ε≪1 in the original vector field of the Hodgkin-Huxley equations (57a), (57b), (57c), and (57d) may be identified while factorizing the right-hand side of (57a) by g-Na and set:(58)g-Na⟶g-Nag-Na=1,g-K⟶g-Kg-Na=0.3,g-L⟶g-Lg-Na=0.0025.Other parameters are kept as for the original Hodgkin-Huxley model [31]:(59)CM=1.0μF/cm2,V-Na=-115mV,V-K=12mV,V-L=-10.613mV.Then, by posingI-→I-/g-Na, ε=CM/g-Na=1/120, and (n,h,m,V) = (x1,x2,y1,y2) to be consistent with the notations of Section 3, we obtain(60a)dx1dt=x-1-x1τ1=f1x1,x2,y1,y2,(60b)dx2dt=x-2-x2τ2=f1x1,x2,y1,y2,(60c)εdy1dt=y-1-y1τ3=g1x1,x2,y1,y2,(60d)εdy2dt=I--g-Kx14y2-VK-g-Nay13x2y2-VNa-g-Ly2-VL=g2x1,x2,y1,y2,where (x-1,x-2,y-1) = (n-,h-,m-) and τ1,2,3=τn,h,m.Let us notice that the multiplicative parameterε has been introduced artificially in (60c). This is due to the fact that it has been stated above that the timescale of variable m, that is, y1, is tenth times greater than the timescale of variables n and h, that is, of variables x1 and x2. Moreover, this parameter is identical to those used in (60d) since it has been also considered that m and V, that is, y1 and y2, evolve on the samefast timescale.According to theGeometric Singular Perturbation Theory, the zero-order approximation in ε of theslow manifold associated with the Hodgkin-Huxley model (60a), (60b), (60c), and (60d) is obtained by posing ε=0 in (60c) and (60d). So, theslow manifold is given by(61a)x2=I--g-Kx14y2-VK-g-Ly2-VLg-Nay-13y2-VNa,(61b)y1=y-1y2.Then, thefast foliation is within the planes x1=constant and x2=constant.Thefold curve is defined as the location of the points, where g1x1,x2,y1,y2=0, g2x1,x2,y1,y2=0, and det⁡J(g1,g2)=0. For the Hodgkin-Huxley model (60a), (60b), (60c), and (60d), thefold curve is thus given by (61a) and (61b) and by the determinant of the Jacobian matrix of the following fast foliation:(62a)dy1dt=y-1-y1τ3=g1x1,x2,y1,y2,(62b)dy2dt=I--g-Kx14y2-VK-g-Nay13x2y2-VNa-g-Ly2-VL=g2x1,x2,y1,y2.The Jacobian matrix of thefast foliation (62a) and (62b) reads(63)Jg1,g2=y1-′τ3-τ3′y1--y1τ32-1τ3-g-Kx4+g-Nay13x2+g-L-3g-Nay12x2y2-VNa,where (′) denotes the derivative with respect to y2. Then, taking into account (61b), that is, y1=y-1, we have(64)Jg1,g2=y1-′τ3-1τ3-g-Kx14+g-Nay-13x2+g-L-3g-Nay-12x2y2-VNa.So, the determinant of the Jacobian matrix of thefast foliation (62a) and (62b) is(65)det⁡Jg1,g2=-1τ3g-Kx14+g-Nay-13x2+g-L+3g-Nay1-′y-12x2y2-VNa.Thus, the condition for thefold curve is det⁡J(g1,g2)=0, which gives(66)g-Kx14+g-Nay-13x2+g-L+3g-Nay1-′y-12x2y2-VNa=0.Therefore(67)x2=-g-Kx14+g-Lg-Nay-12y-1+3y1-′y2-VNa.By subtracting (61a) from (67) we obtain x1(68)x1=x1f=-I-y-1+3y1-′y2-VNa+g-LVNa-VLy-1+3y1-′y2-VNay2-VLg-KVK-VNay-1-3y1-′y2-VNay2-VK1/4.Plugging this value ofx1 (68) into (67) provides(69)x2=x2f=I-+g-LVK-VLg-Nay-12VNa-VKy-1+3y1-′y2-VNay2-VK.So, thefold curve is given by the set of parametric equations ((68)-(69)) in terms of y2.Thepseudo singular points are given by (22a), (22b), (22c), (22d), and (22e) which read for the Hodgkin-Huxley model (60a), (60b), (60c), and (60d):(70a)y-1-y1τ3=0,(70b)I--g-Kx14y2-VK-g-Nay13x2y2-VNa-g-Ly2-VL=0,(70c)4g-Kx13y2-VKx1-x-1τ1+g-Nay13y2-VNax2-x-2τ2=0,(70d)4g-Kx13y2-VKx1-x-1τ1+g-Nay13y2-VNax2-x-2τ21τ3=0,(70e)τ3g-Kx14+g-Nay13x2+g-L+3g-Nay12x2y2-VNaτ3y1′+y1-y-1τ3′=0.Let us notice that (70c) and (70d) are identical. Moreover, the definition of τ3 (56f) enables simplifying system (70a), (70b), (70c), (70d), and (70e). Thus, we have(71a)I--g-Kx14y2-VK-g-Nay-13x2y2-VNa-g-Ly2-VL=0,(71b)4g-Kx13y2-VKx1-x-1τ1+g-Nay-13y2-VNax2-x-2τ2=0,(71c)g-Kx14+g-Nay-13x2+g-L+3g-Nay-12y-1′x2y2-VNa=0.Moreover, (71a) and (71c) indicate that thepseudo singular point belongs to theslow manifold and to thefold curve. So, let us replace in (71b) the variables x1 and x2 by the variables x1f and x2f given by (68) and (69), respectively, which represent the parametric equations offold curve(72)4gKx1f3y2-VKx1f-x-1τ1+y-13y2-VNax2f-x-2τ2=0.Thus, it appears that (72) depends on the variable y2, on the functions gates dynamics timescales τ1(y2) and τ2(y2), and on the bifurcation parameter I-. According to Rubin and Wechselberger [26], the function y2(I-), solution of (72), is independent of time multiplicative constants k1 and k2 that one could set in factor of τ1(y2) and τ2(y2).So, following their works, let us plot the functiony2(I-) solution of (72) for various values of these time constants by posing successively in (72) k1=1, 3, 4.75 and 7 and while fixing k2=1. The result is presented in Figure 2.Figure 2
Functiony2(I-) for various values of parameter kn=1,3,4.75,7 exhibiting the the bifurcation parameter value  I-C≈-4.8.Let us notice that this plot (the functiony2(I-) solution of (71a), (71b), and (71c) has been plotted with Mathematica© while using the Contour Plot Function used for representing implicit function since such function cannot be expressed explicitly) is exactly the same as those presented by Rubin and Wechselberger [26] (Figure 8-9) for a nondimensionalized three-dimensional Hodgkin-Huxley singularly perturbed system which had been obtained after the following variable changes: V→-V and I-→-I; then V→V+65 and finally V→V/100.We observe from Figure2 that the bifurcation parameter value I-C≈-4.8 is exactly identical (in absolute value) to those obtained by Rubin and Wechselberger [26]. Numerical resolution (this resolution has been made while using the function FindRoot in Mathematica©) of (72) provides a better approximation of the bifurcation parameter value:(73)I-C=-4.82988μA.This value corresponds to a voltagey2=-3.18136mV.ForI-≈-4.1, the coordinate of thepseudo singular point can be computed numerically:(74)x1,x2,y1,y2=0.362513,0.521793,0.0733782,-2.81908.According to Proposition3 we can state that the eigenpolynomial of the Jacobian matrix associated with the “normalized slow dynamics” of the Hodgkin-Huxley model (60a), (60b), (60c), and (60d) reads(75)λ4-σ1λ3+σ2λ2-σ3λ+σ4=0for which it is easy to prove that σ4=σ3=0. So, this eigenpolynomial reduces to(76)λ2λ2-σ1λ+σ2=0.According to (31) we find that(77)p=TrJ=144.933,q=σ2=-362.924.Thus, according to Proposition3, thepseudo singular point is of saddle-type. Moreover, numerical computation of the eigenvalues of this Jacobian matrix evaluated at thepseudo singular point provides(78)λ1,λ2,λ3,λ4=-2.46224,147.396,0,0.So, according to Proposition3, thispseudo singular point is of saddle-type and canard solution may occur in the four-dimensional Hodgkin-Huxley singularly perturbed system (60a), (60b), (60c), and (60d) for the original set of parameter values.In Figures3, 4, and 5 canard solution of the four-dimensional Hodgkin-Huxley singularly perturbed system for the “canard value” of I-≈-4.1 has been plotted in the (x1,x2,y2) phase space and then in the (x1,y2) phase plane. The green point represents thepseudo singular point. The trajectory curve, that is, the canard solution, has been plotted in red while thefold curve is in yellow. We observe in Figure 3 that when the trajectory curve reaches thefold at thepseudo singular point it “jumps” suddenly to the other part of theslow manifold before being reinjected towards thepseudo singular point.Figure 3
Phase portrait, canard solution, andslow manifold of the Hodgkin-Huxley system (60a), (60b), (60c), and (60d) in the (n,h,V) phase space.Figure 4
Phase portrait, canard solution, andslow manifold of the Hodgkin-Huxley system (60a), (60b), (60c), and (60d) in the (V,n) phase plane.Figure 5
Phase portrait, canard solution, andslow manifold of the Hodgkin-Huxley system (60a), (60b), (60c), and (60d) in the (V,h) phase plane.


## 6. Discussion

In a previous paper entitled “Canards Existence in Memristor’s Circuits” (see Ginoux and Llibre [12]) we have proposed a new method for proving the existence of “canard solutions” for three- and four-dimensional singularly perturbed systems with only onefast variable which improves the methods used until now. This method enabled stating a unique “generic” condition for the existence of “canard solutions” for such three- and four-dimensional singularly perturbed systems which is based on the stability offolded singularities of thenormalized slow dynamics deduced from a well-known property of linear algebra. This unique condition which is completely identical to that provided by Benoît [2] and then by Szmolyan and Wechselberger [5] and finally by Wechselberger [11] was considered as “generic” since it was exactly the same for singularly perturbed systems of dimensions three and four with only onefast variable. In this work we have extended this new method to the case of four-dimensional singularly perturbed systems with twoslow and twofast 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); we have highlighted and provide a simple and efficient tool for testing the occurrence of “canard solutions” in any three- or four-dimensional singularly perturbed systems with one or twofast variables. Applications of this method to the famous coupled FitzHugh-Nagumo equations and to the Hodgkin-Huxley model has enabled showing the existence of “canard solutions” in such systems. However, in this paper, only the case ofpseudo singular points orfolded singularities of saddle-type has been analyzed. Of course, the case of ofpseudo singular points orfolded singularities of node-type and focus-type could be also studied with the same method.



---

*Source: 342010-2015-12-30.xml*