# Three-Dimensional Memristive Hindmarsh–Rose Neuron Model with Hidden Coexisting Asymmetric Behaviors

Authors: Bocheng Bao; Aihuang Hu; Han Bao; Quan Xu; Mo Chen; Huagan Wu
Journal: Complexity
 (2018)

---

## Abstract

Since the electrical activities of neurons are closely related to complex electrophysiological environment in neuronal system, a novel three-dimensional memristive Hindmarsh–Rose (HR) neuron model is presented in this paper to describe complex dynamics of neuronal activities with electromagnetic induction. The proposed memristive HR neuron model has no equilibrium point but can show hidden dynamical behaviors of coexisting asymmetric attractors, which has not been reported in the previous references for the HR neuron model. Mathematical model based numerical simulations for hidden coexisting asymmetric attractors are performed by bifurcation analyses, phase portraits, attraction basins, and dynamical maps, which just demonstrate the occurrence of complex dynamical behaviors of electrical activities in neuron with electromagnetic induction. Additionally, circuit breadboard based experimental results well confirm the numerical simulations.


---

## Body

## 1. Introduction

In the past three decades, numerous simplified neuron models had been fantastically extended from the classical Hodgkin–Huxley model [1] to reconstruct the main dynamical characteristics of neuronal electrical activities [2–8], among which the two- and three-dimensional Hindmarsh–Rose (HR) neuron models are effective and available for dynamical analysis in electrical activities of biological neurons [9, 10]. In the last few years, a wide variety of the HR neuron models, such as original three-dimensional HR models [10–16], extended or nonlinear feedback coupled HR models [17–20], time delayed HR models [20–22], fractional-order HR models [23, 24], and memristor based HR models under electromagnetic radiations [9, 25–27], have been proposed and further studied by bifurcation analysis methods for understanding the dynamics of electrical activities among neurons [8]. For this reason, bifurcation analysis theory plays an essential role in describing mode transitions between spiking and bursting in the neuronal electrical activities [9–27].Inspired by the constructing approach of the three-dimensional HR neuron model [3, 24], a novel three-dimensional memristive HR neuron model is presented in this paper, which could be used to better describe complex dynamical characteristics of neuronal electrical activities with electromagnetic induction or further exhibit some undiscovered complex dynamical behaviors in neuronal electrical activities. Interestingly, our proposed memristive HR neuron model has no equilibrium point, which can be classified as a particular dynamical system with hidden oscillating patterns [28–31]. Furthermore, coexisting asymmetric attractors’ behavior can also be observed in such memristive HR neuron model as well, indicating the emergence of bistability dynamics, which has been found in some specified neuron models [32–35]. However, the phenomenon of hidden coexisting asymmetric attractors has not been previously reported for the HR neuron model.This paper is organized as follows. In Section2, based on the brief reviews on the HR neuron model, a three-dimensional memristive HR neuron model is presented, upon which hidden coexisting asymmetric attractors are numerically revealed by phase portraits and time series and its bistability dynamics are confirmed by the attraction basins related to the initial values. In Section 3, hidden coexisting asymmetric attractors’ behaviors are demonstrated by bifurcation diagrams, Lyapunov exponents, and dynamical maps, from which numerous types of coexisting asymmetric attractors are easily observed. In addition, a physical implementation circuit is fabricated and breadboard experiments are carried out to confirm the hidden coexisting asymmetric attractors in Section 4. The conclusions are summarized in Section 5.


## 2. Three-Dimensional Memristive HR Neuron Model

### 2.1. Brief Reviews on the HR Neuron Model

Through simplifying the classical Hodgkin–Huxley model [1], the two-dimensional Hindmarsh–Rose (HR) neuron model was proposed by Hindmarsh and Rose [2] in 1982, which is described by two first-order ordinary differential equations as(1)x˙=y-ax3+bx2+Iy˙=c-dx2-y,in which two variables x and y are the membrane potential and recovery variable (also called spiking variable), respectively, and a term I is the externally applied current. The parameters a, b, c, and d are four positive constants, which are often assumed as a=1, b=3, c=1, and d=5, respectively [2, 11–14].To permit numerous dynamical behaviors, for example, chaotic dynamics, for the membrane potential, an extra third equation was introduced by Hindmarsh and Rose [3] in 1984 to improve the two-dimensional neuron model (1), which is expressed by three first-order ordinary differential equations as(2)x˙=y-ax3+bx2+I-zy˙=c-dx2-yz˙=rsx-x1-z,where the variable z is the bursting variable and the constant x1 is the resting potential of the model. The newly added parameters r and s are two positive constants but r is very small. Thus, a new variable z, a slowly evolving current, is coupled into the first equation of the two-dimensional model (1) to tune the externally applied current I. If the three-dimensional neuron model (2) is in its firing state, the value of z increases [24].


### 2.2. Constructed Memristive HR Neuron Model

Motivated by the above constructing approach of the model (2), a three-dimensional memristive HR neuron model with electromagnetic induction effect is proposed through introducing a flux-controlled ideal memristor into the first equation of the model (1), which can be mathematically modeled by(3)x˙=y-ax3+bx2+I+kφxy˙=c-dx2-yφ˙=x,where the new variable φ is the magnetic flux indicating the time integral of the membrane potential x. The newly adding term kφx represents the externally applied electromagnetic induction and k is the strength of the electromagnetic induction.It is important to stress that the introduced memristor in (3) is ideal and flux-controlled. According to the definitions of ideal memristor by state-dependent Ohm’s laws between the terminal voltage v and terminal current i [36, 37], a flux-controlled ideal memristor W is thereby given as(4)i=Wφv,φ˙=v,where the memductance W(φ) can be interpreted as the flux-dependent rate of change of charge. Therefore, the memductance W(φ) utilized in (3) can be written by(5)Wφ=kφ,in which the coefficient k is positive.To exhibit three characteristic fingerprints of pinched hysteresis loop of the memristor modeled by (4) and (5) [38], a sinusoidal voltage source v=Asin⁡(2πFτ) is connected at the input terminals of the memristor, where A and F are the amplitude and frequency, respectively. Let k=1. When A=4 is maintained unchanged and F is assigned as 0.1, 0.2, and 0.5, respectively, the v-i plots are shown in Figure 1(a), while when F=0.1 is fixed and A is determined as 3, 4, and 5, respectively, the v-i plots are shown in Figure 1(b). It is seen from Figure 1 that the v-i plots are the hysteresis loops pinched at the origin. The hysteresis loop is pinched regardless of the stimulus amplitudes but shrinks into a linear function at infinite frequency and its lobe area decreases with increasing the frequency. The numerical results in Figure 1 indicate that the memristor modeled by (4) can behave three fingerprints for distinguishing memristors [38].Figure 1
Pinched hysteresis loops of the flux-controlled ideal memristor. (a)A=4 with F=0.1, 0.2, and 0.5. (b) F=0.1 with A=3, 4, and 5.
(a) (b)In the next work, the three-dimensional memristive HR neuron model given in (3) is considered. It should be remarked that the adjustable parameters of interest are I and k, and their regions are correspondent to the first quadrant of the parameter space (I>0 and k>0). For any uncertain parameter I, the existence of any equilibrium point is not allowed in the three-dimensional memristive HR neuron model, neither stable nor unstable. Only if the applied current I=-1 will the model show an equilibrium point, which is not in the considered parameter region. This case is often encountered in various kinds of nonlinear dynamical systems that are known to generate the specified hidden attractors [28–31].


### 2.3. Coexisting Asymmetric Attractors

When the original parameters are selected asa=1, b=3, c=1, and d=5, respectively [2], an example for model (3) with I=1 and k=0.9 is given as shown in Figure 2, where the orbits marked by the red and blue colors emerge from the initial values (0, 0, −2) and (0, 0, 2), respectively. In Figure 2(a), the phase portraits in the x-φ plane display the bistability phenomenon of hidden coexisting asymmetric attractors consisting of chaotic attractor and limit cycle in the memristive HR neuron model, whereas in Figure 2(b), the time series of the membrane potential x demonstrate the coexistence of chaotic and periodic spikes in the memristive HR neuron model as well. Correspondingly, three Lyapunov exponents for the initial values (0, 0, −2) are 0.0782, 0, and −3.0684, respectively, while those for (0, 0, 2) are 0, −0.2717, and −2.8556, respectively. Remark that Wolf et al.’s method [39] with MATLAB ODE113 algorithm is here used to calculate three Lyapunov exponents.Figure 2
Hidden coexisting asymmetric attractors emerged from the initial values (0, 0, −2) and (0, 0, 2). (a) Phase portraits in thex-φ plane. (b) Time series of the variable x.
(a) (b)For the coexisting asymmetric attractors shown in Figure2(a), the corresponding attraction basins in the x(0)-φ(0) and x(0)-y(0) planes of the initial values are drawn in Figures 3(a) and 3(b), where the attraction basins for chaotic attractors and periodic limit cycles are colored in the fuchsia and cyan regions, respectively. The results effectively indicate the emergence of bistability phenomenon in the memristive HR neuron model.Figure 3
Attraction basins in two different planes forI=1 and k=0.9, indicating the emergence of bistability phenomenon. (a) The x(0)-φ(0) plane with y(0)=0. (b) The x(0)-y(0) plane with φ(0)=0.
(a) (b)Particularly, the emerging coexisting asymmetric attractors do not associate with any equilibrium point, indicating that the memristive HR neuron model always operates in hidden oscillating patterns [28–31]. Additionally, it is interesting to note that, just like the self-excited coexisting asymmetric attractors in hyperbolic-type memristor based Hopfield neural network [32], such hidden coexisting asymmetric attractors in the memristive HR neuron model are induced by electromagnetic induction also, which illustrates the occurrence of complex dynamical behaviors of electrical activities in neuron with electromagnetic induction.


## 2.1. Brief Reviews on the HR Neuron Model

Through simplifying the classical Hodgkin–Huxley model [1], the two-dimensional Hindmarsh–Rose (HR) neuron model was proposed by Hindmarsh and Rose [2] in 1982, which is described by two first-order ordinary differential equations as(1)x˙=y-ax3+bx2+Iy˙=c-dx2-y,in which two variables x and y are the membrane potential and recovery variable (also called spiking variable), respectively, and a term I is the externally applied current. The parameters a, b, c, and d are four positive constants, which are often assumed as a=1, b=3, c=1, and d=5, respectively [2, 11–14].To permit numerous dynamical behaviors, for example, chaotic dynamics, for the membrane potential, an extra third equation was introduced by Hindmarsh and Rose [3] in 1984 to improve the two-dimensional neuron model (1), which is expressed by three first-order ordinary differential equations as(2)x˙=y-ax3+bx2+I-zy˙=c-dx2-yz˙=rsx-x1-z,where the variable z is the bursting variable and the constant x1 is the resting potential of the model. The newly added parameters r and s are two positive constants but r is very small. Thus, a new variable z, a slowly evolving current, is coupled into the first equation of the two-dimensional model (1) to tune the externally applied current I. If the three-dimensional neuron model (2) is in its firing state, the value of z increases [24].


## 2.2. Constructed Memristive HR Neuron Model

Motivated by the above constructing approach of the model (2), a three-dimensional memristive HR neuron model with electromagnetic induction effect is proposed through introducing a flux-controlled ideal memristor into the first equation of the model (1), which can be mathematically modeled by(3)x˙=y-ax3+bx2+I+kφxy˙=c-dx2-yφ˙=x,where the new variable φ is the magnetic flux indicating the time integral of the membrane potential x. The newly adding term kφx represents the externally applied electromagnetic induction and k is the strength of the electromagnetic induction.It is important to stress that the introduced memristor in (3) is ideal and flux-controlled. According to the definitions of ideal memristor by state-dependent Ohm’s laws between the terminal voltage v and terminal current i [36, 37], a flux-controlled ideal memristor W is thereby given as(4)i=Wφv,φ˙=v,where the memductance W(φ) can be interpreted as the flux-dependent rate of change of charge. Therefore, the memductance W(φ) utilized in (3) can be written by(5)Wφ=kφ,in which the coefficient k is positive.To exhibit three characteristic fingerprints of pinched hysteresis loop of the memristor modeled by (4) and (5) [38], a sinusoidal voltage source v=Asin⁡(2πFτ) is connected at the input terminals of the memristor, where A and F are the amplitude and frequency, respectively. Let k=1. When A=4 is maintained unchanged and F is assigned as 0.1, 0.2, and 0.5, respectively, the v-i plots are shown in Figure 1(a), while when F=0.1 is fixed and A is determined as 3, 4, and 5, respectively, the v-i plots are shown in Figure 1(b). It is seen from Figure 1 that the v-i plots are the hysteresis loops pinched at the origin. The hysteresis loop is pinched regardless of the stimulus amplitudes but shrinks into a linear function at infinite frequency and its lobe area decreases with increasing the frequency. The numerical results in Figure 1 indicate that the memristor modeled by (4) can behave three fingerprints for distinguishing memristors [38].Figure 1
Pinched hysteresis loops of the flux-controlled ideal memristor. (a)A=4 with F=0.1, 0.2, and 0.5. (b) F=0.1 with A=3, 4, and 5.
(a) (b)In the next work, the three-dimensional memristive HR neuron model given in (3) is considered. It should be remarked that the adjustable parameters of interest are I and k, and their regions are correspondent to the first quadrant of the parameter space (I>0 and k>0). For any uncertain parameter I, the existence of any equilibrium point is not allowed in the three-dimensional memristive HR neuron model, neither stable nor unstable. Only if the applied current I=-1 will the model show an equilibrium point, which is not in the considered parameter region. This case is often encountered in various kinds of nonlinear dynamical systems that are known to generate the specified hidden attractors [28–31].


## 2.3. Coexisting Asymmetric Attractors

When the original parameters are selected asa=1, b=3, c=1, and d=5, respectively [2], an example for model (3) with I=1 and k=0.9 is given as shown in Figure 2, where the orbits marked by the red and blue colors emerge from the initial values (0, 0, −2) and (0, 0, 2), respectively. In Figure 2(a), the phase portraits in the x-φ plane display the bistability phenomenon of hidden coexisting asymmetric attractors consisting of chaotic attractor and limit cycle in the memristive HR neuron model, whereas in Figure 2(b), the time series of the membrane potential x demonstrate the coexistence of chaotic and periodic spikes in the memristive HR neuron model as well. Correspondingly, three Lyapunov exponents for the initial values (0, 0, −2) are 0.0782, 0, and −3.0684, respectively, while those for (0, 0, 2) are 0, −0.2717, and −2.8556, respectively. Remark that Wolf et al.’s method [39] with MATLAB ODE113 algorithm is here used to calculate three Lyapunov exponents.Figure 2
Hidden coexisting asymmetric attractors emerged from the initial values (0, 0, −2) and (0, 0, 2). (a) Phase portraits in thex-φ plane. (b) Time series of the variable x.
(a) (b)For the coexisting asymmetric attractors shown in Figure2(a), the corresponding attraction basins in the x(0)-φ(0) and x(0)-y(0) planes of the initial values are drawn in Figures 3(a) and 3(b), where the attraction basins for chaotic attractors and periodic limit cycles are colored in the fuchsia and cyan regions, respectively. The results effectively indicate the emergence of bistability phenomenon in the memristive HR neuron model.Figure 3
Attraction basins in two different planes forI=1 and k=0.9, indicating the emergence of bistability phenomenon. (a) The x(0)-φ(0) plane with y(0)=0. (b) The x(0)-y(0) plane with φ(0)=0.
(a) (b)Particularly, the emerging coexisting asymmetric attractors do not associate with any equilibrium point, indicating that the memristive HR neuron model always operates in hidden oscillating patterns [28–31]. Additionally, it is interesting to note that, just like the self-excited coexisting asymmetric attractors in hyperbolic-type memristor based Hopfield neural network [32], such hidden coexisting asymmetric attractors in the memristive HR neuron model are induced by electromagnetic induction also, which illustrates the occurrence of complex dynamical behaviors of electrical activities in neuron with electromagnetic induction.


## 3. Hidden Coexisting Asymmetric Attractors’ Behavior

When the applied currentI and electromagnetic induction strength k are considered as two bifurcation parameters, hidden coexisting asymmetric behaviors of the memristive HR neuron model are numerically studied by MATLAB ODE45 algorithm under two sets of the initial values (0, 0, −2) and (0, 0, 2).


### 3.1. Bifurcation Behaviors with Increasingk

Figure4 gives the bifurcation diagrams of x and the first two Lyapunov exponents as I = 1 and k = 0.5~1.4, where in Figure 4(a) the orbits marked by the red and blue colors emerge from the initial values (0, 0, −2) and (0, 0, 2), respectively, and in Figure 4(b) the Lyapunov exponents marked by the red and fuchsia colors correspond to the initial values (0, 0, −2) and those marked by blue and green colors correspond to (0, 0, 2). It can be seen from Figure 4 that complex dynamics appear in the memristive HR neuron model, in which chaotic attractors with different topologies, limit cycles with different periodicities, period-doubling bifurcation routes, tangent bifurcation routes, crisis scenarios, coexisting bifurcation modes, and so on can be found. Therefore, the electromagnetic induction by the introduced memristor induces numerous complex dynamics for the membrane potential x, especially including hidden coexisting asymmetric behaviors.Figure 4
For the initial values (0, 0, −2) and (0, 0, 2), hidden coexisting asymmetric behaviors withI = 1 and k = 0.5~1.4. (a) Bifurcation diagrams of x. (b) First two Lyapunov exponents.
(a) (b)The concernedly coexisting behaviors of asymmetric attractors mainly locate in two parameter regions[0.718,0.736] and [0.754,0.909], in which some different types of hidden coexisting asymmetric attractors occur. When I=1 and k=0.735 and 0.81, respectively, the phase portraits in the x-φ plane for the other two types of hidden coexisting asymmetric attractors are depicted in Figure 5. In detail, Figure 5(a) exhibits the coexistence of hidden chaotic attractor and hidden limit cycle, and Figure 5(b) displays the coexistence of two hidden limit cycles with different periodicities.Figure 5
Phase portraits in thex-φ plane for the other two types of hidden coexisting asymmetric attractors. (a) I = 1 and k = 0.735. (b) I = 1 and k = 0.81.
(a) (b)


### 3.2. Bifurcation Behaviors with IncreasingI

Figure6 demonstrates the bifurcation diagrams of x and the first two Lyapunov exponents as I = 0~2.4 and k = 0.9, where the initial values for different colored orbits in Figure 6(a) and different colored Lyapunov exponents in Figure 6(b) are consistent with those used in Figures 4(a) and 4(b). In the same way, it can be observed from Figure 6 that complex dynamics are coined in the memristive HR neuron model, reflecting the dynamical effect of the externally applied current I in neurons.Figure 6
For the initial values (0, 0, −2) and (0, 0, 2), hidden coexisting asymmetric behaviors withI = 0~2.4 and k = 0.9. (a) Bifurcation diagrams of x. (b) First two Lyapunov exponents.
(a) (b)The parameter region[0.98,1.65] has the benefit for the coexisting behaviors of asymmetric attractors, in which several different types of hidden coexisting asymmetric attractors can be clearly found. When k = 0.9 is fixed and I is set to 1.15 and to 1.62, respectively, the phase portraits in the x-φ plane for the two types of hidden coexisting asymmetric attractors are plotted in Figure 7, where in Figure 7(a) the coexistence of hidden chaotic attractor and hidden period 1 limit cycle is displayed, and in Figure 7(b) the coexistence of hidden period 2 limit cycle and hidden chaotic attractor with large size is presented.Figure 7
Phase portraits in thex-φ plane for the other two types of hidden coexisting asymmetric attractors. (a) I = 1.15 and k = 0.9. (b) I = 1.62 and k = 0.9.
(a) (b)


### 3.3. Coexisting Asymmetric Behaviors in the Parameter Space

For intuitively manifesting the coexisting behaviors of asymmetric attractors in the memristive HR neuron model, dynamical maps depicted by the largest Lyapunov exponent under two sets of the initial values are numerically plotted in thek-I parameter space [40], as shown in Figures 8(a) and 8(b), where the luminous yellow, red, and black colored regions stand for the chaotic, periodic, and divergent behaviors, respectively. Figure 8 indicates how dynamical evolution in the electromagnetic induction strength k and applied current I affects the coexisting behaviors under different initial values being considered. When the two parameters k and I are evolved, some chaotic regions are embedded in the periodic regions; however, different chaotic regions appear on the parameter spaces of Figures 8(a) and 8(b), which are caused by the coexisting asymmetric attractors’ behaviors under different initial values. The numerical results in Figure 8 illustrate that the dynamical behaviors depicted by the largest Lyapunov exponent based dynamical maps are well agreed with those revealed by the bifurcation behaviors in Figures 4 and 6.Figure 8
In thek-I parameter space, dynamical maps depicted by the largest Lyapunov exponent under different initial values. (a) Dynamical map for the initial values (0, 0, −2). (b) Dynamical map for the initial values (0, 0, 2).
(a) (b)It should be remarkable for the dynamical maps in Figure8 that the chaotic attractors in different locations of the chaotic region have different topologies and the limit cycles in different locations of the periodic region have different periodicities. Specifically, except for several types of coexisting asymmetric behaviors shown in Figures 2, 5, and 7, another type of coexisting asymmetric behaviors of chaotic attractor and divergent orbit can also be uncovered, which means that another form of bistability exists in the memristive HR neuron model.


## 3.1. Bifurcation Behaviors with Increasingk

Figure4 gives the bifurcation diagrams of x and the first two Lyapunov exponents as I = 1 and k = 0.5~1.4, where in Figure 4(a) the orbits marked by the red and blue colors emerge from the initial values (0, 0, −2) and (0, 0, 2), respectively, and in Figure 4(b) the Lyapunov exponents marked by the red and fuchsia colors correspond to the initial values (0, 0, −2) and those marked by blue and green colors correspond to (0, 0, 2). It can be seen from Figure 4 that complex dynamics appear in the memristive HR neuron model, in which chaotic attractors with different topologies, limit cycles with different periodicities, period-doubling bifurcation routes, tangent bifurcation routes, crisis scenarios, coexisting bifurcation modes, and so on can be found. Therefore, the electromagnetic induction by the introduced memristor induces numerous complex dynamics for the membrane potential x, especially including hidden coexisting asymmetric behaviors.Figure 4
For the initial values (0, 0, −2) and (0, 0, 2), hidden coexisting asymmetric behaviors withI = 1 and k = 0.5~1.4. (a) Bifurcation diagrams of x. (b) First two Lyapunov exponents.
(a) (b)The concernedly coexisting behaviors of asymmetric attractors mainly locate in two parameter regions[0.718,0.736] and [0.754,0.909], in which some different types of hidden coexisting asymmetric attractors occur. When I=1 and k=0.735 and 0.81, respectively, the phase portraits in the x-φ plane for the other two types of hidden coexisting asymmetric attractors are depicted in Figure 5. In detail, Figure 5(a) exhibits the coexistence of hidden chaotic attractor and hidden limit cycle, and Figure 5(b) displays the coexistence of two hidden limit cycles with different periodicities.Figure 5
Phase portraits in thex-φ plane for the other two types of hidden coexisting asymmetric attractors. (a) I = 1 and k = 0.735. (b) I = 1 and k = 0.81.
(a) (b)


## 3.2. Bifurcation Behaviors with IncreasingI

Figure6 demonstrates the bifurcation diagrams of x and the first two Lyapunov exponents as I = 0~2.4 and k = 0.9, where the initial values for different colored orbits in Figure 6(a) and different colored Lyapunov exponents in Figure 6(b) are consistent with those used in Figures 4(a) and 4(b). In the same way, it can be observed from Figure 6 that complex dynamics are coined in the memristive HR neuron model, reflecting the dynamical effect of the externally applied current I in neurons.Figure 6
For the initial values (0, 0, −2) and (0, 0, 2), hidden coexisting asymmetric behaviors withI = 0~2.4 and k = 0.9. (a) Bifurcation diagrams of x. (b) First two Lyapunov exponents.
(a) (b)The parameter region[0.98,1.65] has the benefit for the coexisting behaviors of asymmetric attractors, in which several different types of hidden coexisting asymmetric attractors can be clearly found. When k = 0.9 is fixed and I is set to 1.15 and to 1.62, respectively, the phase portraits in the x-φ plane for the two types of hidden coexisting asymmetric attractors are plotted in Figure 7, where in Figure 7(a) the coexistence of hidden chaotic attractor and hidden period 1 limit cycle is displayed, and in Figure 7(b) the coexistence of hidden period 2 limit cycle and hidden chaotic attractor with large size is presented.Figure 7
Phase portraits in thex-φ plane for the other two types of hidden coexisting asymmetric attractors. (a) I = 1.15 and k = 0.9. (b) I = 1.62 and k = 0.9.
(a) (b)


## 3.3. Coexisting Asymmetric Behaviors in the Parameter Space

For intuitively manifesting the coexisting behaviors of asymmetric attractors in the memristive HR neuron model, dynamical maps depicted by the largest Lyapunov exponent under two sets of the initial values are numerically plotted in thek-I parameter space [40], as shown in Figures 8(a) and 8(b), where the luminous yellow, red, and black colored regions stand for the chaotic, periodic, and divergent behaviors, respectively. Figure 8 indicates how dynamical evolution in the electromagnetic induction strength k and applied current I affects the coexisting behaviors under different initial values being considered. When the two parameters k and I are evolved, some chaotic regions are embedded in the periodic regions; however, different chaotic regions appear on the parameter spaces of Figures 8(a) and 8(b), which are caused by the coexisting asymmetric attractors’ behaviors under different initial values. The numerical results in Figure 8 illustrate that the dynamical behaviors depicted by the largest Lyapunov exponent based dynamical maps are well agreed with those revealed by the bifurcation behaviors in Figures 4 and 6.Figure 8
In thek-I parameter space, dynamical maps depicted by the largest Lyapunov exponent under different initial values. (a) Dynamical map for the initial values (0, 0, −2). (b) Dynamical map for the initial values (0, 0, 2).
(a) (b)It should be remarkable for the dynamical maps in Figure8 that the chaotic attractors in different locations of the chaotic region have different topologies and the limit cycles in different locations of the periodic region have different periodicities. Specifically, except for several types of coexisting asymmetric behaviors shown in Figures 2, 5, and 7, another type of coexisting asymmetric behaviors of chaotic attractor and divergent orbit can also be uncovered, which means that another form of bistability exists in the memristive HR neuron model.


## 4. Circuit Design and Breadboard Experiments

### 4.1. Physical Circuit Designs and Parameter Selections

The flux-controlled ideal memristorW characterized by (4) and (5) and its constructing three-dimensional memristive HR neuron model expressed by (3) can be physically realized by using an electronic circuit via analog multipliers and operational amplifiers connected with resistors and/or capacitors [41–43], as drawn in Figures 9(a) and 9(b), respectively. Of course, this three-dimensional memristive HR neuron model also can be digitally implemented in field-programmable gate arrays (FPGA) as well [44, 45].Figure 9
Physical electronic circuit implementation for the memristive HR neuron model. (a) Implementation circuit of the flux-controlled ideal memristor. (b) Main circuit of the memristive HR neuron model.
(a) (b)The implementation circuit of the flux-controlled ideal memristor in Figure9(a) contains an integrator with time constant RC, an inverter, a multiplier M0, and a resistor Rk. For the input voltage v and output current i, the mathematical model for the memristor emulator can be easily given as(6)i=Wvφv=g0RkCvφv=1RC·kvφvv˙φ=1RC·v,where vφ is the inner variable of the memristor emulator and g0 is the gain of the multipliers M0, k = g0R/Rk, and W(vφ)=kvφ/RC.The main circuit of the memristive HR neuron model in Figure9(b) has two integrating channels for implementing the first and second equations of (3). According to Kirchhoff’s circuit laws and electrical properties of the circuit components, the circuit equations of Figure 9(b) are written as(7)RCdvxdt=RR2vy-g1g2RR4vx3+g2RR3vx2+RR1VI+g0RRkvφvxRCdvydt=RR6Vy0-g3RR5vx2-RR7vyRCdvφdt=vx,where vx and vy are two circuit variables, VI and Vy0 are two applied voltages, and g1, g2, and g3 are the gains of the multipliers M1, M2, and M3, respectively.Considering that the dynamic amplitude of the recovery variabley in the numerical simulations exceeds the linear operation ranges of operational amplifier and multiplier, the following linear transformation(8)vx,vy,vφ⟶vx,2.5vy,vφshould be utilized to reduce the dynamic voltage amplitude of vy in the circuit equations of the memristive HR neuron model. Thus, by comparing (8) with (3), there yields(9)R1=RI,Rk=g0Rk,R2=R2.5,R3=g2R3,R4=g1g2R,R5=g3R2,R6=2.5R,R7=R.Let the time constant RC = 10 kΩ  ×  33 nF = 330 μs; that is, R = 10 kΩ and C = 33 nF and the multiplier gains g0 = g1 = 0.1 and g2 = g3 = 1. According to (9), the circuit parameters for the breadboard experiments of the memristive HR neuron model can be calculated, as listed in Table 1.Table 1
Circuit parameters of the memristive HR neuron model for breadboard experiments.
Parameters Significations Values R, R7 Resistance 10 kΩ R 1 Resistance 10 kΩ (adjustable) R k Resistance 1.11 kΩ (adjustable) R 2 Resistance 4 kΩ R 3 Resistance 3.3 kΩ R 4 Resistance 1 kΩ R 5 Resistance 5 kΩ R 6 Resistance 25 kΩ C Capacitance 33 nF V I, Vy0 DC voltage 1 V g 0, g1 Multiplier gain 0.1 g 2, g3 Multiplier gain 1


### 4.2. Results Captured from Breadboard Experiments

According to the circuit diagrams in Figure9 and circuit parameters in Table 1, a hardware circuit using commercially discrete components can be welded on a breadboard. The operational amplifiers AD711JN and analog multipliers AD633JN supplied by ±15 V voltage modules are chosen. The DC voltages VI and Vy0 are provided by Tektronix PWS 2326 DC Power Supply and the experimental results are measured by Tektronix TDS 3054C Digital Phosphor Oscilloscope. The photograph of the connectedly experimental prototype for the memristive HR neuron model is displayed in Figure 10.Figure 10
Photograph of the experimental breadboard and typical chaotic attractor captured by the digital oscilloscope.For experimentally measuring the pinched hysteresis loops of the memristor emulator given in Figure9(a), a sinusoidal voltage source v=Asin⁡(2πft) generated by Tektronix AFG 3102C Function Generator is linked to the input terminals of the memristor emulator, where the physical frequency is calculated by f=F/RC. For k = 1, the adjustable resistance Rk = 10 kΩ. When the amplitudes and frequencies of the sinusoidal voltage source used during numerical simulations in Figure 1 are employed, the pinched hysteresis loops for the corresponding amplitudes and frequencies are captured, as shown in Figures 11(a) and 11(b), respectively, which experimentally validate the characteristic fingerprints of the memristor emulator. It should be addressed that for better observing the experimental results, all the output currents sensed by the current probe are magnified by ten times.Figure 11
Experimentally measured pinched hysteresis loops of the memristor emulator. (a)A = 4 V with different frequencies. (b) f = 303.03 Hz with different amplitudes.
(a) (b)The circuit parameters listed in Table1 are used and the different initial voltages of three capacitors are randomly sensed by repeatedly switching on and off the experimental power supply [46]. For the typical circuit parameters in Table 1, two adjustable circuit parameters of R1 and Rk correspond to the adjustable model parameters of I = 1 and k = 0.9. Corresponding to Figure 2, the phase portraits in the vx-vφ plane and time series of the variable vx that emerged from different initial voltages are experimentally obtained, as shown in Figure 12. The experimental results indicate that hidden coexisting asymmetric attractors also can be measured from the breadboard experiments of the memristive HR neuron model as well.Figure 12
Experimentally measured hidden coexisting asymmetric attractors while repeatedly switching on and off the experimental power supply. (a) Phase portraits in thevx-vφ plane. (b) Time series of the variable vx.
(a) (b)When the applied currentI = 1, that is, the resistance R1 is fixed as 10 kΩ and the resistance Rk is set to 1.36 kΩ and to 1.23 kΩ, respectively, the phase portraits in the vx-vφ plane are captured, as shown in Figures 13(a) and 13(b). Furthermore, when the electromagnetic induction strength k = 0.9, that is, Rk = 1.11 kΩ and R1 is set to 8.70 kΩ and to 6.17 kΩ, respectively, the phase portraits in the vx-vφ plane are captured, as shown in Figures 13(c) and 13(d). Ignoring some tiny differences between numerical simulations and breadboard experiments due to the computational errors and parasitic circuit parameters, the experimental results are almost the same as the numerical simulations, which imply that the coexisting asymmetric attractors’ behaviors that emerged from the memristive HR neuron model can be validated experimentally.Figure 13
Experimentally measured phase portraits in thevx-vφ plane for hidden coexisting asymmetric attractors under different circuit parameters. (a) R1 = 10 kΩ and Rk = 1.36 kΩ. (b) R1 = 10 kΩ and Rk = 1.23 kΩ. (c) R1 = 8.70 kΩ and Rk = 1.11 kΩ. (d) R1 = 6.17 kΩ and Rk = 1.11 kΩ.
(a) (b) (c) (d)


## 4.1. Physical Circuit Designs and Parameter Selections

The flux-controlled ideal memristorW characterized by (4) and (5) and its constructing three-dimensional memristive HR neuron model expressed by (3) can be physically realized by using an electronic circuit via analog multipliers and operational amplifiers connected with resistors and/or capacitors [41–43], as drawn in Figures 9(a) and 9(b), respectively. Of course, this three-dimensional memristive HR neuron model also can be digitally implemented in field-programmable gate arrays (FPGA) as well [44, 45].Figure 9
Physical electronic circuit implementation for the memristive HR neuron model. (a) Implementation circuit of the flux-controlled ideal memristor. (b) Main circuit of the memristive HR neuron model.
(a) (b)The implementation circuit of the flux-controlled ideal memristor in Figure9(a) contains an integrator with time constant RC, an inverter, a multiplier M0, and a resistor Rk. For the input voltage v and output current i, the mathematical model for the memristor emulator can be easily given as(6)i=Wvφv=g0RkCvφv=1RC·kvφvv˙φ=1RC·v,where vφ is the inner variable of the memristor emulator and g0 is the gain of the multipliers M0, k = g0R/Rk, and W(vφ)=kvφ/RC.The main circuit of the memristive HR neuron model in Figure9(b) has two integrating channels for implementing the first and second equations of (3). According to Kirchhoff’s circuit laws and electrical properties of the circuit components, the circuit equations of Figure 9(b) are written as(7)RCdvxdt=RR2vy-g1g2RR4vx3+g2RR3vx2+RR1VI+g0RRkvφvxRCdvydt=RR6Vy0-g3RR5vx2-RR7vyRCdvφdt=vx,where vx and vy are two circuit variables, VI and Vy0 are two applied voltages, and g1, g2, and g3 are the gains of the multipliers M1, M2, and M3, respectively.Considering that the dynamic amplitude of the recovery variabley in the numerical simulations exceeds the linear operation ranges of operational amplifier and multiplier, the following linear transformation(8)vx,vy,vφ⟶vx,2.5vy,vφshould be utilized to reduce the dynamic voltage amplitude of vy in the circuit equations of the memristive HR neuron model. Thus, by comparing (8) with (3), there yields(9)R1=RI,Rk=g0Rk,R2=R2.5,R3=g2R3,R4=g1g2R,R5=g3R2,R6=2.5R,R7=R.Let the time constant RC = 10 kΩ  ×  33 nF = 330 μs; that is, R = 10 kΩ and C = 33 nF and the multiplier gains g0 = g1 = 0.1 and g2 = g3 = 1. According to (9), the circuit parameters for the breadboard experiments of the memristive HR neuron model can be calculated, as listed in Table 1.Table 1
Circuit parameters of the memristive HR neuron model for breadboard experiments.
Parameters Significations Values R, R7 Resistance 10 kΩ R 1 Resistance 10 kΩ (adjustable) R k Resistance 1.11 kΩ (adjustable) R 2 Resistance 4 kΩ R 3 Resistance 3.3 kΩ R 4 Resistance 1 kΩ R 5 Resistance 5 kΩ R 6 Resistance 25 kΩ C Capacitance 33 nF V I, Vy0 DC voltage 1 V g 0, g1 Multiplier gain 0.1 g 2, g3 Multiplier gain 1


## 4.2. Results Captured from Breadboard Experiments

According to the circuit diagrams in Figure9 and circuit parameters in Table 1, a hardware circuit using commercially discrete components can be welded on a breadboard. The operational amplifiers AD711JN and analog multipliers AD633JN supplied by ±15 V voltage modules are chosen. The DC voltages VI and Vy0 are provided by Tektronix PWS 2326 DC Power Supply and the experimental results are measured by Tektronix TDS 3054C Digital Phosphor Oscilloscope. The photograph of the connectedly experimental prototype for the memristive HR neuron model is displayed in Figure 10.Figure 10
Photograph of the experimental breadboard and typical chaotic attractor captured by the digital oscilloscope.For experimentally measuring the pinched hysteresis loops of the memristor emulator given in Figure9(a), a sinusoidal voltage source v=Asin⁡(2πft) generated by Tektronix AFG 3102C Function Generator is linked to the input terminals of the memristor emulator, where the physical frequency is calculated by f=F/RC. For k = 1, the adjustable resistance Rk = 10 kΩ. When the amplitudes and frequencies of the sinusoidal voltage source used during numerical simulations in Figure 1 are employed, the pinched hysteresis loops for the corresponding amplitudes and frequencies are captured, as shown in Figures 11(a) and 11(b), respectively, which experimentally validate the characteristic fingerprints of the memristor emulator. It should be addressed that for better observing the experimental results, all the output currents sensed by the current probe are magnified by ten times.Figure 11
Experimentally measured pinched hysteresis loops of the memristor emulator. (a)A = 4 V with different frequencies. (b) f = 303.03 Hz with different amplitudes.
(a) (b)The circuit parameters listed in Table1 are used and the different initial voltages of three capacitors are randomly sensed by repeatedly switching on and off the experimental power supply [46]. For the typical circuit parameters in Table 1, two adjustable circuit parameters of R1 and Rk correspond to the adjustable model parameters of I = 1 and k = 0.9. Corresponding to Figure 2, the phase portraits in the vx-vφ plane and time series of the variable vx that emerged from different initial voltages are experimentally obtained, as shown in Figure 12. The experimental results indicate that hidden coexisting asymmetric attractors also can be measured from the breadboard experiments of the memristive HR neuron model as well.Figure 12
Experimentally measured hidden coexisting asymmetric attractors while repeatedly switching on and off the experimental power supply. (a) Phase portraits in thevx-vφ plane. (b) Time series of the variable vx.
(a) (b)When the applied currentI = 1, that is, the resistance R1 is fixed as 10 kΩ and the resistance Rk is set to 1.36 kΩ and to 1.23 kΩ, respectively, the phase portraits in the vx-vφ plane are captured, as shown in Figures 13(a) and 13(b). Furthermore, when the electromagnetic induction strength k = 0.9, that is, Rk = 1.11 kΩ and R1 is set to 8.70 kΩ and to 6.17 kΩ, respectively, the phase portraits in the vx-vφ plane are captured, as shown in Figures 13(c) and 13(d). Ignoring some tiny differences between numerical simulations and breadboard experiments due to the computational errors and parasitic circuit parameters, the experimental results are almost the same as the numerical simulations, which imply that the coexisting asymmetric attractors’ behaviors that emerged from the memristive HR neuron model can be validated experimentally.Figure 13
Experimentally measured phase portraits in thevx-vφ plane for hidden coexisting asymmetric attractors under different circuit parameters. (a) R1 = 10 kΩ and Rk = 1.36 kΩ. (b) R1 = 10 kΩ and Rk = 1.23 kΩ. (c) R1 = 8.70 kΩ and Rk = 1.11 kΩ. (d) R1 = 6.17 kΩ and Rk = 1.11 kΩ.
(a) (b) (c) (d)


## 5. Conclusions

This paper presents a novel three-dimensional memristive HR neuron model to describe complex dynamics of neuronal activities with electromagnetic induction. The most prominent feature of this neuron model is that it does not contain any equilibrium point but can exhibit hidden coexisting behaviors of asymmetric attractors. Through executing bifurcation analyses, phase portraits, attraction basins, and dynamical maps, hidden coexisting asymmetric attractors are uncovered from the mathematical model and verified from the corresponding breadboard experiments. Thus, the proposed memristive HR neuron model can imitate the complex dynamical behaviors of electrical activities in neuron with electromagnetic induction. Further investigations will be performed in our future works.



---

*Source: 3872573-2018-02-28.xml*