SYNCHRONIZATION CONDITIONS IN NETWORKS OF HINDMARSH-ROSE SYSTEMS

The algebraic connectivity is crucial parameter in studying of synchronization of diffusively coupled networks. This paper studies the synchronization in networks of Hindmarsh-Rose systems, which is one of the most used neuron models. It presents sufﬁcient condition for synchronization in these networks using the Lyapunov function method. This is a simple condition which depends on the algebraic connectivity and on the parameters of the individual system. Numerical examples are presented to illustrate the obtained results.


Introduction
Synchronization is a fundamental phenomenon and plays an important role in various fields of science and technology [Blekhman, 1988;Osipov et al., 2007]. It is often observed in biological systems from the level of individual cells, physiological subsystems to organisms and their populations [Glass, 2001;Pikovsky et al., 2003]. This phenomenon is necessary for the normal functioning of various body systems, for example, for the functioning of the lungs or the heart, when cells oscillate synchronously and develop a macroscopic rhythm by this ways, which controls respiration, heartbeat, etc [Peskin, 1975]. It also plays a crucial role in the functioning of the brain, where it is associated with various cognitive abilities [Fries, 2005]. On the other hand, the appearance of areas of synchronous activity of neurons can lead to several pathologies, for example, in the case of Parkinson's disease, blocking of many neurons leads to tremor activity [Hammond et al., 2007]. The important role of synchronization in information processing in the brain makes it an important research topic in medical and computational neuroscience.
From dynamical-systems viewpoint, synchronization depends on some key factors which include the dynamics of individual systems and the type, strength and the topology of the interconnection between the nodes. The simplest case of the network is a network of linear systems with linear couplings. Such networks are well studied, and there are a lot of different results about linear networks synchronization, just to mention a few [Dzhunusov and Fradkov, 2011;Furtat et al., 2014]. However, the presence of nonlinearities in networks makes a problem harder. The real world systems are always nonlinear, which also applies to networks of neural models (see eg. [Guzenko et al., 2013;Andreev and Maksimenko, 2019]), on which this paper is focused. For the networks under diffusive coupling the algebraic connectivity, i.e. the second minimal eigenvalue of Laplace matrix, plays the crucial role in networks dynamics [Steur et al., 2009;Panteley and Loría, 2017]. In particular, the nonlinear networks of diffusively coupled semi-passive systems synchronize for large enough values of algebraic connectivity [Steur et al., 2009]. For various network nodes the value of algebraic connectivity needed for synchronization can be different. In [Plotnikov and Fradkov, 2019] the approach to obtain synchronization conditions of heterogeneous FitzHugh-Nagumo systems [FitzHugh, 1961;Nagumo et al., 1962] is presented. It is based on the Lyapunov function methods and allows to obtain simple sufficient condition which guarantees the network synchronization. This paper follows the same approach but for the networks of Hindmarsh-Rose (HR) systems [Hindmarsh and Rose, 1984]. This model is also describes the dynamics of individual neuron, but it is described by a nonlinear system of third-order differential equations, which makes the problem more difficult. Note that the synchroniza-tion condition for directly coupled HR networks was obtained in [Checco et al., 2009], while some results about bifurcations in two coupled HR systems can be found in [Jaimes-Reátegui et al., 2017;Raznoglazova and Plotnikov, 2020]. Sufficient synchronization conditions can help easily conclude will the network synchronized or not. In the case of synchronization absence one can use control algorithms to ensure it [Semenov and Fradkov, 2021;Rehak and Lynnyk, 2021].
The rest of the paper is organized as follows. Section 2 introduces HR model. In Sec. 3 the conditions of HR network synchronization are obtained. Section 4 provides the simulation results of HR network dynamics. Finally, the conclusion is given in Sec. 5.

Model
The HR was initially described by a nonlinear system of second order differential equations and was a simplified version of the Hodgkin-Huxley model [Hodgkin and Huxley, 1952]. However, after a series of experiments, new modes of neuron functioning were discovered that were not taken into account in the model [Hindmarsh and Rose, 1982]. Therefore, the authors J. Hindmarsh and M. Rose included the third equation, which made it possible to take into account most of the new regimes of behavior of the biological neuron. The HR system is a model of a biological neuron and is described by a nonlinear system of third-order differential equations [Hindmarsh and Rose, 1984]: where x describes the dynamics of the membrane potential, while y and z describe the dynamics of ionic currents. Parameter 0 < r 1 separates fast and slow dynamics, meaning that z describes the dynamics of a slow potassium current, while y describes the dynamics a fast sodium current. I is an external current, and a, b, c, d, s, w are constant positive parameters.

Main Result
Consider the network of N diffusively coupled HR systems: (2) where where C = (c ij ), i, j = 1, . . . , N is the adjacency matrix of graph G which corresponds to the considered network. Suppose that the graph G is connected and undirected.
To study the network (2) synchronization find the mean-field dynamics by addition of all equations (2) and dividing the sum by N : Introduce the coordinate synchronization errors then the synchronization error system dynamics can be described by the following equations: Thus, the problem of studying network (2) synchronization can be reduced to studying of synchronization error system (6) stability. For this purpose introduce the following Lyapunov function: where ∆ = col(δ 11 , δ 12 , δ 13 , . . . , δ 1N , δ 2N , δ 3N ), and µ is a positive constant, which will be defined later. Find the derivative of the Lyapunov function with respect to the system (6) and make some transformations: where β is some constant parameter, which will be also defined later. Since N i=1 δ ki = 0, k = 1, 2, 3, these transformations are valid.
Consider the difference of cubes , where x 2 i + x ix +x 2 is nonnegative function. Let find the constant α, for which the following inequality holds: The inequality (9) can be presented in a matrix form: This inequality is fulfilled if and only if the corresponding matrix is nonnegative definite. Using Sylvester criterion we obtain that α should lie in the interval α ∈ [0; 0.25] to ensure the fulfillment of the inequality (9). Also consider the difference of squares x 2 i − βx 2 . Let find the constant γ, for which the following inequality holds: The inequality (10) can be also presented in a matrix form: This inequality is fulfilled if and only if the corresponding matrix is nonnegative definite. Using Sylvester criterion we obtain that γ > 1 and β > γ/(γ − 1). Choosing α = 0.25 the Lyapunov function derivative (8) can be estimated as: Now apply the inequality 2|uv| ≤ λu 2 + v 2 /λ for any variables u and v and λ > 0 to the cross-terms of (11): where ν 1 > 0, ν 2 > 0, ν 3 > 0 are positive constants, which will be defined later. Using obtained inequalities (12) the Lyapunov function derivative (11) can be estimated as follows: The coefficients before terms δ 4 1i and δ 2 2i in (13) should be negative to make the Lyapunov function derivative be negative. The coefficient before the term δ 2 3i is negative, while the coefficient before the term δ 2 1i is always positive. However, there is also a coupling term, which is a quadratic term depending on δ 1i , i = 1, . . . , N . Therefore one can make the sum of quadratic terms depending on δ 1i , i = 1, . . . , N be negative. Thus, one obtains the following minimization problem: which can be rewritten as: The minimal solution can be reached if ν 2 = ν 1 /(bγ), therefore the estimates of ν 2 and ν 1 /(bγ) should be equal and maximal. Thus, one gets: where µ can be expressed to obtain: Meaning that ν 3 > 0, one can find the optimal value of ν 3 , which is a local maximum of the function (16): Thus, other unknown parameters can also be calculated: . (17)  Now consider the sum of the second and the third terms in (13) meaning that ν 2 = ν 1 /(bγ) and formulas (5): where d i is an in-degree of ith node, L(G 0 ) is the Laplace matrix of the complete graph G 0 , and L(G) is the Laplace matrix of the connectivity graph G.

Simulation
For simulation let consider diffusively coupled network of N = 50 HR systems (2). Parameters of each node are the following: a = 1, b = 2, c = 1, d = 0.5, r = 0.01, s = 4, w = 1.5 and I = 5. Consider symmetric sparse matrix with density 0.9 as the adjacency matrix of the network (2). This matrix has approximately 0.3N 2 nonzero entries, which are uniformly distributed on the interval (0; 12), i.e. the corresponding graph is weighted and undirected. The algebraic connectivity of presented graph λ 2 (G) is equal to 15.4172. Theorem 1 guarantees the HR network synchronization, if λ 2 (G) > 8.5 for systems with considered parameters. This means that for the considered graph one obtains synchronization. Figure 1 presents presents the results of the simulation. For Fig. 1 x-axis corresponds to the time, while y-axis corresponds to the node number. The values of state variables of the nodes are marked by the color. One can see that for chosen parameters there is synchronization between all the state coordinates. Now consider the graph with algebraic connectivity λ 2 (G) = 0.001. Theorem 1 does not guarantee the HR network synchronization for such λ 2 (G). One can see the results of simulation in Fig. 2: there is no synchronization between the state coordinates. Thus, the inequality (20) can be used to check if the HR network is synchronized.

Conclusion
Here the synchronization problem in diffusively coupled HR network is studied. It is known that the algebraic connectivity is crucial parameter for studying the synchronization of diffusively coupled networks. For large enough value of algebraic connectivity the network of semi-passive diffusively coupled system will synchronize. The interesting question is to find the estimate of the algebraic connectivity ensuring network synchronization. Here the Lyapunov function method is applied to study the stability of synchronization error system. The simple condition guaranteeing network synchronization is obtained. This condition depends on the algebraic connectivity and on the parameters of the individual system. Also the simulations have been made to illustrate the obtained condition. Some features presented here which are used for the estimation of the derivative of the Lyaponov function can also be used to study the stability of other nonlinear systems.