About this Journal Submit a Manuscript Table of Contents
Mathematical Problems in Engineering
Volume 2015 (2015), Article ID 342010, 17 pages
http://dx.doi.org/10.1155/2015/342010
Research Article

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

1Laboratoire LSIS, CNRS, UMR 7296, Université de Toulon, BP 20132, 83957 La Garde Cedex, France
2Departament de Matemàtiques, Universitat Autònoma de Barcelona, 08193 Bellaterra, Barcelona, Spain

Received 5 August 2015; Accepted 3 December 2015

Academic Editor: Matjaz Perc

Copyright © 2015 Jean-Marc Ginoux and Jaume Llibre. This is an open access article distributed under the Creative Commons Attribution License, which permits unrestricted use, distribution, and reproduction in any medium, provided the original work is properly cited.

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 one fast 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 two slow and two fast 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 of folded singularities (pseudo singular points in this case) of the normalized 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.

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 . In the article entitled “Systèmes Lents-Rapides dans et Leurs Canards,” Benoît [2, p. 170] proved the existence of canards solution for three-dimensional singularly perturbed systems with two slow variables and one fast variable while using “nonstandard analysis” according to a theorem which stated that canard solutions exist in such systems provided that the pseudo singular point (this concept has been originally introduced by Argémi [4]; see Section 2.8) of the slow dynamics, that is, of the reduced vector field, is of saddle-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 and provided a “standard version” of Benoît’s theorem [2]. Very recently, Wechselberger [11] generalized this theorem for -dimensional singularly perturbed systems with slow variables and fast (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 the normal form of such singularly perturbed systems which is expressed according to some coefficients ( and for dimension three and , , and 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 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 Argémi [4, p. 336]). This leads Szmolyan and Wechselberger [5, p. 427] and Wechselberger [11, p. 3298] to a “classification of folded singularities (pseudo singular points).” Thus, they showed that for three-dimensional singularly perturbed systems such folded singularity is of saddl-type if the following condition is satisfied, while for four-dimensional singularly perturbed systems such folded singularity is of saddle-type if . Then, Szmolyan and Wechselberger [5, p. 439] and Wechselberger [11, p. 3304] established their Theorem . 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 ( and ) 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 the normalized slow dynamics and not from the projection of the “desingularized vector field” of the normal 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 of folded singularities of the normalized 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 one fast variable.

The aim of this work is to extend this method to the case of four-dimensional singularly perturbed systems with slow and fast variables. Since the dimension of the system is , such problem is known as “canards existence in ” Moreover, in this particular case, where , the folded singularities of Wechselberger [11, p. 3298] are nothing else but the pseudo 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 systems pseudo singular points are of saddle-type if . Then, according to Theorem . of Wechselberger [11, p. 3304] we provide the expression of this constant which is necessary to establish the existence of canard solutions in such systems. So, we can state that the condition for existence of canards in such is “generic” since it is exactly the same for singularly perturbed systems of dimensions three and four with only one fast variable.

The outline of this paper is as follows. In Section 1, 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 two fast 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 aswhere , , , and the prime denotes differentiation with respect to the independent variable . The functions and are assumed to be functions (in certain applications these functions will be supposed to be , ) of , , and in , where is an open subset of and is an open interval containing .

In the case when , that is, is a small positive number, the variable is called slow variable and is called fast variable. Using Landau’s notation, represents a function of and such that is bounded for positive going to zero, uniformly for in the given domain.

In general we consider that evolves at an rate while evolves at an slow rate. Reformulating system (1) in terms of the rescaled variable , we obtain

The dot represents the derivative with respect to the new independent variable .

The independent variables and are referred to the fast and slow times, respectively, and (1) and (2) are called the fast and slow systems, respectively. These systems are equivalent whenever , and they are labeled singular perturbation problems when . The label “singular” stems in part from the discontinuous limiting behavior in system (1) as .

2.2. Reduced Slow System

In such case system (2) leads to a differential-algebraic system (DAE) called reduced slow system whose dimension decreases from to . Then, the slow variable partially evolves in the submanifold called the critical manifold (it represents the approximation of the slow invariant manifold, with an error of ). The reduced slow system is

2.3. Slow Invariant Manifold

The critical manifold is defined by

Such a normally hyperbolic invariant manifold (4) of the reduced slow system (3) persists as a locally invariant slow manifold of the full problem (1) for sufficiently small. This locally slow invariant manifold is close to the critical manifold.

When is invertible, thanks to the Implicit Function Theorem, is given by the graph of a function for , where is a compact, simply connected domain and the boundary of is a -dimensional submanifold (the set is overflowing invariant with respect to (2) when ; see Kaper [10] and Jones [9]).

According to Fenichel [6, 7] theory if is sufficiently small, then there exists a function defined on such that the manifold,is locally invariant under the flow of system (1). Moreover, there exist perturbed local stable (or attracting) and unstable (or repelling) branches of the slow invariant manifold . Thus, normal hyperbolicity of is lost via a saddle-node bifurcation of the reduced slow system (3). Then, it gives rise to solutions of “canard” type.

2.4. Canards: Singular Canards and Maximal Canards

A canard is a solution of singularly perturbed dynamical system (1) following the attracting branch 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 of the slow invariant manifold.

A singular canard is a solution of reduced slow system (3) following the attracting branch of the critical manifold, passing near a bifurcation point located on the fold of this critical manifold, and then following the repelling branch of the critical manifold.

A maximal canard corresponds to the intersection of the attracting and repelling branches 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 of the slow manifold and then follow, rather surprisingly, a repelling/saddle branch 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 and . However, in the singular limit , 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 ,    and one fast variable , 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 the reduced slow system (3) (obtained by setting in (2)), Takens [16] introduced the “constrained system” defined as follows:

Since, according to Fenichel [6, 7], the critical manifold may be considered as locally invariant under the flow of system (1), we have

By replacing by leads to

This justifies the introduction of the constrained system.

Now, let denote the adjoint of the matrix which is the transpose of the cofactor matrix ; then while multiplying the left-hand side of (6) by the inverse matrix obtained by the adjoint method we have

2.6. Normalized Slow Dynamics

Then, by rescaling the time by setting we obtain the following system which has been called by Benoît [2, p. 166] “normalized slow dynamics”:where the overdot now denotes the time derivation with respect to .

Let us notice that Argémi [4] proposed to rescale time by setting 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, as a function of the other variables. This implies that is locally the graph of a function over the base , where . Thus, we can span the “normalized slow dynamics” on the tangent bundle at the critical manifold at the pseudo singular point. This leads to the so-called desingularized vector field:

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 such pseudo singular points folded singularities. In a recent publication entitled “A Propos de Canards” Wechselberger [11, p. 3295] proposed to define such singularities for -dimensional singularly perturbed systems with slow variables and fast as the solutions of the following system:

Thus, for dimensions higher than three, his concept encompasses that of Argémi. Moreover, Wechselberger [11, p. 3296] proved that folded singularities form a -dimensional manifold. Thus, for the folded singularities are nothing else than the pseudo singular points defined by Argémi [4], while for the folded singularities are no more points but a -dimensional manifold. Moreover, let us notice on one hand that the original system (1) includes variables and on the other hand that system (12) comprises equations. However, in the particular case , two equations of system (12) are linearly dependent. So, such system only comprises equations. So, all the variables (unknowns) of system (12) can be determined. The solutions of this system are called pseudo singular points. We will see in Section 3 that the stability analysis of these pseudo 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-dimensional singularly perturbed dynamical system (2) with slow variables and fast may be written aswhere , , , and the functions and are assumed to be functions of .

3.1. Critical Manifold

The critical manifold equation of system (13a), (13b), (13c), and (13d) is defined by setting in ((13c) and (13d)). Thus, we obtain

By application of the Implicit Function Theorem, let us suppose that we can explicitly express from ((14a) and (14b)), say, without loss of generality, and as functions of the other variables:

3.2. Constrained System

The constrained system is obtained by equating with zero the time derivative of :

Equations (16a) and (16b) may be written as

By solving the system of (17a) and (17b) with two unknowns we find

So, we have the following constrained system:

3.3. Normalized Slow Dynamics

By rescaling the time by setting we obtain the “normalized slow dynamics”:where the overdot now denotes the time derivation with respect to .

3.4. Desingularized System on the Critical Manifold

Then, since we have supposed that and 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

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

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

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” the pseudo singular point can be shifted at the origin and that by a “standard rotation” of -axis that the slow manifold is tangent to ()-hyperplane, so we have

Then, let us make the following assumptions for the nondegeneracy of the folded singularity:

According to these generic hypotheses ((24)-(25)), the Jacobian matrix (23) readswhere

Thus, we have the following Cayley-Hamilton eigenpolynomial associated with such Jacobian matrix (26) evaluated at the pseudo singular point, that is, at the origin:where is the sum of all first-order diagonal minors of , that is, the trace of the Jacobian matrix , represents the sum of all second-order diagonal minors of , and represents the sum of all third-order diagonal minors of . It appears that since one row of the Jacobian matrix (26) is null. So, the eigenpolynomial reduces to

But, according to Wechselberger [11], vanishes at a pseudo singular point as it is easy to prove it. So, the eigenpolynomial (29) is reduced to

Let be the eigenvalues of eigenpolynomial (30) and let us denote by the obvious double root of this polynomial. We havewhere is is the sum of all first-order diagonal minors of , that is, the trace of the Jacobian matrix , and represents the sum of all second-order diagonal minors of . Thus, the pseudo singular point is of saddle-type iff the following conditions and are verified:

Condition is systematically satisfied provided that condition is verified. Thus, the pseudo singular point is of saddle-type iff .

3.7. Canard Existence in

Following the works of Wechselberger [11] it can be stated, while using a standard polynomial change of variables, that any -dimensional singularly perturbed systems with slow variables () and fast () (1) can be transformed into the following “normal form”:

We establish in the Appendix for any four-dimensional singularly perturbed systems (13a), (13b), (13c), and (13d) with slow and fast variables that

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 -dimensional singularly perturbed systems with slow variables () and fast (). 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 .

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 a pseudo singular point of saddle-type, that is, if the sum of all second-order diagonal minors of the Jacobian matrix of the normalized slow dynamics (20) evaluated at the pseudo singular point is negative, that is, if 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”:Then, we deduce that the condition for the pseudo singular point to be of saddle-type is . According to (32) it is easy to verify thatSo, the condition for which the pseudo singular point is of saddle-type, that is, , is identical to that proposed by Wechselberger [11, p. 3298] in his theorem; that is, .

So, Proposition 3 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 two slow and two fast variables will enable proving, as many previous works such as those of Tchizawa and Campbell [18] and Tchizawa [1924], 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 and on the other hand that the fast system can give attractive, repulsive, or attractive-repulsive part at each pseudo 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 two slow and two fast variables:where and is the “canard parameter” or “duck parameter” while 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 in (37c) and (37d). Thus, we obtain

4.2. Normalized Slow Dynamics

Then, by rescaling the time by setting we obtain the “normalized slow dynamics”:

4.3. Pseudo Singular Points

From (22a), (22b), (22c), (22d), and (22e), the pseudo singular points of system (37a), (37b), (37c), and (37d) are defined by

According to Tchizawa and Campbell [18] and Tchizawa [19, 20], there are six pseudo singular points with the last four depending on the parameter :

4.4. Canard Existence in Coupled FitzHugh-Nagumo Equations

The Jacobian matrix of system (39) evaluated at the pseudo singular points (41a) reads

Remark 4. Although the pseudo 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 .

According to (31) we find that

Thus, according to Proposition 3, the pseudo singular points are of saddle-type if and only if:

So, we have conditions and :

Let us choose arbitrarily as the “canard parameter” or “duck parameter.” Obviously, it appears that the condition is still satisfied. Finally, the pseudo singular points are of saddle-type if and only if we have

Remark 5. Let us notice that the pseudo singular points are of node-type if as stated by Tchizawa and Campbell [18] and Tchizawa [19, 20].

The Jacobian matrix of system (39) evaluated at the pseudo singular points (41b) reads

Remark 6. Although the pseudo 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 .

According to (31) we find that

Thus, according to Proposition 3, the pseudo singular points are of saddle-type if and only if

So, we have conditions and :

Let us choose arbitrarily as the “canard parameter” or “duck parameter.” Obviously, it appears that if the condition is verified, then the condition is de facto satisfied. Finally, the pseudo singular points are of saddle-type if and only if we have

Remark 7. Because of the symmetry of these coupled FitzHugh-Nagumo equations, the Jacobian matrix of system (39) evaluated at the pseudo 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:

where:

The first equation (52a) results from the application of Kirchhoff’s law to the space clamped squid giant axon. Thus, the total membrane current for which represents the specific membrane capacity and the displacement of the membrane potential from its resting value is equal to the sum of the following intrinsic currents:where is a delayed rectifier potassium current, is fast sodium current, and is the “leakage current.” The parameter is the total membrane current density, inward positive, that is, the total current injected into the space clamped squid giant axon, and , , and are the equilibrium potentials of potassium, sodium, and “leakage current,” respectively. The maximal specific conductance degrees of the ionic currents are denoted by , , and , respectively. Functions and are gates’ opening and closing rates depending on . Variable denotes the activation of the sodium current, variable the inactivation of the sodium current, and variable the activation of the potassium current. These dimensionless gating variables vary in the range .

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 and the reference direction of the current 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:

According to Suckley and Biktashev [32] and Suckley [33], dimensionless functions , , and called gates’ instant equilibrium values, that is, steady-state relation for gating variables , , and , respectively, as well as , , and called gates dynamics timescales in , that is, time constant for gating variables , , and , respectively, may be defined as follows:

By using (56a), (56b), (56c), (56d), (56e), and (56f), the original Hodgkin-Huxley model [31] reads

Now, in order to apply the singular perturbation method to the Hodgkin-Huxley model, two small multiplicative parameters 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 () and () enables justifying such an introduction. So, in order to differentiate slow variables from fast variables, Suckley and Biktashev [32], Suckley [33], and Rubin and Wechselberger [26] have plotted the inverse of “time constant for gating variable ,” that is, according to with . In Figure 1, they have been plotted for the original functions and (53a). However, let us notice that this plot is exactly the same as those presented by Rubin and Wechselberger [26] (Figure ) for a nondimensionalized three-dimensional Hodgkin-Huxley singularly perturbed system obtained after the following variable changes: and ; then and finally .

Figure 1: Graph of (ms−1) against ().

Figure 1 shows a plot of the functions according to with over the physiological range. We observe that is of an order of magnitude bigger than and , which are of comparable size. Indeed, we can deduce that the values of times scales are approximately while . Then, it appears that corresponds to the fast variable while and correspond to slow variables. Moreover, since the activation of the sodium channel is directly related to the dynamics of the membrane (action) potential , Rubin and Wechselberger [26] consider that and evolve on the same fast timescale. So, the Hodgkin-Huxley model may be transformed into a singularly perturbed system with two timescales in which the slow variables are and the fast variables are .

So, according to Awiszus et al. [34], Suckley and Biktashev [32], Suckley [33], and Rubin and Wechselberger [26] small multiplicative parameters 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 and set:Other parameters are kept as for the original Hodgkin-Huxley model [31]:

Then, by posing , , and () () to be consistent with the notations of Section 3, we obtainwhere () () and .

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 , that is, , is tenth times greater than the timescale of variables and , that is, of variables and . Moreover, this parameter is identical to those used in (60d) since it has been also considered that and , that is, and , evolve on the same fast timescale.

According to the Geometric Singular Perturbation Theory, the zero-order approximation in of the slow manifold associated with the Hodgkin-Huxley model (60a), (60b), (60c), and (60d) is obtained by posing in (60c) and (60d). So, the slow manifold is given by

Then, the fast foliation is within the planes and .

The fold curve is defined as the location of the points, where , , and . For the Hodgkin-Huxley model (60a), (60b), (60c), and (60d), the fold curve is thus given by (61a) and (61b) and by the determinant of the Jacobian matrix of the following  fast foliation:

The Jacobian matrix of the fast foliation (62a) and (62b) readswhere denotes the derivative with respect to . Then, taking into account (61b), that is, , we have

So, the determinant of the Jacobian matrix of the fast foliation (62a) and (62b) is

Thus, the condition for the fold curve is , which gives

Therefore

By subtracting (61a) from (67) we obtain

Plugging this value of (68) into (67) provides

So, the fold curve is given by the set of parametric equations ((68)-(69)) in terms of .

The pseudo singular points are given by (22a), (22b), (22c), (22d), and (22e) which read for the Hodgkin-Huxley model (60a), (60b), (60c), and (60d):