Multistability and Phase Synchronization of Rulkov Neurons Coupled with a Locally Active Discrete Memristor

: In order to enrich the dynamic behaviors of discrete neuron models and more effectively mimic biological neural networks, this paper proposes a bistable locally active discrete memristor (LADM) model to mimic synapses. We explored the dynamic behaviors of neural networks by introducing the LADM into two identical Rulkov neurons. Based on numerical simulation, the neural network manifested multistability and new ﬁring behaviors under different system parameters and initial values. In addition, the phase synchronization between the neurons was explored. Additionally, it is worth mentioning that the Rulkov neurons showed synchronization transition behavior; that is, anti-phase synchronization changed to in-phase synchronization with the change in the coupling strength. In particular, the anti-phase synchronization of different ﬁring patterns in the neural network was investigated. This can characterize the different ﬁring behaviors of coupled homogeneous neurons in the different functional areas of the brain, which is helpful to understand the formation of functional areas. This paper has a potential research value and lays the foundation for biological neuron experiments and neuron-based engineering applications.


Introduction
The brain is regarded as a complex neural network composed of many neurons.Information transmission [1][2][3] and conscious perception [4] are closely related to the dynamical behaviors of biological neural networks.At the same time, synchronization plays a role in some pathologies in neuroscience, such as epilepsy.Experiments [5,6] showed that the nervous system exhibited high synchronization during seizures.In order to better understand the structure of the brain, researchers have studied a large number of dynamic characteristics of neural networks [7][8][9][10][11][12].However, continuous neuron models have the characteristics of high dimension and multiscale, which are not conducive to the study of complex and large-scale neural networks.In recent years, discrete neuron models have drawn more attention due to their simple form and flexible calculation.Extensive studies [13][14][15][16][17][18] proved that discrete neuron models, especially the Rulkov neuron, could more effectively mimic the real behaviors of biological neurons.On this foundation, Cao et al. [19][20][21][22][23] disclosed two identical Rulkov neurons coupled with chemical and electrical synapses, respectively, and reported the influence of the coupling strength on synchronization in detail.By constructing a small-world network made up of Rulkov neurons, Ferrari et al. [24] proposed the idea of using delayed feedback to suppress synchronization.Rakshit S et al. [25] revealed several types of firing patterns and synchronous behavior of Rulkov neurons under the joint effect of internal connections and a chemical synapse.
In 1971, Chua first put forward the concept of a memristor as a device to describe the relationship between the charge and the flux [26].Since HP Laboratory developed the entity of memristor in 2008 [27], continuous memristors have been widely used in neural networks [28][29][30][31][32][33][34] and neural morphological circuits [35][36][37][38][39]. Ding et al. [40] studied the hidden coexisting firing patterns of two heterogeneous fractional-order HR neurons coupled with a memristor and applied them to image encryption.However, the synchronization between the neurons was not analyzed.Because of the properties of memorability, nonvolatility, nanoscale, and local activity, memristors are regarded as the best choice to mimic synapses, which have been proved in continuous neuron models.Additionally, phenomena such as coexisting firing patterns [41][42][43][44] and synchronization [45][46][47] have been found.It was found in [48] that anti-phase synchronization helps to distinguish the different functional areas of the brain.However, the anti-phase synchronization of different firing patterns was not found.Recently, researchers have shifted their focus to discrete memristors.Introducing discrete memristors to different maps, it was found that discrete memristive chaotic systems have the advantages of high complexity [49][50][51], coexisting attractors [52][53][54], and offset boosting [55].Furthermore, Peng et al. [56] also studied parameter identification for discrete memristive chaotic maps.In addition, it was found in [57] that the integer-order discrete memristor has an ideal memory effect, whereas the fractional-order discrete memristor does not.
As far as we know, there is very little discussion about discrete memristors in the field of neural networks, and currently, it basically stays at the single-neuron coupling stage.For example, it was revealed in [58][59][60] that the introduction of discrete memristors could enrich the firing patterns of neurons, and phenomena such as hyperchaotic firing and regime transition were surveyed.Lu et al. [61] studied the dynamic characteristics of Rulkov neurons under electromagnetic induction by introducing a discrete memristor and found that the dynamic behaviors of fractional-order systems are more complex than that of integer-order systems.However, compared with integer-order systems, the calculation cost of fractional-order systems is higher.Xu et al. [62] discovered coexisting firing behaviors by introducing a memristor to a single continuous nonautonomous Rulkov neuron.However, the firing pattern of this system was of a relatively single type.Li et al. [63] reported hidden attractors and coexisting firing patterns in a discrete memristive synapse-coupled single HR neuron model.It can be seen that it is very meaningful to introduce discrete memristors into discrete neuron models.
Inspired by the above literature, we investigate the firing behaviors and phase synchronization of two Rulkov neurons coupled with the LADM in this paper.The main contributions are summarized as follows: (1) Compared with the original neurons, the Rulkov neurons coupled with the LADM produce many novel firing patterns.Meanwhile, we show the coexisting firing behaviors of the neural network; (2) The synchronization transition behavior is explored; (3) The neural network exhibits the anti-phase synchronization of different firing patterns of homogeneous neurons, which has not been found in the previous literature.
The remaining structure of this paper is as follows: Section 2 shows the new bistable LADM model, and its characteristics of nonvolatility, bistability, and local activity are verified through numerical simulation.In Section 3, the Rulkov neuron model coupled with LADM is established.Then, we discuss the influence of the LADM on the firing behavior of discrete neurons in Section 4. The synchronization behaviors between the neurons are reported in Section 5, and finally, the conclusions are given in Section 6.
A sinusoidal voltage signal v(n) = Asin(wT(n)) is added, and the volt-ampere characteristic curves of the LADM are shown in Figure 1.Changing the amplitude, it can be seen from Figure 1a that the pinched hysteresis loops through the origin present an oblique "8" shape, and the area of the loops gradually becomes larger as the amplitude increases.The corresponding pinched hysteresis loops can be observed in Figure 1b as the frequency varies from 0 to 1.When the frequency is small, the pinched hysteresis loop is located in all quadrants.When the frequency increases, the pinched hysteresis loop area decreases and is gradually limited to the second and fourth quadrants.Finally, it is a negative single-valued function, namely a negative resistance.Thus, the LADM meets the definition of a memristor [64].
A sinusoidal voltage signal v(n) = Asin(wT(n)) is added, and the volt-ampere characteristic curves of the LADM are shown in Figure 1.Changing the amplitude, it can be seen from Figure 1a that the pinched hysteresis loops through the origin present an oblique "8" shape, and the area of the loops gradually becomes larger as the amplitude increases.The corresponding pinched hysteresis loops can be observed in Figure 1b as the frequency varies from 0 to 1.When the frequency is small, the pinched hysteresis loop is located in all quadrants.When the frequency increases, the pinched hysteresis loop area decreases and is gradually limited to the second and fourth quadrants.Finally, it is a negative singlevalued function, namely a negative resistance.Thus, the LADM meets the definition of a memristor [64].A multistable memristor refers to a memristor with multiple coexisting stable pinched hysteresis loops under different initial values [65].As shown in Figure 2, with frequency w = 0.2 Hz and amplitude A = 1 V, when the two internal initial states are set as ϕ(0) = 0 and ϕ(0) = 2, respectively, there are two different pinched hysteresis loops on the V-I phase plane, which means that the LADM is bistable.A multistable memristor refers to a memristor with multiple coexisting stable pinched hysteresis loops under different initial values [65].As shown in Figure 2, with frequency w = 0.2 Hz and amplitude A = 1 V, when the two internal initial states are set as φ(0) = 0 and φ(0) = 2, respectively, there are two different pinched hysteresis loops on the V-I phase plane, which means that the LADM is bistable.

Nonvolatility and Local Activity
The nonvolatile property of the memristor allows the memristor to remain up-to-date when the power is off, which can be characterized by the power-off plot (POP) [66], i.e.,

Nonvolatility and Local Activity
The nonvolatile property of the memristor allows the memristor to remain up-todate when the power is off, which can be characterized by the power-off plot (POP) [66], i.e., the existence of two or more negative slopes intersecting the φ-axis.Letting v(n) = 0, the equilibrium equation can be obtained as follows: ( The dynamic curve of Equation ( 2) is drawn in Figure 3.According to the direction of the arrows shown in Figure 3, we can judge what state the internal state variable will eventually stabilize in different situations.It shows that the intersection points of the curve and φ-axis are Q 1 (−1, 0), Q 2 (0, 0), and Q 3 (1, 0), respectively.Q 1 and Q 3 are the two intersection points with negative slopes.Therefore, the nonvolatile memristor has two stable states under different initial states, that is to say,

Nonvolatility and Local Activity
The nonvolatile property of the memristor allows the memristor to re when the power is off, which can be characterized by the power-off plo the existence of two or more negative slopes intersecting the ϕ-axis.Let equilibrium equation can be obtained as follows: The dynamic curve of Equation ( 2) is drawn in Figure 3.According of the arrows shown in Figure 3, we can judge what state the internal st eventually stabilize in different situations.It shows that the intersecti curve and ϕ-axis are Q1 (−1, 0), Q2 (0, 0), and Q3 (1, 0), respectively.Q1 an intersection points with negative slopes.Therefore, the nonvolatile me stable states under different initial states, that is to say, According to the above analysis, it can be concluded that the memr tile.According to the above analysis, it can be concluded that the memristor is nonvolatile.Local activity is the root of complexity [67].The locally active nature of the memristor can be judged by the DC V-I diagram, i.e., the existence of at least one negative slope curve of the diagram [64].For a better illustration, let φ(n + 1) − φ(n) = 0, which gives the equilibrium equation In Equation (4), V is the DC voltage, I is the DC current, and Φ implies a variable equilibrium state of the equation and satisfies φ(n + 1) − φ(n) = 0, φ(n) = Φ.Taking the range of Φ as [−2, 2], the curve of the V-I plane is drawn in Figure 4, and it is clear that the memristor has two locally active regions.Then, the LADM has the characteristic of local activity.

 
In Equation (4), V is the DC voltage, I is the DC current, and Φ implies a variable equilibrium state of the equation and satisfies ϕ(n + 1) − ϕ(n) = 0, ϕ(n) = Φ.Taking the range of Φ as [−2, 2], the curve of the V-I plane is drawn in Figure 4, and it is clear that the memristor has two locally active regions.Then, the LADM has the characteristic of local activity.

The Neural Network Model and Equilibrium Point Analysis
The Rulkov neuron has the following main firing patterns: period-1 bursting firing, triangular bursting firing, square-wave bursting firing, etc.By coupling two Rulkov neurons with the proposed LADM, we can obtain a new LADM-based integer-order Rulkov neural network (IORNN) model as x

n y n k x n x n n xn y n y n x x n y n k x n x n n xn
where x1 and x2 denote the membrane potentials of neurons; y1 and y2 are the neuronal recovery variables; (x1(n) − x2(n)) means the membrane potential difference between two Rulkov neurons; αi (i = 1, 2), αi, and σ are parameters controlling the firing patterns of neurons, 0 < μ ≪ 1; and k represents the coupling strength.What needs to be noted is that the control parameters of the new neural network are fixed with μ = 0.001 and σ = −1 in this paper [68].

The Neural Network Model and Equilibrium Point Analysis
The Rulkov neuron has the following main firing patterns: period-1 bursting firing, triangular bursting firing, square-wave bursting firing, etc.By coupling two Rulkov neurons with the proposed LADM, we can obtain a new LADM-based integer-order Rulkov neural network (IORNN) model as where x 1 and x 2 denote the membrane potentials of neurons; y 1 and y 2 are the neuronal recovery variables; (x 1 (n) − x 2 (n)) means the membrane potential difference between two Rulkov neurons; α i (i = 1, 2), α i , and σ are parameters controlling the firing patterns of neurons, 0 < µ 1; and k represents the coupling strength.What needs to be noted is that the control parameters of the new neural network are fixed with µ = 0.001 and σ = −1 in this paper [68].
It is known that the equilibrium points for calculating Equation ( 5) are The eigenvalues can be obtained from P(λ) = λI − J, where I is an identity matrix.For the convenience of calculation, we assume that α 1 = α 2 .
It is obvious that there is a constant λ = 1.1 > 1, so E 3 is always unstable.In the following numerical simulation, it can be inferred that unstable points will lead to chaos.Similarly, for E 1 and E 2 , it can be calculated that the eigenvalues are the same, both of which are (0.8, 1, −0.5α, 1, −0.5α).Additionally, it is valid when the control parameter α is positive for the neural network.Therefore, the equilibrium points E 1 and E 2 are stable.
Discrete fractional calculus is suggested to describe neural networks with memory effects [69,70].In order to present the fractional-order model, some definitions and lemmas are given below.Definition 1.Following [61], for a given fractional-order q > 0, q / ∈ N, when x(t) is defined in N t 0 , its Caputo difference is defined by where , and ρ(s) denotes the next point in the time scale after s, namely ρ(s) = s + 1 for s ∈ N t 0 .Obviously, when q = 1, the fractional difference becomes ∆x( Theorem 1.Following [61], for the Caputo-like fractional-order system, The definition can continue to be equivalent to the equation Furthermore, let s + q = j; then, the equation can be changed to where x(t) is the system variable, and g(•) is the system equation.
Based on the Caputo fractional-order difference operator and Equation ( 5), we can obtain a new LADM-based fractional-order Rulkov neural network (FORNN) model as where x 1 (0) and x 2 (0) denote the initial membrane potentials of neurons, y 1 (0) and y 2 (0) are the initial neuronal recovery variables, and q is the fractional order.Other parameters are the same as those in Equation ( 5).

Multistability and Novel Firing Patterns
Some simulation results are presented in this section.These results indicate that the introduction of the LADM improves the complexity of neurons, and the firing states of the neural network can be selected by adjusting the initial values or the coupling strength k, which has a potential application value in engineering applications.

Transition of Firing Patterns
According to [25], it is assumed that the two neurons have the same firing patterns before coupling, which means α = α 1 = α 2 .As the control variable α of a neuron is set in the range of [2,5], and the coupling strength k is varied from 0 to 1, the spectral entropy complexity [71] of the IORNN and FORNN are shown in Figure 5. Different colors represent different degrees of complexity, and the darker the color, the more complex it will be.Brown and black imply chaotic or hyperchaotic firing patterns, and yellow indicates periodic firing behaviors.The diversity of colors indicates that neurons have rich firing behaviors, and the transition of firing patterns can be expressed by the variation in the colors with the change in the coupling strength.It is obvious that the types of firing patterns of the FORNN at q = 0.6 are less than those of the IORNN.Here, a few examples are presented to illustrate this point.For α = α1 = α2 = 3, the time sequences of the membrane potential x1 under diff coupling strengths are depicted in Figures 6-8.The firing patterns of the FORNN shown in Figures 6 and 7, while the firing patterns of the IORNN are presented in Figu When the fractional order q = 0.6, and the coupling strength k varies from 0 to 1, a tr tion occurs from a silent state to silent firing and periodic spike firing.When the fract order is 0.95, it can be seen that the firing patterns of the FORNN are basically the s as that of the IORNN.The periodic triangular bursting firing pattern is drawn in Fi 8a.However, in Figure 8b, the neurons exhibit a mixed firing pattern with chaotic t gular bursting firing and square-wave bursting firing.When k = 0.8, the neurons ex chaotic square-wave bursting firing, as presented in Figure 8c.In addition, the trans process corresponds to the SE complexity in Figure 5.For α = α 1 = α 2 = 3, the time sequences of the membrane potential x 1 under different coupling strengths are depicted in Figures 6-8.The firing patterns of the FORNN are shown in Figures 6 and 7, while the firing patterns of the IORNN are presented in Figure 8.When the fractional order q = 0.6, and the coupling strength k varies from 0 to 1, a transition occurs from a silent state to silent firing and periodic spike firing.When the fractional order is 0.95, it can be seen that the firing patterns of the FORNN are basically the same as that of the IORNN.The periodic triangular bursting firing pattern is drawn in Figure 8a.However, in Figure 8b, the neurons exhibit a mixed firing pattern with chaotic triangular bursting firing and square-wave bursting firing.When k = 0.8, the neurons exhibit chaotic square-wave bursting firing, as presented in Figure 8c.In addition, the transition process corresponds to the SE complexity in Figure 5.
order is 0.95, it can be seen that the firing patterns of the FORNN are basically the same as that of the IORNN.The periodic triangular bursting firing pattern is drawn in Figure 8a.However, in Figure 8b, the neurons exhibit a mixed firing pattern with chaotic triangular bursting firing and square-wave bursting firing.When k = 0.8, the neurons exhibi chaotic square-wave bursting firing, as presented in Figure 8c.In addition, the transition process corresponds to the SE complexity in Figure 5.  order is 0.95, it can be seen that the firing patterns of the FORNN are basically the sam as that of the IORNN.The periodic triangular bursting firing pattern is drawn in Figur 8a.However, in Figure 8b, the neurons exhibit a mixed firing pattern with chaotic trian gular bursting firing and square-wave bursting firing.When k = 0.8, the neurons exhibi chaotic square-wave bursting firing, as presented in Figure 8c.In addition, the transition process corresponds to the SE complexity in Figure 5.   9c.Similarly, when q = 0.95, the firing patterns of the FORNN are basically the same as that of the IORNN shown in Figures 10 and 11. Figure 11 shows the transition process of the firing patterns of the IORNN with α1 = 3 and α2 = 4.As shown in Figure 11a, the membrane potential x2 is in a silent state, while the membrane potential x1 is in a triangular bursting firing state at k = 0.With a slight increase in the coupling strength, the membrane potential x2 shown in Figure 11b starts to excite as a mixed firing pattern, which indicates that a synapse acts as an intermediary to successfully transmit information between neurons.In Figure 11c, the two Rulkov neurons reach synchronization at k = 0.8.It  9a,b, the membrane potentials x 1 and x 2 are in a silent state.When k = 0.8, the firing patterns of the neurons are periodic firing, as presented in Figure 9c.Similarly, when q = 0.95, the firing patterns of the FORNN are basically the same as that of the IORNN shown in Figures 10 and 11. Figure 11 shows the transition process of the firing patterns of the IORNN with α 1 = 3 and α 2 = 4.As shown in Figure 11a, the membrane potential x 2 is in a silent state, while the membrane potential x 1 is in a triangular bursting firing state at k = 0.With a slight increase in the coupling strength, the membrane potential x 2 shown in Figure 11b starts to excite as a mixed firing pattern, Fractal Fract.2023, 7, 82 9 of 18 which indicates that a synapse acts as an intermediary to successfully transmit information between neurons.In Figure 11c, the two Rulkov neurons reach synchronization at k = 0.8.It is clear that the coupling strength plays a decisive role in regulating the firing states of neurons, which also indicates that synaptic modulation plays an important role in the synchronization process.
in Figure 9c.Similarly, when q = 0.95, the firing patterns of the FORNN are basically the same as that of the IORNN shown in Figures 10 and 11. Figure 11 shows the transition process of the firing patterns of the IORNN with α1 = 3 and α2 = 4.As shown in Figure 11a, the membrane potential x2 is in a silent state, while the membrane potential x1 is in a triangular bursting firing state at k = 0.With a slight increase in the coupling strength, the membrane potential x2 shown in Figure 11b starts to excite as a mixed firing pattern, which indicates that a synapse acts as an intermediary to successfully transmit information between neurons.In Figure 11c, the two Rulkov neurons reach synchronization at k = 0.8.It is clear that the coupling strength plays a decisive role in regulating the firing states of neurons, which also indicates that synaptic modulation plays an important role in the synchronization process.  in Figure 9c.Similarly, when q = 0.95, the firing patterns of the FORNN are basically the same as that of the IORNN shown in Figures 10 and 11. Figure 11 shows the transition process of the firing patterns of the IORNN with α1 = 3 and α2 = 4.As shown in Figure 11a, the membrane potential x2 is in a silent state, while the membrane potential x1 is in a triangular bursting firing state at k = 0.With a slight increase in the coupling strength, the membrane potential x2 shown in Figure 11b starts to excite as a mixed firing pattern, which indicates that a synapse acts as an intermediary to successfully transmit information between neurons.In Figure 11c, the two Rulkov neurons reach synchronization at k = 0.8.It is clear that the coupling strength plays a decisive role in regulating the firing states of neurons, which also indicates that synaptic modulation plays an important role in the synchronization process.It is found that the firing patterns of the LADM-based FORNN are more complex than that of the LADM-based IORNN, but the dynamic behaviors of the IORNN are more convenient to characterize the information transmission process between biological neu rons with lower computation costs.Therefore, we will specify the dynamics of the IORNN model in the next section.It is found that the firing patterns of the LADM-based FORNN are more complex than that of the LADM-based IORNN, but the dynamic behaviors of the IORNN are more convenient to characterize the information transmission process between biological neurons with lower computation costs.Therefore, we will specify the dynamics of the IORNN model in the next section.

Novel Firing Patterns
Due to the nonvolatile and bistable advantages of the proposed LADM, the bidirectional memristive synapse-coupled neural network leads to novel firing patterns as well as mixed firing patterns, enriching the dynamical behaviors of neurons.
Figure 12 depicts the time sequences of the two membrane potentials x 1 and x 2 with different coupling strengths, respectively, and blue represents the membrane potential x 1 , while magenta represents the membrane potential x 2 .When α 1 = 2.5, α 2 = 2, and k = 0.1, a new periodic bursting firing pattern appears in Figure 12a, and its topological structure is different from the original Rulkov neuron.With the increase in the coupling strength to 0.2, a crossed firing pattern with different topologies of the bursting firing behavior is shown in Figure 12c.In Figure 12b,d, the membrane potential x 2 shows a mixed firing behavior of spike firing and bursting firing.With the increase in the coupling strength, the frequency of the spike firing pattern decreases.It is obvious that the introduction of the LADM enriches the firing patterns of neurons.
It is found that the firing patterns of the LADM-based FORNN are more comple than that of the LADM-based IORNN, but the dynamic behaviors of the IORNN are mor convenient to characterize the information transmission process between biological neu rons with lower computation costs.Therefore, we will specify the dynamics of the IORNN model in the next section.

Novel Firing Patterns
Due to the nonvolatile and bistable advantages of the proposed LADM, the bidirec tional memristive synapse-coupled neural network leads to novel firing patterns as we as mixed firing patterns, enriching the dynamical behaviors of neurons.
Figure 12 depicts the time sequences of the two membrane potentials x1 and x2 wit different coupling strengths, respectively, and blue represents the membrane potential x while magenta represents the membrane potential x2.When α1 = 2.5, α2 = 2, and k = 0.1, new periodic bursting firing pattern appears in Figure 12a, and its topological structure i different from the original Rulkov neuron.With the increase in the coupling strength t 0.2, a crossed firing pattern with different topologies of the bursting firing behavior i shown in Figure 12c.In Figure 12b,d  With different control variables, the time series of the neuronal membrane potentials x 1 and x 2 are depicted in Figure 13.Along with the change in the coupling strength, the membrane potential x 1 not only produces novel triangular firing behavior but also generates a mixed firing behavior of triangular bursting firing with different topologies, as shown in Figure 13a-c.More interestingly, Figure 13d-f indicate that the membrane potential x 2 has a mixed firing pattern with spike firing, period-1 bursting firing, and triangular bursting firing.The frequency of the triangular bursting firing increases with the increase in the coupling strength.These results lead to a similar conclusion that the introduction of the LADM makes the firing behaviors of Rulkov neurons more complicated.
(c) (d) With different control variables, the time series of the neuronal membrane potential x1 and x2 are depicted in Figure 13.Along with the change in the coupling strength, th membrane potential x1 not only produces novel triangular firing behavior but also gener ates a mixed firing behavior of triangular bursting firing with different topologies, a shown in Figure 13a-c.More interestingly, Figure 13d-f indicate the membrane po tential x2 has a mixed firing pattern with spike firing, period-1 bursting firing, and trian gular bursting firing.The frequency of the triangular bursting firing increases with th increase in the coupling strength.These results lead to a similar conclusion that the intro duction of the LADM makes the firing behaviors of Rulkov neurons more complicated.

Transient Chaotic Firing Behavior
The phase diagrams and time sequences are shown in Figure 14

Coexisting Firing Patterns
Coexisting behavior is a hot topic in the field of neural networks.The parameters of the neural network are relatively fixed, and the coexisting firing pattern can be obtained only by changing the initial values.The corresponding parameters are assigned as follows:

Coexisting Firing Patterns
Coexisting behavior is a hot topic in the field of neural networks.The parameters of the neural network are relatively fixed, and the coexisting firing pattern can be obtained only by changing the initial values.The corresponding parameters are assigned as follows:  When α = 2.2, and the coupling strength k is set to 0.5, the neuron exhibits the coexisting pattern of periodic and transient chaotic firing.The phase diagrams of the membrane potential x 1 under two sets of initial values (x 1 (0), y 1 (0), x 2 (0), y 2 (0), φ(0)) = (−1.1,−3, −1, −3, 0.8) and (−1, 0.5, −1, 0, 0) are drawn in Figure 15.The membrane potential x 1, presenting a periodic firing pattern, is shown in Figure 15a.In Figure 15b, the membrane potential x 1 presents a transient chaotic firing pattern, and the blue track represents periodic firing, while the green track represents chaotic firing.When α = 2.5, and the coupling strength k is set to 0.1, there is the coexisting behavior of period-1 bursting firing and triangular bursting firing in Rulkov neurons.Figure 16 shows the time series and phase diagrams with different initial values, in which red and blue correspond to the initial values of (1, 0.5, 1, 0.51, 0.3) and (1, 0.5, 1, 0, 1), respectively.It can be seen that when the initial values are chosen as (1, 0.5, 1, 0, 1), the neuron exhibits period-1 bursting firing.When the initial values are changed to (1, 0.5, 1, 0.51, 0.3), the firing pattern switches to triangular bursting firing.It is clear that the multiple firing behaviors of the neural network are very sensitive to the initial values.

Phase Synchronization and Synchronization Transition
Each time the neuron begins to burst, the neuronal recovery variable y reaches its maximum value.For a better analysis of phase synchronization, the phase of the recovery variable y is discussed here, giving the definition of phase θ [24] where nk is the time of the appearance of the kth bursting.
The phase difference can be expressed as

Phase Synchronization and Synchronization Transition
Each time the neuron begins to burst, the neuronal recovery variable y reaches its maximum value.For a better analysis of phase synchronization, the phase of the recovery variable y is discussed here, giving the definition of phase θ [24] where nk is the time of the appearance of the kth bursting.
The phase difference can be expressed as

Phase Synchronization and Synchronization Transition
Each time the neuron begins to burst, the neuronal recovery variable y reaches its maximum value.For a better analysis of phase synchronization, the phase of the recovery variable y is discussed here, giving the definition of phase θ [24] θ where n k is the time of the appearance of the kth bursting.
The phase difference can be expressed as When α 1 = 2 and α 2 = 3, with the coupling strength k increasing from 0.1 to 0.5, the phenomenon of synchronization transition occurs, which is different from all the previous synchronization transitions.Figure 17a,c,e,g show firing patterns of two neurons and Figure 17b,d,f,h present phase difference between two neurons at k = 0.1, 0.25, 0.3, and 0.5 respectively.According to the time sequences shown in Figure 17a, it can be realized that the two neurons exhibit anti-phase synchronization at k = 0.1, and the firing patterns of the two Rulkov neurons are period-1 bursting firing and triangular bursting firing, respectively.The phase difference is maintained at a constant value in Figure 17b, which reflects the anti-phase synchronization behavior.When the coupling strength increases to 0.25, the firing patterns of the Rulkov neurons are mixed firing patterns, as presented in Figure 17c.As shown in Figure 17d, the phase difference is similar to a period function and fluctuates within the interval of [0, 3], which indicates that the Rulkov neurons are partially synchronized.Comparing Figure 17f with Figure 17d, it is obvious that the overall trend of the phase difference is similar, but it has a short stay at [4800, 5600].When k = 0.5, the firing patterns of the two Rulkov neurons are triangular bursting firing, as presented in Figure 17g.Furthermore, it can be observed from Figure 17h that the phase difference between the two neurons remains constant, demonstrating the phenomenon of synchronization.Through the whole process, it is found that the firing patterns of the two Rulkov neurons change.
phenomenon of synchronization transition occurs, which is different from all the previous synchronization transitions.Figure 17a,c,e,g show firing patterns of two neurons and Figure 17b,d,f,h present phase difference between two neurons at k = 0.1, 0.25, 0.3, and 0.5 respectively.According to the time sequences shown in Figure 17a, it can be realized that the two neurons exhibit anti-phase synchronization at k = 0.1, and the firing patterns of the two Rulkov neurons are period-1 bursting firing and triangular bursting firing, respectively.The phase is maintained at a constant value in Figure 17b, which reflects the anti-phase synchronization behavior.When the coupling strength increases to 0.25, the firing patterns of the Rulkov neurons are mixed firing patterns, as presented in Figure 17c.As shown in Figure 17d, the phase difference is similar to a period function and fluctuates within the interval of [0, 3], which indicates that the Rulkov neurons are partially synchronized.Comparing Figure 17f with Figure 17d, it is obvious that the overall trend of the phase difference is similar, but it has a short stay at [4800, 5600].When k = 0.5, the firing patterns of the two Rulkov neurons are triangular bursting firing, as presented in Figure 17g.Furthermore, it can be observed from Figure 17h that the phase difference between the two neurons remains constant, demonstrating the phenomenon of synchronization.Through the whole process, it is found that the firing patterns of the two Rulkov neurons change.In addition, with α1 = 3 and α2 = 3.5, the synchronization transition behavior is depicted in Figure 18.When k increases from 0.1 to 0.15, the transition from anti-phase synchronization to in-phase synchronization takes place.Figure 18a,c  In addition, with α 1 = 3 and α 2 = 3.5, the synchronization transition behavior is depicted in Figure 18.When k increases from 0.1 to 0.15, the transition from anti-phase synchronization to in-phase synchronization takes place.Figure 18a,c  In addition, with α1 = 3 and α2 = 3.5, the synchronization transition behavior is depicted in Figure 18.When k increases from 0.1 to 0.15, the transition from anti-phase synchronization to in-phase synchronization takes place.Figure 18a,c

Conclusions
In this paper, a simple neural network was proposed, in which an LADM was used to mimic the synapse connecting two Rulkov neurons.Because of its bistable, nonlinear, and biomimetic characteristics, the LADM enriched the dynamic behaviors of the neural

Figure 4 .
Figure 4. DC V-I diagram of the bistable LADM.

Figure 4 .
Figure 4. DC V-I diagram of the bistable LADM.

Figure 8 .
Figure 8. Transition of firing patterns of the IORNN with α = 3: (a) k = 0.1; (b) k = 0.8; (c) k = 0.9.In Figures 9-11, the blue dots and lines indicate the membrane potential x 1 , while the magenta dots and lines indicate the membrane potential x 2 .Figures 9 and 10 show the transition process of the firing patterns of the FORNN with α 1 = 3 and α 2 = 4.With q = 0.6, k = 0, and k = 0.01, as plotted in Figure9a,b, the membrane potentials x 1 and x 2 are in a silent state.When k = 0.8, the firing patterns of the neurons are periodic firing, as presented in Figure9c.Similarly, when q = 0.95, the firing patterns of the FORNN are basically the same as that of the IORNN shown in Figures10 and 11.Figure11shows the transition process of the firing patterns of the IORNN with α 1 = 3 and α 2 = 4.As shown in Figure11a, the membrane potential x 2 is in a silent state, while the membrane potential x 1 is in a triangular bursting firing state at k = 0.With a slight increase in the coupling strength, the membrane potential x 2 shown in Figure11bstarts to excite as a mixed firing pattern,
2 and α1 = α2 = α= 2.5 as examples, respectively, to reveal two different coexisting firing patterns.When α = 2.2, and the coupling strength k is set to 0.5, the neuron exhibits the coexisting pattern of periodic and transient chaotic firing.The phase diagrams of the membrane potential x1 under two sets of initial values (x1(0), y1(0), x2(0), y2(0), ϕ(0)) = (−1.1,−3, −1, −3, 0.8) and (−1, 0.5, −1, 0, 0) are drawn in Figure 15.The membrane potential x1, presenting a periodic firing pattern, is shown in Figure15a.In Figure15b, the membrane potential x1 presents a transient chaotic firing pattern, and the blue track represents periodic firing, while the green track represents chaotic firing.When α = 2.5, and the coupling strength k is set to 0.1, there is the coexisting behavior of period-1 bursting firing and triangular bursting firing in Rulkov neurons.Figure16shows the time series and phase diagrams with different initial values, in which red and blue correspond to the initial values
,e show firing patterns of two neurons and Figure 18b,d,f present phase difference between two neurons with k = 0.1, 0.108, and 0.15 respectively.When the coupled strength k = 0.1, the two Rulkov neu-
,e show firing patterns of two neurons and Figure18b,d,f present phase difference between two neurons with k = 0.1, 0.108, and 0.15 respectively.When the coupled strength k = 0.1, the two Rulkov neurons exhibit anti-phase synchronization with different triangular bursting firing patterns, as presented in Figure18a,b.Figure18cpresents the intermediate transition process.Under different iteration times, the neurons present different synchronization behaviors.It can also be characterized by the phase difference shown in Figure18d.When the coupled strength k = 0.15, the two neurons are exactly synchronized, as shown in Figure18e,f.In the whole process, it is clear that the two Rulkov neurons do not produce new firing patterns as the coupling strength increases.