Synchronization of Memristive FitzHugh-Nagumo Neural Networks
A new mathematical model of neural networks described by diffusive
FitzHugh-Nagumo equations with memristors and linear synaptic coupling is
proposed and investigated. The existence of absorbing set for the solution
semiflow in the energy space is proved and global dynamics of the memristive
neural networks are dissipative. Through uniform estimates and maneuver of
integral inequalities on the interneuron difference equations, it is shown that
exponential synchronization of the neural network at a uniform convergence rate
occurs if the coupling strength satisfies a threshold condition explicitly
expressed by the system parameters, which is illustrated by an example and
numerical simulation experiments.
Yuncheng You, Jing Tian, Junyi Tu
# Synchronization of memristive FitzHugh-Nagumo neural networks

###### Abstract.

A new mathematical model of neural networks described by diffusive FitzHugh-Nagumo equations with memristors and linear synaptic coupling is proposed and investigated. The existence of absorbing set for the solution semiflow in the energy space is proved and global dynamics of the memristive neural networks are dissipative. Through uniform estimates and maneuver of integral inequalities on the interneuron difference equations, it is shown that exponential synchronization of the neural network at a uniform convergence rate occurs if the coupling strength satisfies a threshold condition explicitly expressed by the system parameters, which is illustrated by an example and numerical simulation experiments.

Key words and phrases:Memristive FitzHugh-Nagumo equations, dissipative dynamics, exponential synchronization, coupling strength, neural network 2010 Mathematics Subject Classification: 35B40, 35B41, 35K55, 37L30, 92B20

## 1. Introduction

Recently the global dynamics and exponential synchronization of the neural networks modeled by the diffusive Hindmarsh-Rose equations with memristors were proposed and studied by the first author in [37, 38]. In this paper, we shall consider a new mathematical model of neural networks described by the diffusive FitzHugh-Nagumo equations with memristors and linear synaptic interneuron coupling.

Let a network of \(m\) fully coupled memristive neuron cells be denoted by \(\mathcal{NW}=\{\mathcal{N}_{i}:i=1,2,\cdots,m\}\), where \(m\geq 2\) is a positive integer, which is described by the following model of memristive and diffusive FitzHugh-Nagumo (FHN) equations. Each neuron \(\mathcal{N}_{i},1\leq i\leq m\), in this network is presented by the differential equations:

\[\begin{split}&\frac{\partial u_{i}}{\partial t}=\eta\Delta u_{i}+f(u_{ i},x)-\sigma w_{i}+J-k\tanh(\rho_{i})u_{i}+\sum_{j=1}^{m}P(u_{j}-u_{i}),\\ &\frac{\partial w_{i}}{\partial t}=au_{i}+c-bw_{i},\\ &\frac{\partial\rho_{i}}{\partial t}=qu_{i}-r\rho_{i},\end{split} \tag{1.1}\]
for \(t>0\), \(x\in\Omega\subset\mathbb{R}^{n}\) (\(n\leq 3\)), where \(\Omega\) is a bounded domain with locally Lipschitz continuous boundary \(\partial\Omega\).

In the membrane potential \(u_{i}\)-equations, the nonlinear term \(k\tanh(\rho_{i})u_{i}\) presents the memristive coupling effect [24, 33, 36], where \(\rho_{i}(t,x)\) stands for the memductance of the memristor and \(\tanh(\rho_{i})\) is the the electromagnetic induction flux with its coupling strength coefficient \(k\). In this system, the fast excitatory variable \(u_{i}(t,x)\) refers to the transmembrane electrical potential of a neuron cell and the slow recovering variable \(w_{i}(t,x)\) represents the integrated ionic current across the neuron membrane. The network neuron coupling terms are assumed to be linear with a common strength coefficient \(P\) in the membrane potential equation.

We impose the homogeneous Neumann boundary condition

\[\frac{\partial u_{i}}{\partial\nu}(t,x)=0,\quad\text{for}\;\;t>0,\;x\in \partial\Omega,\quad 1\leq i\leq m. \tag{1.2}\]

The initial states of the system (1.1) will be denoted by

\[u_{i}^{0}(x)=u_{i}(0,x),\;w_{i}^{0}(x)=w_{i}(0,x),\;\rho_{i}^{0}=\rho_{i}(0,x),\;\;1\leq i\leq m. \tag{1.3}\]

The following Assumption is made on the scalar function \(f\in C^{1}(\mathbb{R}\times\Omega)\):

\[\begin{split}& f(s,x)s\leq-\lambda|s|^{4}+\varphi(x),\quad s \in\mathbb{R},\;x\in\Omega,\\ &\frac{\partial f}{\partial s}(s,x)\leq\beta,\quad s\in\mathbb{R },\;x\in\Omega,\end{split} \tag{1.4}\]

where \(\lambda\) and \(\beta\) are positive constants, \(\varphi\in L^{2}(\Omega)\) is a given function. Note that the prototype nonlinear homogeneous function \(f(s)=s(s-\kappa)(1-s)\) with \(\kappa>0\) in the original FitzHugh-Nagumo equations [9] satisfies the properties (1.4), in particular,

\[\begin{split}& f(s)s=-s^{4}+(1+\kappa)s^{3}-\kappa s^{2}\leq- \frac{1}{4}s^{4}+\frac{1}{4}(1+\kappa)^{4},\\ & f^{\,\prime}(s)=-3s^{2}+2(1+\kappa)s-\kappa\leq\frac{1}{3}(1+ \kappa)^{2}-\kappa\leq\beta=\frac{1}{3}(1+\kappa)^{2}.\end{split}\]

All the parameters \(\eta,\sigma,k,a,c,b,q,r\), and \(P\) can be any positive constants, while the reference membrane potential \(J\) can be any real number constant.

Typical models of neuron dynamics are four-dimensional Hodgkin-Huxley equations [12], two-dimensional FitzHugh-Nagumo equations [9], and three-dimensional Hindmarsh-Rose equations [11], which originally consist of ordinary differential equations without memristors and characterize the periodic firing-bursting dynamics for neurons and nerve systems. Analysis through Hopf bifurcations and energy or Hamiltonian methods with semi-numerical simulations are the main approaches to show many solution patterns and collective synchronization behavior [6, 14, 20, 39].
Global dynamics and synchronization of ensemble neurons and neural networks modeled by partly diffusive Hindmarsh-Rose equations and FitzHugh-Nagumo equations have been studied by the authors' group in recent years [20, 21, 25, 26]. These models are hybrid differential equations and reflect the structural feature of neuron cells, which contain short-branch dendrites receiving incoming signals and long-branch axons propagating and transmitting outgoing signals through synapses.

The concept of memristor was coined by Chua [4] to describe the effect of electromagnetic flux on moving electric charges. Physical and generic memristive systems initially reported in [5, 30] attracted broad scientific interests in the recent decade, especially recognized in biological neuron models and artificial intelligence computing [1, 7, 15, 17, 29, 31] as a new type (other than electrical and chemical) synapsis or as an ideal component which has the nonvolatile properties and can process dynamically memorized signal information to exhibit more complex or chaotic dynamics in neural networks. Memristor-based mathematical models now penetrate many fields with applications to image encryption, DNA sequences operation, brain criticality, cell physiology, cybersecurity, and quantum computers, cf. [17, 18, 23, 28, 29, 33, 34].

The researches on memristive FitzHugh-Nagumo and Hindmarsh-Rose neural networks in ordinary differential equations have been expanding in the recent decade, cf. [1, 8, 16, 24, 40] and many references therein. Various synchronization results with memristive effect of these models are achieved [10, 13, 19, 22, 31, 32, 35] mainly by the methods of generalized Hamiltonian functions, Lyapunov exponents, and the computational algebra with numerical simulations.

In this work we shall rigorously prove a threshold condition on the neuron coupling strength \(P\) to ensure an exponential synchronization of the memristive neural networks (1.1) through the approach of dissipative dynamical analysis and sharp uniform estimates, which can be extended to study complex or artificial neural networks.

## 2. **Formulation and Preliminaries**

Define two Hilbert spaces of functions:

\[E=[L^{2}(\Omega,\mathbb{R}^{3})]^{m}\quad\text{and}\quad\Pi=[H^{1}(\Omega) \times L^{2}(\Omega,\mathbb{R}^{2})]^{m}\]

where \(H^{1}(\Omega)\) is a Sobolev space. Call \(E\) the energy space and \(\Pi\) the regular space. The norm and inner-product of \(L^{2}(\Omega)\) or \(E\) will be denoted by \(\|\,\cdot\,\|\) and \(\langle\,\cdot,\cdot\,\rangle\), respectively. We use \(|\,\cdot\,|\) to denote a vector norm or a set measure in \(\mathbb{R}^{n}\). The initial-boundary value problem (1.1)-(1.3) can be formulated into an initial value problem of the evolutionary equation:

\[\begin{split}&\frac{\partial g}{\partial t}=Ag+F(g),\,\,\,t>0,\\ & g(0)=g^{0}\in E.\end{split} \tag{2.1}\]
The unknown function is a column vector \(g(t)=\operatorname{col}\,(g_{1}(t),g_{2}(t),\cdots,g_{m}(t))\), where

\[g_{i}(t)=\operatorname{col}\,(u_{i}(t,\cdot),\,w_{i}(t,\cdot),\,\rho_{i}(t, \cdot)),\quad 1\leq i\leq m,\]

characterizes the dynamics of the neuron \(\mathcal{N}_{i}\). The initial data function in (2.1) is

\[g(0)=g^{0}=\operatorname{col}\,(g_{1}^{0},\,g_{2}^{0},\cdots,g_{m}^{0})\quad \text{where}\;\;g_{i}^{0}=\operatorname{col}\,(u_{i}^{0},\,w_{i}^{0},\,\rho_{i} ^{0}),\;1\leq i\leq m.\]

The energy norm \(\|g(t)\|\) of the solution for the evolutionary equation (2.1) in the space \(E\) is given by

\[\|g(t)\|^{2}=\sum_{i=1}^{m}\|g_{i}(t)\|^{2}=\sum_{i=1}^{m}\left(\|u_{i}(t)\|^{ 2}+\|w_{i}(t)\|^{2}+\|\rho_{i}(t)\|^{2}\right).\]

The closed linear operator \(A\) in (2.1) is defined by \(A=\operatorname{diag}\,(A_{1},A_{2},\cdots,A_{m})\), where

\[A_{i}=\begin{pmatrix}\eta\Delta&0&0\\ 0&-bI&0\\ 0&0&-rI\end{pmatrix}:\mathcal{D}(A)\to E,\quad i=1,2,\cdots,m, \tag{2.2}\]

with the domain \(\mathcal{D}(A)=\{g\in[H^{2}(\Omega)\times L^{2}(\Omega,\mathbb{R}^{2})]^{m}: \partial u_{i}/\partial\nu=0,\,1\leq i\leq m\}\), is the generator of the \(C_{0}\)-semigroup \(\{e^{At}\}_{t\geq 0}\) on the space \(E\) and \(I\) is the identity operator. By the fact that \(H^{1}(\Omega)\hookrightarrow L^{6}(\Omega)\) is a continuous imbedding for space dimension \(n\leq 3\) and by the Assumption (1.4), the nonlinear mapping

\[F(g)=\begin{pmatrix}f(u_{1},x)-\sigma w_{1}+J-k\tanh(\rho_{1})u_{1}+\sum_{j=1 }^{m}P(u_{j}-u_{1})\\ au_{1}+c\\ qu_{1}\\ \vdots\\ f(u_{m},x)-\sigma w_{m}+J-k\tanh(\rho_{m})u_{m}+\sum_{j=1}^{m}P(u_{j}-u_{m}) \\ au_{m}+c\\ qu_{m}\end{pmatrix}:\Pi\longrightarrow E \tag{2.3}\]

is a locally Lipschitz continuous mapping.

In this work we shall consider the weak solutions [3, Section XV.3] of this initial value problem (2.1).

**Definition 2.1**.: A \(3m\)-dimensional vector function \(g(t,x)\), where \((t,x)\in[0,\tau]\times\Omega\), is called a weak solution to the initial value problem of the evolutionary equation (2.1), if the following two conditions are satisfied:

(i) \(\frac{d}{dt}\langle g,\xi\rangle=\langle Ag,\xi\rangle+\langle F(g),\xi\rangle\) is satisfied for a.e. \(t\in[0,\tau]\) and any \(\xi\in E^{*}=E\).

(ii) \(g(t,\cdot)\in C([0,\tau];E)\cap C^{1}((0,\tau);E)\) and \(g(0)=g^{0}\).

Here \(E^{*}\) is the dual space of the Hilbert space \(E\).
The following proposition can be proved by the Galerkin approximation method [3] with the regularity property [27] of the parabolic operator semigroup \(e^{At}\).

**Proposition 2.2**.: _For any initial state \(g^{0}\in E\), there exists a unique weak solution \(g(t;g^{0}),\,t\in[0,\tau]\), for some \(\tau>0\) may depending on \(g^{0}\), of the initial value problem (2.1) formulated from the memristive FitzHugh-Nagumo equations (1.1). The weak solution \(g(t;g^{0})\) continuously depends on the initial data and satisfies_

\[g\in C([0,\tau];E)\cap C^{1}((0,\tau);E)\cap L^{2}((0,\tau);\Pi). \tag{2.4}\]

_Moreover, for any initial state \(g^{0}\in E\), the weak solution \(g(t;g^{0})\) becomes a strong solution for \(t\in(0,\tau)\), which has the regularity_

\[g\in C((0,\tau];\Pi)\cap C^{1}((0,\tau);\Pi). \tag{2.5}\]

An infinite dimensional dynamical systems [3, 27] for time \(t\geq 0\) only is called a semiflow. Absorbing set defined below is the key concept to characterize dissipative dynamics of a semiflow.

**Definition 2.3**.: Let \(\{S(t)\}_{t\geq 0}\) be a semiflow on a Banach space \(\mathscr{X}\). A bounded set \(B^{*}\) of \(\mathscr{X}\) is called an absorbing set of this semiflow, if for any given bounded set \(B\subset\mathscr{X}\) there exists a finite time \(T_{B}\geq 0\) depending on \(B\), such that \(S(t)B\subset B^{*}\) for all \(t>T_{B}\). The semiflow is called dissipative if there exists an absorbing set.

The Young's inequality in a generic form will be used throughout this paper: For any two positive numbers \(x\) and \(y\), if \(\frac{1}{p}+\frac{1}{q}=1\) and \(p>1,q>1\), one has

\[x\,y\leq\frac{1}{p}\varepsilon x^{p}+\frac{1}{q}C(\varepsilon,p)\,y^{q}\leq \varepsilon x^{p}+C(\varepsilon,p)\,y^{q},\quad C(\varepsilon,p)=\varepsilon^{ -q/p}, \tag{2.6}\]

where constant \(\varepsilon>0\) can be arbitrarily small. In Section 4, the Gagliardo-Nirenberg interpolation inequalities [27, Theorem B.3] will be used in a crucial step to prove the main result on exponential synchronization of the memristive neural networks.

## 3. **Dissipative Dynamics of Memristive FitzHugh-Nagumo Semiflow**

In this section, we shall prove the global existence of weak solutions in time for the initial value problem (2.1) to establish a solution semiflow of the memristive FitzHugh-Nagumo neural network modeled by (1.1). Then we show the dissipative dynamics in terms of the existence of an absorbing set of this semiflow in the state spaces \(E\).

**Theorem 3.1**.: _For any initial state \(g^{0}\in E\), there exists a unique global weak solution in time, \(g(t;g^{0})=\operatorname{col}\left(u_{i}(t),w_{i}(t),\rho_{i}(t):1\leq i\leq m\right),\,t\in[0,\infty)\), to the initial value problem (2.1) of the memristive FitzHugh-Nagumo equations (1.1) for the neural network \(\mathcal{N}\mathcal{N}\)._
Proof.: Conduct the \(L^{2}\) inner-products of the \(u_{i}\)-equation with \(C_{1}u_{i}(t,x)\) for \(1\leq i\leq m\), with the scaling constant \(C_{1}>0\) to be chosen later. Then sum them up to get

\[\begin{split}&\frac{C_{1}}{2}\frac{d}{dt}\sum_{i=1}^{m}\|u_{i}\|^{2}+C _{1}\eta\,\sum_{i=1}^{m}\|\nabla u_{i}\|^{2}=-C_{1}P\,\sum_{i=1}^{m}\sum_{j=1}^ {m}\int_{\Omega}(u_{i}-u_{j})^{2}\,dx\\ &+C_{1}\sum_{i=1}^{m}\int_{\Omega}(f(u_{i},x)u_{i}-\sigma u_{i}w_ {i}+Ju_{i}-k\tanh(\rho_{i})u_{i}^{2})\,dx\\ \leq&\,C_{1}\sum_{i=1}^{m}\int_{\Omega}\left[-\lambda |u_{i}|^{4}+|\varphi(x)|-\sigma u_{i}w_{i}+Ju_{i}+k|\tanh(\rho_{i})|u_{i}^{2} \right]dx\\ \leq&\,C_{1}\sum_{i=1}^{m}\int_{\Omega}\left[- \lambda u_{i}^{4}+|\varphi(x)|+\frac{1}{2}\left(\lambda u_{i}^{2}+\frac{\sigma ^{2}}{\lambda}w_{i}^{2}\right)+\frac{1}{2}\left(\frac{J^{2}}{\lambda}+\lambda u _{i}^{2}\right)+ku_{i}^{2}\right]dx\\ =&\,C_{1}\sum_{i=1}^{m}\int_{\Omega}((\lambda+k)u_{ i}^{2}-\lambda u_{i}^{4})\,dx+\frac{C_{1}\sigma^{2}}{2\lambda}\sum_{i=1}^{m} \|w_{i}\|^{2}+C_{1}m\left(\|\varphi\|_{L^{1}}+\frac{J^{2}}{2\lambda}|\Omega| \right)\\ \leq&\,-\frac{1}{2}C_{1}\lambda\sum_{i=1}^{m}\int_{ \Omega}u_{i}^{4}(t,x)\,dx+\frac{C_{1}\sigma^{2}}{2\lambda}\sum_{i=1}^{m}\|w_{ i}\|^{2}\\ &\,+C_{1}m\left(\|\varphi\|\|\Omega|^{1/2}+\frac{1}{2\lambda}(( \lambda+k)^{2}+J^{2})|\Omega|\right),\end{split} \tag{3.1}\]

where the Gauss divergence theorem and the Assumption (1.4) are used. In the last step above, we notice that \((\lambda+k)u_{i}^{2}\leq\frac{(\lambda+k)^{2}}{2\lambda}+\frac{\lambda}{2}u_{i} ^{4}\). Then summing up the \(L^{2}\) inner-products of the \(w_{i}\)-equations with \(w_{i}(t,x)\) for \(1\leq i\leq m\), by Young's inequality (2.6), we get

\[\begin{split}&\frac{1}{2}\frac{d}{dt}\sum_{i=1}^{m}\|w_{i}\|^{2}= \sum_{i=1}^{m}\int_{\Omega}(au_{i}w_{i}+cw_{i}-bw_{i}^{2})\,dx\\ \leq&\,\sum_{i=1}^{m}\int_{\Omega}\left[\left(\frac{ a^{2}}{b}u_{i}^{2}+\frac{1}{4}b\,w_{i}^{2}\right)+\left(\frac{c^{2}}{b}+\frac{1}{4 }b\,w_{i}^{2}\right)-b\,w_{i}^{2}\right]dx\\ =&\,\sum_{i=1}^{m}\frac{a^{2}}{b}\int_{\Omega}u_{i}^{ 2}(t,x)\,dx-\frac{b}{2}\sum_{i=1}^{m}\|w_{i}\|^{2}+\frac{mc^{2}}{b}|\Omega|. \end{split} \tag{3.2}\]

Next, we sum up the \(L^{2}\) inner-products of the \(\rho_{i}\)-equations with \(\rho_{i}(t,x),1\leq i\leq m\), to obtain

\[\frac{1}{2}\frac{d}{dt}\sum_{i=1}^{m}\|\rho_{i}\|^{2}=\sum_{i=1}^{m}\int_{ \Omega}(qu_{i}\rho_{i}-r\rho_{i}^{2})\,dx\leq\frac{q^{2}}{2r}\sum_{i=1}^{m}u_{ i}^{2}(t,x)\,dx-\frac{r}{2}\sum_{i=1}^{m}\|\rho_{i}\|^{2}. \tag{3.3}\]
Now add the above three inequalities. We come up with

\[\begin{split}&\frac{1}{2}\frac{d}{dt}\sum_{i=1}^{m}\left(C_{1}\|u_{i} \|^{2}+\|w_{i}\|^{2}+\|\rho_{i}\|^{2}\right)+C_{1}\eta\sum_{i=1}^{m}\|\nabla u_ {i}\|^{2}+C_{1}P\sum_{i=1}^{m}\sum_{j=1}^{m}\int_{\Omega}(u_{i}-u_{j})^{2}dx\\ \leq&-\sum_{i=1}^{m}\int_{\Omega}\left(\frac{1}{2}C_ {1}\lambda u_{i}^{4}-\left(\frac{a^{2}}{b}+\frac{q^{2}}{2r}\right)u_{i}^{2} \right)dx+\sum_{i=1}^{m}\left(\frac{C_{1}\sigma^{2}}{2\lambda}-\frac{b}{2} \right)\|w_{i}\|^{2}-\frac{r}{2}\sum_{i=1}^{m}\|\rho_{i}\|^{2}\\ &+C_{1}m\left(\|\varphi\|\|\Omega|^{1/2}+\frac{1}{2\lambda}(( \lambda+k)^{2}+J^{2})|\Omega|\right)+\frac{mc^{2}}{b}|\Omega|,\quad t\in I_{ max}=[0,T_{max}),\end{split} \tag{3.4}\]

where \(I_{max}\) is the maximal existence interval of a weak solution. Now we can choose the scaling constant to be

\[C_{1}=\frac{b\lambda}{2\sigma^{2}}\quad\text{so that}\quad\frac{C_{1}\sigma^{2} }{2\lambda}-\frac{b}{2}=-\frac{b}{4}\,. \tag{3.5}\]

With this choice, from (3.4) it follows that

\[\begin{split}&\frac{1}{2}\frac{d}{dt}\sum_{i=1}^{m}\left(C_{1}\|u_{i} \|^{2}+\|w_{i}\|^{2}+\|\rho_{i}\|^{2}\right)+C_{1}\eta\sum_{i=1}^{m}\|\nabla u _{i}\|^{2}+\frac{b}{4}\sum_{i=1}^{m}\|w_{i}\|^{2}+\frac{r}{2}\sum_{i=1}^{m}\| \rho_{i}\|^{2}\\ &+\sum_{i=1}^{m}\int_{\Omega}\left(\frac{1}{2}C_{1}\lambda u_{i}^ {4}-\left(\frac{a^{2}}{b}+\frac{q^{2}}{2r}\right)u_{i}^{2}\right)dx+C_{1}P\sum _{i=1}^{m}\sum_{j=1}^{m}\int_{\Omega}(u_{i}-u_{j})^{2}\,dx\\ \leq& C_{1}m\|\varphi\|^{2}+m\left(C_{1}+\frac{C_{1}} {2\lambda}((\lambda+k)^{2}+J^{2})+\frac{c^{2}}{b}\right)|\Omega|,\quad t\in I _{max}=[0,T_{max}).\end{split} \tag{3.6}\]

By completing square and (3.5), we have

\[\begin{split}&\sum_{i=1}^{m}\int_{\Omega}\left[\frac{1}{2}C_{1} \lambda u_{i}^{4}-\left(\frac{a^{2}}{b}+\frac{q^{2}}{2r}\right)u_{i}^{2}\right] dx\\ =&\sum_{i=1}^{m}\left(\frac{a^{2}}{b}+\frac{q^{2}}{2r }\right)\|u_{i}\|^{2}+\sum_{i=1}^{m}\int_{\Omega}\left[\frac{b\lambda^{2}}{4 \sigma^{2}}u_{i}^{4}-2\left(\frac{a^{2}}{b}+\frac{q^{2}}{2r}\right)u_{i}^{2} \right]dx\\ =&\sum_{i=1}^{m}\left(\frac{a^{2}}{b}+\frac{q^{2}}{2r }\right)\|u_{i}\|^{2}+\sum_{i=1}^{m}\int_{\Omega}\left[\frac{\sqrt{b}\lambda}{ 2\sigma}u_{i}^{2}-\frac{2\sigma}{\sqrt{b}\lambda}\left(\frac{a^{2}}{b}+\frac{ q^{2}}{2r}\right)\right]^{2}dx\\ &-\frac{4m\sigma^{2}}{b\lambda^{2}}\left[\frac{a^{2}}{b}+\frac{q ^{2}}{2r}\right]^{2}|\Omega|\geq\sum_{i=1}^{m}\left(\frac{a^{2}}{b}+\frac{q^{ 2}}{2r}\right)\|u_{i}\|^{2}-\frac{4m\sigma^{2}}{b\lambda^{2}}\left[\frac{a^{2} }{b}+\frac{q^{2}}{2r}\right]^{2}|\Omega|.\end{split} \tag{3.7}\]
Substitute (3.7) in (3.6). It yields the inequality

\[\begin{split}&\frac{1}{2}\frac{d}{dt}\sum_{i=1}^{m}\left(C_{1}\|u_{i} \|^{2}+\|w_{i}\|^{2}+\|\rho_{i}\|^{2}\right)+C_{1}\eta\sum_{i=1}^{m}\|\nabla u_ {i}\|^{2}+C_{1}P\sum_{i=1}^{m}\sum_{j=1}^{m}\int_{\Omega}(u_{i}-u_{j})^{2}\,dx \\ &+\sum_{i=1}^{m}\left(\frac{a^{2}}{b}+\frac{q^{2}}{2r}\right)\|u_ {i}\|^{2}+\frac{b}{4}\sum_{i=1}^{m}\|w_{i}\|^{2}+\frac{r}{2}\sum_{i=1}^{m}\| \rho_{i}\|^{2}\\ \leq&\,C_{1}m\|\varphi\|^{2}+m\left[C_{1}+\frac{C_{1 }}{2\lambda}((\lambda+k)^{2}+J^{2})+\frac{c^{2}}{b}+\frac{4\sigma^{2}}{b \lambda^{2}}\left[\frac{a^{2}}{b}+\frac{q^{2}}{2r}\right]^{2}\right]|\Omega|, \quad t\in I_{max}.\end{split} \tag{3.8}\]

Denote by

\[C_{2}=C_{1}+\frac{C_{1}}{2\lambda}((\lambda+k)^{2}+J^{2})+\frac{c^{2}}{b}+ \frac{4\sigma^{2}}{b\lambda^{2}}\left[\frac{a^{2}}{b}+\frac{q^{2}}{2r}\right] ^{2}. \tag{3.9}\]

Then (3.8) gives rise to the Gronwall-type differential inequality

\[\begin{split}&\frac{d}{dt}\sum_{i=1}^{m}\left[C_{1}\|u_{i}\|^{2}+\|w_{i} \|^{2}+\|\rho_{i}\|^{2}\right]+\mu\sum_{i=1}^{m}\left[C_{1}\|u_{i}\|^{2}+\|w_{i }\|^{2}+\|\rho_{i}\|^{2}\right]\\ \leq&\,\frac{d}{dt}\sum_{i=1}^{m}\left[C_{1}\|u_{i}\| ^{2}+\|w_{i}\|^{2}+\|\rho_{i}\|^{2}\right]+2\sum_{i=1}^{m}\left[\left(\frac{a^ {2}}{b}+\frac{q^{2}}{2r}\right)\|u_{i}\|^{2}+\frac{b}{4}\|w_{i}\|^{2}+\frac{r} {2}\|\rho_{i}\|^{2}\right]\\ \leq&\,2C_{1}m\|\varphi\|^{2}+2C_{2}m|\Omega|,\quad \text{for}\;\;t\in I_{max}=[0,T_{max}),\end{split} \tag{3.10}\]

where

\[\mu=\min\,\left\{\frac{2a^{2}}{b}+\frac{q^{2}}{r},\;\frac{b}{2},\;r\right\}.\]

We can solve the differential inequality (3.10) to obtain the following bounding estimate of all the weak solutions on the maximal existence time interval \(I_{max}\),

\[\begin{split}&\|g(t,g^{0})\|^{2}=\sum_{i=1}^{m}\|g_{i}(t,g_{i}^{0}) \|^{2}=\sum_{i=1}^{m}\left(\|u_{i}(t)\|^{2}+\|w_{i}(t)\|^{2}+\|\rho_{i}(t)\|^{ 2}\right)\\ \leq&\,\frac{\max\{C_{1},1\}}{\min\{C_{1},1\}}e^{- \mu\,t}\|g^{0}\|^{2}+\frac{2m}{\mu\min\{C_{1},1\}}\left(C_{1}\|\varphi\|^{2}+C_ {2}|\Omega|\right),\quad t\in[0,\infty).\end{split} \tag{3.11}\]

Here it is shown that \(I_{max}=[0,\infty)\) for every weak solution \(g(t,g^{0})\) because it will never blow up at any finite time. Therefore, for any initial state \(g^{0}=(g_{1}^{0},\cdots,g_{m}^{0})\in E\), there exists a unique global weak solution in time \(t\in[0,\infty)\) of the initial value problem (2.1) for this memristive neural network model (1.1) in the space \(E\)
Based on the global existence of weak solutions shown in Theorem 3.1, we define the solution semiflow \(\{S(t):E\to E\}_{t\geq 0}\) of the memristive and diffusive FitzHugh-Nagumo equations (1.1) to be

\[S(t):g^{0}\longmapsto g(t;g^{0})=\operatorname{col}\left(u_{i}(t,\cdot),w_{i}(t,\cdot),\rho_{i}(t,\cdot):1\leq i\leq m\right),\quad t\geq 0.\]

We call this semiflow \(\{S(t)\}_{t\geq 0}\) the _memristive FitzHugh-Nagumo neural network semiflow_ generated by the neural network model equations (1.1).

The next theorem shows that the memristive FitzHugh-Nagumo neural network semiflow \(\{S(t)\}_{t\geq 0}\) is a dissipative dynamical system in the state space \(E\).

**Theorem 3.2**.: _There exists a bounded absorbing set for the memristive FitzHugh-Nagumo neural network semiflow \(\{S(t)\}_{t\geq 0}\) in the state space \(E\), which is the bounded ball_

\[B^{*}=\{h\in E:\|h\|^{2}\leq K\}. \tag{3.12}\]

_Here the constant_

\[K=1+\frac{2m}{\mu\min\{C_{1},1\}}\left(C_{1}\|\varphi\|^{2}+C_{2}|\Omega| \right), \tag{3.13}\]

_where the constants \(C_{1}\) and \(C_{2}\) are given in (3.5) and (3.9)._

Proof.: This is the consequence of the global uniform estimate (3.11) shown in the proof of Theorem 3.1, which implies that

\[\limsup_{t\to\infty}\|g(t;g^{0})\|^{2}=\limsup_{t\to\infty}\,\sum_{i=1}^{m}\|g _{i}(t;g_{i}^{0})\|^{2}<K \tag{3.14}\]

for all weak solutions of (2.1) with any initial data \(g^{0}\) in \(E\). Moreover, for any given bounded set \(B=\{h\in E:\|h\|^{2}\leq L\}\) in \(E\), there exists a finite time

\[T_{B}=\frac{1}{\mu}\log^{+}\left(L\,\frac{\max\{C_{1},1\}}{\min\{C_{1},1\}}\right)\]

such that all the solution trajectories started at the initial time \(t=0\) from the set \(B\) will permanently enter the bounded ball \(B^{*}\) shown in (3.12) for \(t>T_{B}\). Therefore, the bounded ball \(B^{*}\) is an absorbing set in \(E\) for the semiflow \(\{S(t)\}_{t\geq 0}\) so that this memristive FitzHugh-Nagumo neural network semiflow is dissipative. 

We shall further prove an ultimate uniform bound of the membrane potential functions \(\{u_{i}(t,x):1\leq i\leq m\}\) for all the weak solutions in the higher-order integrable space \(L^{4}(\Omega)\), which paves the way for the attempt to achieve the memristive neural network synchronization in the next section.
**Theorem 3.3**.: _There exists a constant \(Q>0\) such that for any initial data \(g^{0}\in E\), the \(u_{i}(t)\) components, \(1\leq i\leq m\), of the weak solution \(g(t;g^{0})=(g_{1}(t),\cdots,g_{m}(t))\) of the initial value problem (2.1) for the memristive FitzHugh-Nagumo neural network \(\mathcal{NW}\) satisfies the absorbing property in the space \(L^{4}(\Omega)\),_

\[\limsup_{t\to\infty}\,\sum_{i=1}^{m}\|u_{i}(t)\|_{L^{4}}^{4}<1+Q. \tag{3.15}\]

Proof.: Take the \(L^{2}\) inner-product of the \(u_{i}\)-equation in (1.1) with \(u_{i}^{3}(t),1\leq i\leq m\), and sum them up. By the boundary condition (1.2) and Assumption (1.4), we get

\[\begin{split}&\frac{1}{4}\,\frac{d}{dt}\sum_{i=1}^{m}\|u_{i}(t)\|_{L^ {4}}^{4}+3\eta\sum_{i=1}^{m}\|u_{i}\nabla u_{i}\|_{L^{2}}^{2}\\ &+P\sum_{i=1}^{m}\sum_{j=1}^{m}\int_{\Omega}(u_{i}-u_{j})^{2}(u_ {i}^{2}+u_{i}u_{j}+u_{j}^{2})\,dx\\ &=\sum_{i=1}^{m}\int_{\Omega}(f(u_{i},x)u_{i}^{3}-\sigma u_{i}^{ 3}w_{i}+Ju_{i}^{3}-k\tanh(\rho_{i})u_{i}^{4})\,dx\\ &\leq\sum_{i=1}^{m}\int_{\Omega}(-\lambda u_{i}^{6}+u_{i}^{2} \varphi(x)-\sigma u_{i}^{3}w_{i}+Ju_{i}^{3}+ku_{i}^{4})\,dx,\quad t>0.\end{split} \tag{3.16}\]

By Cauchy inequality, it is seen that

\[u_{i}^{2}\varphi(x)-\sigma u_{i}^{3}w_{i}+Ju_{i}^{3}\leq\frac{1}{6}\lambda u_ {i}^{4}+\frac{1}{3}\lambda u_{i}^{6}+\frac{6}{\lambda}\left(\varphi^{2}(x)+ \sigma^{2}|w_{i}(t,x)|^{2}+J^{2}\right). \tag{3.17}\]

Using Young's inequality (2.6),

\[ku_{i}^{4}\leq\frac{1}{3}\left(\frac{16\,k^{3}}{\lambda^{2}}\right)+\frac{2}{ 3}\left(\frac{\lambda}{4}u_{i}^{6}\right)\leq\frac{6k^{3}}{\lambda^{2}}+\frac {1}{6}\lambda u_{i}^{6}. \tag{3.18}\]

Note that \(\sum_{i=1}^{m}\sum_{j=1}^{m}(u_{i}-u_{j})^{2}(u_{i}^{2}+u_{i}u_{j}+u_{j}^{2})\geq 0\) always holds and

\[u_{i}^{4}\leq\frac{1}{3}+\frac{2}{3}u_{i}^{6}\leq 1+u_{i}^{6},\]

so that

\[-\frac{1}{2}{u_{i}}^{6}\leq\frac{1}{2}-\frac{1}{2}{u_{i}}^{4}, \tag{3.19}\]

From (3.16) wherein we successively use the above inequalities (3.17), (3.18), (3.19) and (3.14), it follows that
\[\frac{1}{4}\,\frac{d}{dt}\sum_{i=1}^{m}\|u_{i}(t)\|_{L^{4}}^{4}+3\eta \sum_{i=1}^{m}\|u_{i}\nabla u_{i}\|^{2}\] \[\leq\sum_{i=1}^{m}\left(\lambda\int_{\Omega}\left(\frac{1}{6}u_{i} ^{4}-\frac{1}{2}u_{i}^{6}\right)dx+\frac{6\sigma^{2}}{\lambda}\|w_{i}(t)\|^{2} \right)+\frac{6m}{\lambda}\left(\|\varphi\|^{2}+\left(J^{2}+\frac{k^{3}}{ \lambda}\right)|\Omega|\right)\] \[\leq\sum_{i=1}^{m}\left(-\frac{\lambda}{3}\int_{\Omega}u_{i}^{4} \,dx+\frac{6\sigma^{2}}{\lambda}\|w_{i}(t)\|^{2}\right)+\frac{m\lambda}{2}| \Omega|+\frac{6m}{\lambda}\left(\|\varphi\|^{2}+\left(J^{2}+\frac{k^{3}}{ \lambda}\right)|\Omega|\right)\] \[\leq-\frac{\lambda}{3}\sum_{i=1}^{m}\|u_{i}(t)\|_{L^{4}}^{4}+ \frac{6\sigma^{2}}{\lambda}K+m\left[\frac{6}{\lambda}\|\varphi\|^{2}+\left( \frac{\lambda}{2}+\frac{6}{\lambda}J^{2}+\frac{6}{\lambda^{2}}k^{3}\right)| \Omega|\right],\;\;t>0. \tag{3.20}\]

Consequently, with the non-negative gradient term removed, (3.20) shows that

\[\frac{d}{dt}\sum_{i=1}^{m}\|u_{i}(t)\|_{L^{4}}^{4}+\frac{4\lambda }{3}\sum_{i=1}^{m}\|u_{i}(t)\|_{L^{4}}^{4} \tag{3.21}\] \[\leq\frac{24\sigma^{2}}{\lambda}K+m\left[\frac{24}{\lambda}\| \varphi\|^{2}+\left(2\lambda+\frac{24}{\lambda}J^{2}+\frac{24}{\lambda^{2}}k^ {3}\right)|\Omega|\right],\quad t>0.\]

By the parabolic regularity stated in Proposition 2.2, for any weak solution \(g(t;g^{0})\) one has \(u_{i}(1)\in H^{1}(\Omega)\subset L^{4}(\Omega)\) for \(1\leq i\leq m\). Then the second statement in Proposition 2.2 shows that any weak solution has the regularity

\[\sum_{i=1}^{m}u_{i}(t)\in C([1,\infty),H^{1}(\Omega))\subset C([1,\infty),L^{ 4}(\Omega)).\]

Apply the Gronwall inequality to (3.21). It results in the bounding estimate of all the \(u_{i}(t)\) components in the space \(L^{4}(\Omega)\) as follows:

\[\sum_{i=1}^{m}\|u_{i}(t)\|_{L^{4}}^{4}\leq e^{-\frac{4\lambda}{3}(t-1)}\sum_{ i=1}^{m}\|u_{i}(1)\|_{L^{4}}^{4}+Q,\quad\text{for}\;\,t\geq 0, \tag{3.22}\]

where the constant \(Q\) is independent of any initial data and given by

\[Q=\frac{18\,\sigma^{2}}{\lambda^{2}}K+m\left[\frac{18}{\lambda^{2}}\|\varphi \|^{2}+\left(\frac{3}{2}+\frac{18}{\lambda^{2}}J^{2}+\frac{18}{\lambda^{3}}k^ {3}\right)|\Omega|\right]. \tag{3.23}\]

Therefore, the claim (3.15) of this theorem is proved.
## 4. **Synchronization of Memristive FitzHugh-Nagumo Neural Networks**

In this section, we shall prove the main result on the exponential synchronization of the memristive FitzHugh-Nagumo neural networks described by (1.1) in the state space \(E\). This result provides a sufficient quantitative threshold condition for the neuron coupling strength to reach the neural network synchronization.

**Definition 4.1**.: For a model evolutionary equation of a neural network called NW such as (2.1) formulated from the memristive and diffusive FitzHugh-Nagumo equations (1.1), we define the asynchronous degree of this neural network in a state space (as a Banach space) \(Z\) to be

\[deg_{s}\left(\mathrm{NW}\right)=\sum_{1\,\leq i\,<j\,\leq\,m}\left\{\sup_{g_{ i}^{0},\,g_{j}^{0}\,\in\,Z}\,\left\{\limsup_{t\to\infty}\,\|g_{i}(t;g_{i}^{0})-g_{j}(t; g_{j}^{0})\|_{Z}\right\}\right\}\]

where \(g_{i}(t)\) and \(g_{j}(t)\) are any two solutions of the model equation with the initial states \(g_{i}^{0}\) and \(g_{j}^{0}\) for two neurons \(\mathcal{N}_{i}\) and \(\mathcal{N}_{j}\) in the network. The neural network is said to be asymptotically synchronized if

\[deg_{s}\left(\mathrm{NW}\right)=0.\]

If the asymptotic convergence to zero of the difference norm above for any two neurons in the network admits a uniform exponential rate, then the neural network is called exponentially synchronized.

Introduce the neuron difference functions: For \(i,j=1,\cdots,m\), we define

\[U_{ij}(t,x)=u_{i}(t,x)-u_{j}(t,x),\,W_{ij}(t,x)=w_{i}(t,x)-w_{j}(t,x),\,R_{ij} (t,x)=\rho_{i}(t,x)-\rho_{j}(t,x).\]

Given any initial state \(g^{0}=\mathrm{col}\left(g_{1}^{0},\cdots,g_{m}^{0}\right)\) in the space \(E\), the difference between any two solutions of (2.1) associated with two neurons \(\mathcal{N}_{i}\) and \(\mathcal{N}_{j}\) in the network is what we consider:

\[g_{i}(t,g_{i}^{0})-g_{j}(t,g_{j}^{0})=\mathrm{col}\left(U_{ij}(t,\cdot),W_{ij }(t,\cdot),R_{ij}(t,\cdot)\right),\quad t\geq 0.\]

By subtraction of the three governing equations for the \(j\)-th neuron from the corresponding equations for the \(i\)-th neuron in (1.1), we obtain the following differencing FitzHugh-Nagumo equations. For \(i,j=1,\cdots,m\),

\[\frac{\partial U}{\partial t}=\eta\Delta U+f(u_{i},x)-f(u_{j},x)-\sigma W-k( \tanh(\rho_{i})u_{i}-\tanh(\rho_{j})u_{j})-mPU,\]

\[\frac{\partial W}{\partial t}=aU-bW, \tag{4.1}\]

\[\frac{\partial R}{\partial t}=qU-rR.\]
Here we can simply write \(U(t,x)=U_{ij}(t,x),W(t,x)=W_{ij}(t,x),R(t,x)=R_{ij}(t,x)\) and further \(U(t)=U(t,\cdot),W(t)=W(t,\cdot),R(t)=R(t,\cdot)\) for notational convenience.

The following exponential synchronization theorem is the main result of this paper.

**Theorem 4.2**.: _For the memristive FitzHugh-Nagumo neural network \(\mathcal{NW}\) with the model (1.1)-(1.4), If the following threshold condition is satisfied by the coupling strength coefficient \(P\),_

\[P>\Gamma. \tag{4.2}\]

_Here the threshold \(\Gamma\) is defined to be_

\[\Gamma=\frac{1}{m}\left(\beta+k+\frac{1}{2b}|a-\sigma|^{2}+\frac{q^{2}}{r}+ \frac{C^{*4}\,k^{8}(1+Q)^{2}}{\eta^{3}\,r^{4}}\right),\]

_where the constant \(Q\) is given in (3.23) and the constant \(C^{*}\) is the coefficient in the Gagliardo-Nirenberg inequality (4.8), then the neural network \(\mathcal{NW}\) is exponentially synchronized in the state space \(E\) at a uniform exponential rate \(\alpha(P)\):_

\[\alpha(P)=\min\left\{b,\,r,\,\left[2mP-2\left(\beta+k+\frac{1}{2b}|a-\sigma|^{ 2}+\frac{q^{2}}{r}+\frac{C^{*4}\,k^{8}(1+Q)^{2}}{\eta^{3}\,r^{4}}\right) \right]\right\}. \tag{4.3}\]

Proof.: The proof will go through two steps.

Step 1. Take the \(L^{2}\) inner-products of the first equation in (4.1) with \(U(t)\), the second equation in (4.1) with \(W(t)\), and the third equation in (4.1) with \(R(t)\). Then sum them altogether and use the Assumption (1.4) to get

\[\begin{split}&\frac{1}{2}\frac{d}{dt}(\|U(t)\|^{2}+\|W(t)\|^{2}+\|R(t) \|^{2})+\eta\|\nabla U(t)\|^{2}+b\,\|W(t)\|^{2}+r\|R(t)\|^{2}\\ =&\int_{\Omega}(f(u_{i},x)-f(u_{j},x))U\,dx+\int_{ \Omega}[(a-\sigma)UW+qUR]\,dx\\ &-\int_{\Omega}k(\tanh(\rho_{i})u_{i}-\tanh(\rho_{j})u_{j})U]\, dx-mP\|U(t)\|^{2}\\ \leq&\int_{\Omega}\frac{\partial f}{\partial s}\left( \ell u_{i}+(1-\ell)u_{j},x\right)U^{2}\,dx+\int_{\Omega}[(a-\sigma)UW+qUR]\,dx \\ &-k\int_{\Omega}\left[\operatorname{sech}^{2}(\xi\rho_{i}+(1-\xi )\rho_{j})R\,u_{i}U+\tanh(\rho_{j})U^{2}\right]dx-mP\|U(t)\|^{2}\\ \leq&\,\beta\|U\|^{2}+\int_{\Omega}[(a-\sigma)UW+qUR ]\,dx-k\int_{\Omega}Ru_{j}U\,dx+(k-mP)\|U\|^{2},\end{split} \tag{4.4}\]

where the properties \(|\tanh(\rho_{j})|\leq 1\) and \(\operatorname{sech}^{2}(\xi\rho_{i}+(1-\xi)\rho_{j})\leq 1\) for the hyperbolic functions and \(\xi,\ell\in[0,1]\) are used.

Next we have to treat the two integral terms on the right-hand side of the differential inequality (4.4). By the Young's inequality (2.6), we have
\[\int_{\Omega}[(a-\sigma)UW+qUR]\,dx\] \[\leq\int_{\Omega}\left[\frac{b}{2}W^{2}(t,x)+\frac{1}{2b}|a-\sigma|^ {2}U^{2}(t,x)\right]dx+\int_{\Omega}\left[\frac{r}{4}R^{2}(t,x)+\frac{q^{2}}{r}U ^{2}(t,x)\right]dx \tag{4.5}\] \[= \,\frac{b}{2}\|W(t)\|^{2}+\frac{r}{4}\|R(t)\|^{2}+\left(\frac{1}{ 2b}|a-\sigma|^{2}+\frac{q^{2}}{r}\right)\|U(t)\|^{2},\quad t>0.\]

For the last integral term in (4.4), by the Holder inequality, we get

\[-k\int_{\Omega}Ru_{j}U\,dx\leq k\int_{\Omega}\left(\frac{r}{4k}R^{2}(t,x)+\frac{k}{r}u_{j}^{2}(t,x)U ^{2}(t,x)\right)dx\] \[\leq \,\frac{r}{4}\|R(t)\|^{2}+\frac{k^{2}}{r}\left[\int_{\Omega}u_{j }^{4}(t,x)\,dx\right]^{1/2}\left[\int_{\Omega}U^{4}(t,x)\,dx\right]^{1/2} \tag{4.6}\] \[= \,\frac{r}{4}\|R(t)\|^{2}+\frac{k^{2}}{r}\|u_{j}(t)\|_{L^{4}}^{2 }\|U(t)\|_{L^{4}}^{2},\quad t>0.\]

Substitute the term estimates (4.5) and (4.6) into the differential inequality (4.4), we obtain

\[\frac{1}{2}\frac{d}{dt}(\|U(t)\|^{2}+\|W(t)\|^{2}+\|R(t)\|^{2})+ \eta\|\nabla U(t)\|^{2}+\frac{b}{2}\|W(t)\|^{2}+\frac{r}{2}\|R(t)\|^{2} \tag{4.7}\] \[\leq \,\left(\beta+k+\frac{1}{2b}|a-\sigma|^{2}+\frac{q^{2}}{r}-mP \right)\|U(t)\|^{2}+\frac{k^{2}}{r}\|u_{j}(t)\|_{L^{4}}^{2}\|U(t)\|_{L^{4}}^{2 },\,\,\,t>0.\]

Step 2. The key challenge here is to handle the last term on the right-hand side of the estimate inequality (4.7). We shall exploit the sharp technique of Gagliardo-Nirenberg interpolation inequalities [27, Theorem B.3]. In view of the Sobolev embedding

\[H^{1}(\Omega)\subset L^{4}(\Omega)\subset L^{2}(\Omega),\]

one has

\[\|U(t)\|_{L^{4}}^{2}\leq C^{*}\|\nabla U(t)\|^{2\theta}\|U(t)\|^{2(1-\theta)} \tag{4.8}\]

where the coefficient \(C^{*}(\Omega)>0\) only depends on the spatial domain \(\Omega\), and the interpolation index \(\theta=3/4\) is determined by

\[-\frac{3}{4}=\theta\left(1-\frac{3}{2}\right)-\left(1-\theta\right)\frac{3}{2}\,.\]

Hence (4.8) shows that

\[\|U(t)\|_{L^{4}}^{2}\leq C^{*}\|\nabla U(t)\|^{3/2}\|U(t)\|^{1/2}. \tag{4.9}\]
According to Theorem 3.3 and (3.15), \(\limsup_{t\to\infty}\sum_{i=1}^{m}\|u_{i}(t)\|_{L^{4}}^{4}<1+Q\), so that there exists a finite time \(T(g^{0})\geq 0\) such that for all \(1\leq i\leq m\),

\[\|u_{i}(t)\|_{L^{4}}^{2}<(1+Q)^{1/2},\quad\text{for all}\;\;t>T(g^{0}).\]

Therefore, for any given initial state \(g^{0}\in E\), by (4.9) and using Young's inequality (2.6), we can estimate

\[\begin{split}&\frac{k^{2}}{r}\|u_{j}(t)\|_{L^{4}}^{2}\|U(t)\|_{L^{4 }}^{2}\leq\frac{k^{2}}{r}(1+Q)^{1/2}\,\|U(t)\|_{L^{4}}^{2}\\ \leq& C^{*}\|\nabla U(t)\|^{3/2}\left[\frac{k^{2}}{ r}(1+Q)^{1/2}\,\|U(t)\|^{1/2}\right]\\ \leq&\eta\|\nabla U(t)\|^{(3/2)\times(4/3)}+\frac{1} {\eta^{3}}\left[\frac{C^{*}k^{2}}{r}(1+Q)^{1/2}\,\|U(t)\|^{1/2}\right]^{4}\\ =&\eta\|\nabla U(t)\|^{2}+\frac{C^{*4}\,k^{8}(1+Q)^ {2}}{\eta^{3}\,r^{4}}\|U(t)\|^{2},\quad t>T(g^{0}).\end{split} \tag{4.10}\]

Substitute (4.10) in (4.7). After cancellation of the gradient terms \(\eta\|\nabla U(t)\|^{2}\) on both sides of the inequality, it follows that for \(t>T(g^{0})\),

\[\begin{split}&\frac{1}{2}\frac{d}{dt}(\|U(t)\|^{2}+\|W(t)\|^{2}+\|R(t)\|^{2})+ \frac{b}{2}\|W(t)\|^{2}+\frac{r}{2}\|R(t)\|^{2}\\ +&\left[mP-\left(\beta+k+\frac{1}{2b}|a-\sigma|^{2}+ \frac{q^{2}}{r}+\frac{C^{*4}\,k^{8}(1+Q)^{2}}{\eta^{3}\,r^{4}}\right)\right] \|U(t)\|^{2}\leq 0.\end{split} \tag{4.11}\]

Therefore, the following Gronwall-type inequality holds:

\[\begin{split}&\frac{d}{dt}(\|U(t)\|^{2}+\|W(t)\|^{2}+\|R(t)\|^{2})+ \alpha(P)(\|U(t)\|^{2}+\|W(t)\|^{2}+\|R(t)\|^{2})\\ \leq&\frac{d}{dt}(\|U(t)\|^{2}+\|W(t)\|^{2}+\|R(t)\|^ {2})+(b\|W(t)\|^{2}+r\|R(t)\|^{2})\\ +&\,2\left[mP-\left(\beta+k+\frac{1}{2b}|a-\sigma|^{2}+ \frac{q^{2}}{r}+\frac{C^{*4}\,k^{8}(1+Q)^{2}}{\eta^{3}\,r^{4}}\right)\right] \|U(t)\|^{2}\leq 0,\end{split} \tag{4.12}\]

for \(t>T(g^{0})\), where the constant

\[\alpha(P)=\min\left\{b,\,r,\,\left[2mP-2\left(\beta+k+\frac{1}{2b}|a-\sigma|^{ 2}+\frac{q^{2}}{r}+\frac{C^{*4}\,k^{8}(1+Q)^{2}}{\eta^{3}\,r^{4}}\right) \right]\right\}\]

as shown in (4.3).

Under the threshold condition (4.2) of this theorem, solving this linear Gronwall inequality (4.12) directly shows the exponential synchronization result: For any initial state \(g^{0}\in E\) and any two neurons \(\mathcal{N}_{i}\) and \(\mathcal{N}_{j}\) in the memristive FitzHugh-Nagumo neural network (1.1), their difference function \(g_{i}(t;g_{i}^{0})-g_{j}(t;g_{j}^{0})\) converges to zero
in the state space \(E\) exponentially at a uniform rate \(\alpha(P)\) shown in (4.3). Namely, for any \(1\leq i<j\leq m\),

\[\begin{split}\|g_{i}(t)-g_{j}(t)\|_{E}^{2}&=\|U_{ij} (t)\|^{2}+\|W_{ij}(t)\|^{2}+\|R_{ij}(t)\|^{2}\\ &\leq e^{-\alpha(P)\,t}\left\|g_{i}^{0}-g_{j}^{0}\right\|^{2} \to 0,\ \ \text{as}\,\ t\to\infty.\end{split} \tag{4.13}\]

Hence it is proved that

\[deg_{s}(\mathcal{NW})=\sum_{1\,\leq\,i\,<\,j\,\leq\,m}\left\{\sup_{g^{0}\,\in \,E}\,\left\{\limsup_{t\to\infty}\|g_{i}(t)-g_{j}(t)\|_{E}^{2}\right\}\right\} =0. \tag{4.14}\]

Thus the exponential synchronization of the memristive and diffusive Hindmarsh-Rose neural network \(\mathcal{NW}\) in the space \(E\) is proved. 

## 5. **Example and Numerical Simulation**

In this section, we test some numerical experiments to verify and illustrate the obtained theoretical result on synchronization stated in Theorem 4.2.

We numerically solve the memristive FitzHugh-Nagumo neural network \(\mathcal{NW}\) with the model (1.1)-(1.4) in a two-dimensional square domain. We use the finite difference method for the numerical scheme and programmed in Python.

Choosing \(f(s)=s(s-1)(1-s)\), we consider the following selection of parameters:

\[\begin{split}& m=4,\ \eta=10;\ \sigma=0.01;\ J=0.5;\ P=1.45;\\ & a=0.35;\ b=0.35;\ c=0.7;\ q=0.35;\ r=10.\end{split}\]

Make the time-step to be \(0.00025\)s and spatial-step to be \(1\) on a \(32*32\) membrane. We compute the \(L^{2}\) norm of neuron potential \(u_{i}\), the recovering variable \(w_{i}\), the memductance \(\rho_{i}\), and also the vector solutions \(g_{i}\) from (4.13) in the energy space \(E\) as showing in Figure 1 to Figure 4.

In Figure 1 to Figure 3, with a comparison between the beginning stages and results after \(10000\) iterations, one can observe the synchronization tendency of the three characterizing variables \((u_{i},w_{i},\rho_{i})\) among the neurons in the simulated mimristive neural network. From Figure 4, we observe that the differences among \(\|g_{i}\|\) tend to \(0\).

We can get the the following constants that used in Theorem 4.2 based on our selection of parameters.

\[\begin{split}&\lambda=0.25,\quad\phi(x)=4,\quad\beta=\frac{4}{3},\\ & C_{1}=437.5\quad C_{2}=876.4\quad\mu=0.175\quad K=41345645.6 \quad 1+Q=1220899.6\quad C^{*}=0.4\\ &\quad P=1.45>\Gamma=0.45,\quad\alpha=0.35.\end{split}\]
Remark: The \(\lambda,\ \phi,\ \beta\) from (1.4) are not unique. The constant \(C^{*}\) from Gargliardo-Nirenberg inequality is chosen to be 0.4 based on [2].

Table 1 to Table 3 list the sampled values of the three components \(u_{i},w_{i},\) and \(\rho_{i}\) of the simulated solution \(g_{i}\) at one same point in the domain at \(t=0\) and at the 10000 time-step. It is seen that with a big difference on the initial values, after a certain time, the curves of \(u_{i}\), \(w_{i},\) and \(\rho_{i}\) tend to be close to each other between various neurons.

The synchronization result rigorously proved in this work is illustrated by the example with sample selections of the system parameters and a randomized set of initial data. Our numerical simulation also exhibits that the neuron potentials \(u_{i}\) seem to be synchronized fastest within a limited time, while it takes much longer time to observe the synchronization on the other two variables \(w_{i}\) and \(\rho_{i}\).

Figure 1. The \(L^{2}\) norm of the neurons \(u_{i}\) at the beginning (the upper figure) and after 10000 iterations (the lower figure)
This observation actually enhances the conjecture that adding a nonlinear memristor coupling in the neuron potential equation would accelerate the synchronization for the main variable of neuron membrane potential. On the other hand, it also hints that the neuron potential equation is a nonlinear function of the neuron potential.

Figure 3. The \(L^{2}\) norm of the neurons \(\rho_{i}\) at the beginning (the upper figure) and after 10000 iterations (the lower figure)

Figure 2. The \(L^{2}\) norm of the neurons \(w_{i}\) at the beginning (the upper figure) and after 10000 iterations (the lower figure)
that although the main result Theorem 4.2 confirmed the exponential synchronization has a uniform but may be small convergence rate, each of the three components may have a different synchronization rate, which turns out to be a new open and interesting problem for further research.

\begin{table}
\begin{tabular}{c|c|c}  & Initial Value & At the 10000 time step \\ \hline \(w_{1}\) & 0.038162808368984426 & 1.4968727659565046 \\ \(w_{2}\) & 0.016109788021287225 & 1.4881936913572227 \\ \(w_{3}\) & 0.028702864538319866 & 1.4927071463421802 \\ \(w_{4}\) & 0.04497009552599465 & 1.4999646128986852 \\ \end{tabular}
\end{table}
Table 2. Comparision of the \(w_{i}\) at the point \(x=10,\ y=10\)

Figure 4. The \(L^{2}\) norm of the neurons from the Energy space at the beginning (the upper figure) and after 10000 iterations (the lower figure)

\begin{table}
\begin{tabular}{c|c|c}  & Initial Value & At the 10000 time step \\ \hline \(\rho_{1}\) & 0.02516293441173103 & 0.03032291951665976 \\ \(\rho_{2}\) & 0.040098450204586085 & 0.030324374854761176 \\ \(\rho_{3}\) & 0.015327673754565209 & 0.030322776082447638 \\ \(\rho_{4}\) & 0.03404629085391771 & 0.03032306624583209 \\ \end{tabular}
\end{table}
Table 3. Comparision of the \(\rho_{i}\) at the point \(x=10,\ y=10\)
## 6. **Conclusions**

We summarize the new contribution of results in this paper.

1. We propose a new mathematical model in (1.1) of memristive neural networks in terms of the partly diffusive FitzHugh-Nagumo equations with a nonlinear memristor and linear synaptic coupling in the membrane potential equations. This model as a hybrid system of partial-ordinary differential equations features the full synaptic coupling among all the neurons. In comparison with the extensively studied ODE models, this model is more meaningful to capture the structure of biological neuron cells with long-branch axons in neurodynamics.

2. The nonlinear memristive coupling across the neurons membrane for this work is in the form of \(k\tanh(\rho_{i})u_{i}\), which appeared in quite a few researches as cited in the references. In this paper we proved the exponential synchronization of this model of memristive neural networks. There is another type of coupled memristors in the quadratic form \(k(c+\gamma\rho_{i}+\delta\rho_{i}^{2})u_{i}\), which have also been actively studied most with the ODE models of Hindmarsh-Rose equations on various topics of neuromorphic patterns and chaotic dynamics. But synchronization of the hybrid PDE neural network models such as the diffusive FitzHugh-Nagumo equations with quadratic memristors and linear membrane potential coupling seems still an open problem.

3. In this work we take the analytic approach of global dynamics for the weak solutions to pursue the synchronization investigation. Through the uniform _a priori_ estimates of grouped component solutions with adjustable scaling and maneuvering the integral inequalities, we are able to show the existence of absorbing set in the \(L^{2}\)-energy space of the solution semiflow and the existence of asymptotic ultimate bound in the higher-order integrable space of the key component solutions. This signifies our methodology from dissipative global dynamics to synchronization for such a complex neural network system.

4. The spirit of the entire mathematical proof is to tackle and control the nonlinear memductance-potential effect by the linear network coupling in different integrable spaces. Many steps of sharp analysis including the crucial Gagliardo-Nirenberg interpolation are carried out and cohesively managed.

5. The main result Theorem 4.2 of this paper provides a sufficient threshold condition for achieving the exponential synchronization of the memristive FitzHugh-Nagumo neural networks described by this hybrid system. Importantly this quantitative threshold condition (4.2)-(4.3) is simply on the linear coupling coefficient \(P\) and the exponentially decaying rate is explicitly expressed by the given biological and mathematical parameters.

It is expected that modeling of biological and artificial neural networks by hybrid differential equations or other types of PDE in neuroscience and in deep learning field can be generalized with new features such as memristors and time delays. The
mathematical approach presented in this work can be further explored and extended with applications in a broad scope.

6. Numerical simulation of our model is presented and aligned with the theoretical proof. All the important constants in the proof are calculated and the threshold of coupling strength P is estimated based on the theoretical result. Synchronization of our example can be seen from the \(L^{2}\) norm figures. Visualization of our model as well as other hybrid equation models of network systems and further investigation of time steps and computational estimate of suboptimal thresholds needed to achieve some required synchronization is another interesting research problem.

**Acknowledgment**

Jing Tian's work is supported in part by the AMS Simons Travel Grant and the Jess and Mildred Fisher Endowed Professor fund of Mathematics from the Fisher College of Science and Mathematics at Towson University.
