Multiple Nonlinear Oscillations in a 𝔻₃×𝔻₃-Symmetrical Coupled System of Identical Cells with Delays

Abstract

A coupled system of nine identical cells with delays and 𝔻₃×𝔻₃-symmetry is considered. The individual cells are modelled by a scalar delay differential equation which includes linear decay and nonlinear delayed feedback. By analyzing the corresponding characteristic equations, the linear stability of the equilibrium is given. We also investigate the simultaneous occurrence of multiple periodic solutions and spatiotemporal patterns of the bifurcating periodic oscillations by using the equivariant bifurcation theory of delay differential equations combined with representation theory of Lie groups. Numerical simulations are carried out to illustrate our theoretical results.

1. Introduction

Over the past decades, symmetry has become a topic of considerable attention in the research of nonlinear dynamical systems [1–11]. In general, the symmetry reflects a certain spatial invariant of the dynamical systems. The work of Golubitsky et al. [1] shows that systems with symmetry can lead to multiple patterns of oscillation, which are predictable based on the theory of equivariant bifurcations. It is well known that the introduction of time delays into some systems is more reasonable and realistic [12]. Wu and coworkers [3, 13, 14] extended the theory of equivariant Hopf bifurcation to delay differential equations.

An artificial neural network is an information processing device that is inspired by the way biological nervous systems, such as the brain, process information simultaneously. It has many applications in different areas including pattern recognition, associative memory, signal processing, and combinatorial optimization. There has been an increasing interest in the investigation of neural networks (see, e.g., [4–6, 9–11, 15]) since Hopfield [16] constructed a simplified neural network model. Ring networks have been found in a variety of neural structures such as cerebellum [17] and even in chemistry and electrical engineering. In the field of neural networks, rings are studied to gain insight into the mechanisms underlying the behavior of recurrent networks [15, 18]. The dynamical behavior of ring networks has been investigated in more detail. For example, in order to understand which patterns occur in a particular system, Huang and Wu [4] studied the following ring neural network of three identical neurons with delayed feedback: 𝑥̇_𝑖(𝑡)=−𝑥_𝑖(𝑡)+𝑓𝑥_𝑖()−𝑔𝑥𝑡−𝜏_𝑖−1(𝑡−𝜏)+𝑔𝑥_𝑖+1(𝑡−𝜏),𝑖(mod3), (1) where 𝑥̇=d𝑥/d𝑡, 𝑥_𝑖(𝑡) represents the state of the 𝑖th neuron at time 𝑡, 𝑓 represents the nonlinear self-feedback function, 𝑔 is the connection function between neurons, and 𝜏≥0 is the time delay. Afterward, some researchers have been studying many ring networks with 𝔻_𝑛-symmetry (see [5–8]). However, previous work just has considered the individual network but not investigated the interactions between multiple networks.

In fact, a wide variety of natural and artificial systems possess a hierarchic structure or functioning and can naturally be modeled by coupled subnetwork. For example, the brain may be conceived as a dynamic network of coupled neurons. In order to describe the complicated interaction between billions of neurons in large neural network, the neurons are often lumped into highly connected subnetworks [19]. Coupled networks of nonlinear dynamical systems can exhibit rich dynamics, such as synchronization, symmetric bifurcation, and chaos. The spatio-temporal dynamics of systems of several coupled nonlinear oscillators is presently receiving great attention and a significant body of research has been carried out [9–11, 20, 21]. It must be pointed out that the hierarchical network of neuronal oscillators with 𝔻₃×𝔻₃-symmetry investigated in [20, 21] is described by a system of ordinary differential equations (ODEs), and the effect of time delays is not considered.

Motivated by the above ideas, in this paper, we consider the two-level hierarchical system which is composed of three coupled modules of interacting nonlinear neuron oscillators with time delays, modeled by the following system of delay differential equations (DDEs): 𝑥̇_𝑖,𝑗(𝑡)=−𝑥_𝑖,𝑗(𝑡)+𝑓𝑥_𝑖,𝑗(𝑡−𝜏)+𝑔𝑥_{𝑖,𝑗−1}(𝑥𝑡−𝜏)+𝑔_𝑖,𝑗+1(𝑡−𝜏)+ℎ𝑥_{𝑖−1,𝑗}()𝑥𝑡−𝜏+ℎ_𝑖+1,𝑗(𝑡−𝜏),𝑖,𝑗=0,1,2(mod3), (2) where 𝑓,𝑔,ℎ∈𝐶¹(ℝ,ℝ), 𝑓(0)=𝑔(0)=ℎ(0)=0, and ℎ represents the connection function between different modules. The individual elements are represented by a scalar equation, consisting of a linear decay term and a nonlinear time-delayed self-feedback. The architecture of the model is given in Figure 1.

Figure 1 

Architecture of model (2).

It is easy to see that all cells are identical, all couplings within each group are identical, and all groups are identically coupled to each other in this model. Therefore, these lead to a 𝔻₃×𝔻₃-symmetry of the associated system. the model is a natural extension of system (1) and is a particularly simple example of a symmetric system exhibiting a hierarchical structure with two levels: a “macro” level concerning the interactions between the groups and a “micro” level concerning the interactions within the groups. On the other hand, system (2) can be regarded as a special example of the general Hopfield neural networks with delays [16].

Although model (2) is a little simple, it would be of great significance for applications to have a detailed analysis and then to understand possible mechanisms behind the observed behaviour. In this paper, our main purpose is to reveal how the time delay can affect the stability of system (2), the simultaneous occurrence of multiple periodic solutions, and spatio-temporal patterns of the bifurcating periodic oscillations depending on the 𝔻₃×𝔻₃-symmetry.

The rest of the paper is organized as follows. In Section 2, the associated characteristic equation is analyzed and the linear stability of the equilibrium is given. In Section 3, we discuss the existence of multiple branches of periodic oscillations and their spatio-temporal patterns with the help of symmetric bifurcation theory of delay differential equations coupled with representation theory of Lie groups. An example and numerical simulations are presented to illustrate the results in Section 4. In Section 5, a brief discussion is drawn to conclude this paper.

2. Distribution of Characteristic Roots and Linear Stability

It is clear that (2) admits the trivial solution, 𝑥=0. The linearization of (2) at this equilibrium point is given by 𝑥̇_𝑖,𝑗(𝑡)=−𝑥_𝑖,𝑗(𝑡)+𝑎𝑥_𝑖,𝑗(𝑡−𝜏)+𝑏𝑥_{𝑖,𝑗−1}(𝑡−𝜏)+𝑏𝑥_𝑖,𝑗+1(𝑡−𝜏)+𝑐𝑥_{𝑖−1,𝑗}(𝑡−𝜏)+𝑐𝑥_𝑖+1,𝑗(𝑡−𝜏),𝑖,𝑗=0,1,2(mod3), (3) where 𝑎=𝑓^′(0),𝑏=𝑔^′(0), and 𝑐=ℎ^′(0). Regarding 𝜏 as the parameter, let 𝐴(𝜏) denote the infinitesimal generator of the semigroup generated by linear system (3). We first determine when 𝐴(𝜏) has a pair of purely imaginary eigenvalues.

The characteristic matrix of (3) is Δ(𝜏,𝜆)=(𝜆+1)Id₉−𝑀e^−𝜆𝜏,𝜆∈ℂ, (4) where Id_𝑛 denotes the identity matrix of order 𝑛, 𝑀=circ(𝑀₁,𝑐Id₃,𝑐Id₃) is a circle block matrix, and 𝑀₁=circ(𝑎,𝑏,𝑏) is a circulant matrix of order 3. Then we have the following lemma.

Lemma 1. The associated characteristic equation of (3) is detΔ(𝜏,𝜆)=²²_𝑞=0𝑝=0Δ_𝑝𝑞=0, (5) where Δ_𝑝𝑞=𝜆+1−(𝑎+2𝑏cos(2𝑞𝜋/3)+2𝑐cos(2𝑝𝜋/3))e^−𝜆𝜏.

Proof. Let 𝜒=e^i(2𝜋/3),𝑣_𝑞=(1,𝜒^𝑞,𝜒^2𝑞)^T, and 𝑣_𝑝𝑞=(𝑣_𝑞,𝜒^𝑝𝑣_𝑞,𝜒^2𝑝𝑣_𝑞)^T. Then 𝑀𝑣_𝑝𝑞=𝑎+2𝑏cos2𝑞𝜋3+2𝑐cos2𝑝𝜋3𝑣_𝑝𝑞. (6) Hence Δ(𝜏,𝜆)𝑣_𝑝𝑞=𝜆+1−𝑎+𝑐𝜒^𝑝+𝑐𝜒^−𝑝+𝑏𝜒^𝑞+𝑏𝜒^−𝑞e^−𝜆𝜏𝑣_𝑝𝑞=𝜆+1−𝑎+2𝑏cos2𝑞𝜋3+2𝑐cos2𝑝𝜋3e^−𝜆𝜏𝑣_𝑝𝑞. (7) The conclusion is obtained.

Note that Δ₀₀=𝜆+1−(𝑎+2𝑏+2𝑐)e^−𝜆𝜏,Δ₀₁=Δ₀₂=𝜆+1−(𝑎+2𝑐−𝑏)e^−𝜆𝜏,Δ₁₀=Δ₂₀=𝜆+1−(𝑎+2𝑏−𝑐)e^−𝜆𝜏,Δ₁₁=Δ₁₂=Δ₂₁=Δ₂₂=𝜆+1−(𝑎−𝑏−𝑐)e^−𝜆𝜏, (8) we have detΔ(𝜏,𝜆)=Δ₀₀Δ₀₁²Δ₁₀²Δ₁₁⁴. (9) In this paper, for the sake of simplicity, we consider only the case where characteristic equation (5) may have a pair of purely imaginary roots of multiplicity 4, that is, focus on the distribution of zeros of the factor Δ₁₁. The other cases are simpler and can be handled in a similar way, and we omit them. Therefore, throughout this paper we suppose(𝐇)|𝑎+2𝑏+2𝑐|<1, |𝑎+2𝑐−𝑏|<1, |𝑎+2𝑏−𝑐|<1, 𝑎−𝑏−𝑐<−1.The following result about the distribution of the characteristic roots is obtained.

Lemma 2. Assume that (𝐇) holds. Define 𝜏_𝑠=arccos(1/(𝑎−𝑏−𝑐))+2𝑠𝜋(𝑎−𝑏−𝑐)²−1,𝑠=0,1,2,…,𝛽=(𝑎−𝑏−𝑐)²−1. (10)(11) Then(i)for all 𝜏⩾0, all zeros of the factors Δ₀₀,Δ₁₀, and Δ₀₁ have negative real parts,(ii)when 𝜏∈[0,𝜏₀), all roots of the characteristic equation (5) have negative real parts; when 𝜏∈(𝜏_𝑠,𝜏_𝑠+1), the characteristic equation (5) has exactly 2𝑠+2 roots with positive real parts; the other roots have negative real parts; at (and only at) 𝜏=𝜏_𝑠, 𝐴(𝜏) has a pair of purely imaginary eigenvalues ±i𝛽 of multiplicity 4, and all other eigenvalues of 𝐴(𝜏) are not integer multiples of ±i𝛽,(iii)for each fixed 𝑠∈ℕ, there exist 𝛿>0 and 𝐶¹-smooth mapping 𝜆∶(𝜏_𝑠−𝛿,𝜏_𝑠+𝛿)→ℂ such that 𝜆(𝜏_𝑠)=i𝛽 and 𝜆(𝜏)+1−(𝑎−𝑏−𝑐)e^{−𝜏𝜆(𝜏)}=0 for all 𝜏∈(𝜏_𝑠−𝛿,𝜏_𝑠+𝛿). Moreover, (d/d𝜏)Re𝜆(𝜏)|_{𝜏=𝜏𝑠}>0,(iv)the generalized eigenspace 𝑈_±i𝛽(𝐴(𝜏_𝑠)) consists of eigenvector of 𝐴(𝜏_𝑠) associated with ±i𝛽. Moreover, 𝑈_±i𝛽𝐴𝜏_𝑠=²²_𝑞=1𝑝=1𝑦_𝑝𝑞𝜖^𝑝𝑞+𝑧_𝑝𝑞𝜍^𝑝𝑞,𝑦_𝑝𝑞,𝑧_𝑝𝑞∈ℝ,𝑝,𝑞=1,2, (12)where 𝜖^𝑝𝑞(𝜃)=Ree^i𝛽𝜃𝑣_𝑝𝑞=cos𝛽𝜃Re𝑣_𝑝𝑞−sin𝛽𝜃Im𝑣_𝑝𝑞,𝜍^𝑝𝑞(𝜃)=Ime^i𝛽𝜃𝑣_𝑝𝑞=sin𝛽𝜃Re𝑣_𝑝𝑞+cos𝛽𝜃Im𝑣_𝑝𝑞,𝑓𝑜𝑟𝜃∈[−𝜏,0]. (13)

Proof. For 𝛾∈ℝ, let 𝑞_𝛾(𝜆)=𝜆+1−𝛾e^−𝜆𝜏. If 𝜆=𝛼+i𝛽 is a zero of 𝑞_𝛾(𝜆), then 1+𝛼+i𝛽=𝛾e^{−(𝛼+i𝛽)𝜏}, from which it follows that (1+𝛼)²+𝛽²=𝛾²e^−2𝛼𝜏. (14) Thus, we claim that all the zeros of 𝑞_𝛾(𝜆) have negative real parts provided that |𝛾|<1. Noticing that (𝐇) holds and applying the above discussions to the factors Δ₀₀,Δ₁₀, and Δ₀₁, we obtain conclusion (i).
In what follows, we consider the distribution of zeros of the factor Δ₁₁. For 𝜏=0, Δ₁₁=0 becomes 𝜆=𝑎−𝑏−𝑐−1<−2. Let i𝛽(𝛽>0) be a root of Δ₁₁=0, then i𝛽+1−(𝑏+𝑐−𝑎)e^{i(𝜋−𝜏𝛽)}=1+𝛽²e^{iarccos(1/1+𝛽2)}−(𝑏+𝑐−𝑎)e^{i(𝜋−𝜏𝛽)}=0. (15) Thus, 1+𝛽²=𝑏+𝑐−𝑎,arc1cos1+𝛽²=𝜋−𝜏𝛽+2𝑠𝜋,𝑠∈ℤ. (16) Therefore, Δ₁₁=0 has a pair of purely imaginary roots ±i𝛽 if and only if 𝛽=(𝑎−𝑏−𝑐)²−1,𝜏=arccos(1/(𝑎−𝑏−𝑐))+2𝑠𝜋(𝑎−𝑏−𝑐)²−1∶=𝜏_𝑠,forsome𝑠∈ℕ. (17) Clearly, for each fixed 𝑠, 𝜏_𝑠⩾0. Since dd𝜆Δ₁₁(𝜆)_{𝜆=i𝛽,𝜏=𝜏𝑠}=1+(𝑎−𝑏−𝑐)𝜏e^−𝜆𝜏_{𝜆=i𝛽,𝜏=𝜏𝑠}=1+𝜏_𝑠i𝛽+1≠0, (18) there exist 𝛿>0 and 𝐶¹-smooth mapping 𝜆∶(𝜏_𝑠−𝛿,𝜏_𝑠+𝛿)→ℂ such that 𝜆(𝜏_𝑠)=i𝛽 and 𝜆(𝜏)+1−(𝑎−𝑏−𝑐)e^{−𝜏𝜆(𝜏)}=0 for all 𝜏∈(𝜏_𝑠−𝛿,𝜏_𝑠+𝛿). Differentiating Δ₁₁=0 with respect to 𝜏, we have 𝜆^′𝜏_𝑠=−(𝑎−𝑏−𝑐)𝜆𝜏_𝑠e^{−𝜆(𝜏𝑠)𝜏𝑠}1+(𝑎−𝑏−𝑐)𝜏_𝑠e^{−𝜆(𝜏𝑠)𝜏𝑠}=−𝜆𝜏_𝑠𝜆𝜏_𝑠+11+𝜏_𝑠𝜆𝜏_𝑠.+1 (19) Therefore, Re𝜆^′𝜏_𝑠=𝛽²1+𝜏_𝑠²+𝜏²_𝑠𝛽²>0. (20) So far, the proof of (ii) and (iii) is complete (for more details see XI.2 in [22]). It remains to verify (iv).
Note that 𝑣₁₁=𝑣₂₂,𝑣₁₂=𝑣₂₁; it follows from the proof of Lemma 1 and the above discussions that the eigenspace of 𝐴(𝜏_𝑠) associated with ±i𝛽 is spanned by e^i𝛽𝜃𝑣₁₁,e^i𝛽𝜃𝑣₁₂,e^i𝛽𝜃𝑣₂₁,e^i𝛽𝜃𝑣₂₂,e^−i𝛽𝜃𝑣₁₁, e^−i𝛽𝜃𝑣₁₂, e^−i𝛽𝜃𝑣₂₁, and e^−i𝛽𝜃𝑣₂₂. Hence, the space has the real basis {𝜖¹¹,𝜍¹¹,𝜖¹²,𝜍¹²,𝜖²¹,𝜍²¹,𝜖²²,𝜍²²}. On the other hand, the eigenspace 𝐴(𝜏_𝑠) associated with i𝛽_𝑝𝑞 is of dimension 4, and the multiplicity of the characteristic root 𝜆=i𝛽_𝑝𝑞 is also 4. According to the folk theorem in functional differential equations (see [23]), 𝑈_±i𝛽(𝐴(𝜏_𝑠)) must coincide with the eigenspace of 𝐴(𝜏_𝑠) associated with ±i𝛽. Therefore, conclusion (iv) is correct. This completes the proof.

From Lemma 2 (ii), we can draw the following conclusion about the linear stability of system (2).

Theorem 3. If the assumption (𝐇) is satisfied, then the equilibrium 𝑥=0 of system (2) is asymptotically stable for 𝜏∈[0,𝜏₀); the equilibrium 𝑥=0 is unstable for 𝜏>𝜏₀.

3. Multiple Patterns of Oscillation

It follows from Lemma 2 that 𝐴(𝜏) has a pair of purely imaginary eigenvalues ±i𝛽 of multiplicity 4 and all other eigenvalues of 𝐴(𝜏) are not integer multiples of ±i𝛽 at 𝜏=𝜏_𝑠(𝑠=0,1,…). Thus, the Hopf bifurcation may provide some asynchronous periodic solutions at each 𝜏_𝑠.

The symmetry of a system is important in determining the patterns of oscillation. In order to explore the symmetry of system and analyze the spatio-temporal patterns of the bifurcated periodic solutions, we need the following definition.

Let 𝔻₃ denote the dihedral group of order 6, which is generated by cyclic group ℤ₃ together with the flip 𝜅 of order 2. Define the action of Γ∶=𝔻₃×𝔻₃ on ℝ⁹ by 𝜌,1⋅𝑥_𝑖,𝑗=𝑥_𝑖+1,𝑗,1,𝜌⋅𝑥_𝑖,𝑗=𝑥_𝑖,𝑗+1,((𝜅,1)⋅𝑥)_𝑖,𝑗=𝑥_{3−𝑖,𝑗},((1,𝜅)⋅𝑥)_𝑖,𝑗=𝑥_{𝑖,3−𝑗},𝑖,𝑗mod3, (21) where 𝜌 is the generator of ℤ₃ and 𝜅 denotes the flip.

Lemma 4. System (2) is Γ-equivalent.

Proof. Let mapping ℱ∶𝐶([−𝜏,0],ℝ⁹)→ℝ⁹ be ℱ𝜙_𝑖,𝑗=−𝜙_𝑖,𝑗(0)+𝑓𝜙_𝑖,𝑗(𝜙−𝜏)+𝑔_{𝑖,𝑗−1}(−𝜏)+𝑔𝜙_𝑖,𝑗+1()𝜙−𝜏+ℎ_{𝑖−1,𝑗}()𝜙−𝜏+ℎ_𝑖+1,𝑗(−𝜏),𝑖,𝑗(mod3), (22) where 𝜙=𝜙₀₀,𝜙₀₁,𝜙₀₂,𝜙₁₀,𝜙₁₁,𝜙₁₂,𝜙₂₀,𝜙₂₁,𝜙₂₂^T∈𝐶[−𝜏,0],ℝ⁹,ℱ𝜙=ℱ𝜙₀₀,ℱ𝜙₀₁,ℱ𝜙₀₂,ℱ𝜙₁₀,ℱ𝜙₁₁,ℱ𝜙₁₂,ℱ𝜙₂₀,ℱ𝜙₂₁,ℱ𝜙₂₂^T∈ℝ⁹. (23) Then ℱ𝜌,1𝜙_𝑖,𝑗=−𝜌,1𝜙_𝑖,𝑗(0)+𝑓𝜌,1𝜙_𝑖,𝑗()𝜙−𝜏+𝑔𝜌,1_{𝑖,𝑗−1}(−𝜏)+𝑔𝜌,1𝜙_𝑖,𝑗+1()𝜙−𝜏+ℎ𝜌,1_{𝑖−1,𝑗}(−𝜏)+ℎ𝜌,1𝜙_𝑖+1,𝑗(𝜙−𝜏)=−_𝑖+1,𝑗(0)+𝑓𝜙_𝑖+1,𝑗()𝜙−𝜏+𝑔_{𝑖+1,𝑗−1}()𝜙−𝜏+𝑔_{𝑖+1,𝑗+1}(−𝜏)+ℎ𝜙_𝑖,𝑗()𝜙−𝜏+ℎ_𝑖+2,𝑗()=ℱ𝜙−𝜏𝜌,1_𝑖,𝑗,ℱ(𝜅,1)𝜙_𝑖,𝑗=−(𝜅,1)𝜙_𝑖,𝑗(0)+𝑓(𝜅,1)𝜙_𝑖,𝑗()()𝜙−𝜏+𝑔𝜅,1_{𝑖,𝑗−1}(−𝜏)+𝑔(𝜅,1)𝜙_𝑖,𝑗+1()()𝜙−𝜏+ℎ𝜅,1_{𝑖−1,𝑗}(−𝜏)+ℎ(𝜅,1)𝜙_𝑖+1,𝑗()𝜙−𝜏=−_{3−𝑖,𝑗}(0)+𝑓𝜙_{3−𝑖,𝑗}()𝜙−𝜏+𝑔_{3−𝑖,𝑗−1}(−𝜏)+𝑔𝜙_{3−𝑖,𝑗+1}()𝜙−𝜏+ℎ_{4−𝑖,𝑗}(−𝜏)+ℎ𝜙_{2−𝑖,𝑗}()=()ℱ𝜙−𝜏𝜅,1_𝑖,𝑗. (24) Similarly, we can prove that ℱ1,𝜌𝜙_𝑖,𝑗=1,𝜌ℱ𝜙_𝑖,𝑗,ℱ(1,𝜅)𝜙_𝑖,𝑗=(1,𝜅)ℱ𝜙_𝑖,𝑗. (25) Therefore, ℱ is 𝔻₃×𝔻₃-equivalent. This completes the proof.

Lemma 5. Let Γ act on ℝ⁴ by ⋅𝑥𝜌,1¹𝑥₂𝑥₃𝑥₄=−12√0−32010−2√0−32√3210−200√3210−2𝑥₁𝑥₂𝑥₃𝑥₄,⋅𝑥1,𝜌¹𝑥₂𝑥₃𝑥₄=−12−√32√0032−1210000−2−√32√0032−12𝑥₁𝑥₂𝑥₃𝑥₄,()⋅𝑥𝜅,1¹𝑥₂𝑥₃𝑥₄=𝑥₁𝑥₂−𝑥₃−𝑥₄,()⋅𝑥1,𝜅¹𝑥₂𝑥₃𝑥₄=𝑥₁−𝑥₂𝑥₃−𝑥₄. (26) Then ℝ⁴ is an absolutely irreducible representation of Γ, and the restricted action of Γ on KerΔ(𝜏,i𝛽) is isomorphic to the action of  Γ on ℝ⁴⊕ℝ⁴.

Proof. It is straightforward to verify the absolute irreducibility of the representation of Γ on ℝ⁴ by the definition (see [1]). Note that KerΔ𝜏,i𝛽=²²_𝑝=1𝑞=1𝑎_𝑝𝑞+i𝑏_𝑝𝑞𝑣_𝑝𝑞,𝑎_𝑝𝑞,𝑏_𝑝𝑞∈ℝ,𝑝,𝑞=1,2. (27) Define 𝐽∶KerΔ(𝜏,i𝛽)→ℝ⁴⊕ℝ⁴ by ²²_𝑝=1𝑞=1𝑎_𝑝𝑞+i𝑏_𝑝𝑞𝑣_𝑝𝑞↦𝐵𝑎₁₁,𝑏₁₁,𝑎₁₂,𝑏₁₂,𝑎₂₁,𝑏₂₁,𝑎₂₂,𝑏₂₂^T, (28) where the matrix 𝐵=.10101010010−1010−101010−10−1−101010−1001010101−1010−1010−10−1010100−101010−1 (29) The nonsingularity of the matrix 𝐵 implies that 𝐽∶KerΔ(𝜏,i𝛽)≅ℝ⁴ is a linear isomorphism. It is easy to see that 𝜌,1⋅²²_𝑝=1𝑞=1𝑎_𝑝𝑞+i𝑏_𝑝𝑞𝑣_𝑝𝑞=²²_𝑝=1𝑞=1𝑎_𝑝𝑞+i𝑏_𝑝𝑞e^{i((−1)𝑝+12𝜋/3)}𝑣_𝑝𝑞,1,𝜌⋅²²_𝑝=1𝑞=1𝑎_𝑝𝑞+i𝑏_𝑝𝑞𝑣_𝑝𝑞=²²_𝑝=1𝑞=1𝑎_𝑝𝑞+i𝑏_𝑝𝑞e^{i((−1)𝑞+12𝜋/3)}𝑣_𝑝𝑞,(𝜅,1)⋅²²_𝑝=1𝑞=1𝑎_𝑝𝑞+i𝑏_𝑝𝑞𝑣_𝑝𝑞=𝑎₂₁+i𝑏₂₁𝑣₁₁+𝑎₂₂+i𝑏₂₂𝑣₁₂+𝑎₁₁+i𝑏₁₁𝑣₂₁+𝑎₁₂+i𝑏₁₂𝑣₂₂,(1,𝜅)⋅²²_𝑝=1𝑞=1𝑎_𝑝𝑞+i𝑏_𝑝𝑞𝑣_𝑝𝑞=𝑎₁₂+i𝑏₁₂𝑣₁₁+𝑎₁₁+i𝑏₁₁𝑣₁₂+𝑎₂₂+i𝑏₂₂𝑣₂₁+𝑎₂₁+i𝑏₂₁𝑣₂₂. (30) Therefore, a straightforward calculation shows that 𝐽𝜌,1⋅²²_𝑝=1𝑞=1𝑎_𝑝𝑞+i𝑏_𝑝𝑞𝑣_𝑝𝑞=𝜌,1⋅𝐽²²_𝑝=1𝑞=1𝑎_𝑝𝑞+i𝑏_𝑝𝑞𝑣_𝑝𝑞,𝐽1,𝜌⋅²²_𝑝=1𝑞=1𝑎_𝑝𝑞+i𝑏_𝑝𝑞𝑣_𝑝𝑞=1,𝜌⋅𝐽²²_𝑝=1𝑞=1𝑎_𝑝𝑞+i𝑏_𝑝𝑞𝑣_𝑝𝑞,𝐽(𝜅,1)⋅²²_𝑝=1𝑞=1𝑎_𝑝𝑞+i𝑏_𝑝𝑞𝑣_𝑝𝑞=(𝜅,1)⋅𝐽²²_𝑝=1𝑞=1𝑎_𝑝𝑞+i𝑏_𝑝𝑞𝑣_𝑝𝑞,𝐽(1,𝜅)⋅²²_𝑝=1𝑞=1𝑎_𝑝𝑞+i𝑏_𝑝𝑞𝑣_𝑝𝑞=(1,𝜅)⋅𝐽²²_𝑝=1𝑞=1𝑎_𝑝𝑞+i𝑏_𝑝𝑞𝑣_𝑝𝑞. (31) This concludes the proof.

Let 𝜔=2𝜋/𝛽 and let 𝑃_𝜔 be the Banach space of all continuous 𝜔-periodic functions 𝑥∶ℝ→ℝ⁹. Then for the circle group 𝑆¹,  Γ×𝑆¹ acts on 𝑃_𝜔 by 𝛾,e^i𝜃(𝑡)𝜔⋅𝑥=𝛾⋅𝑥𝑡+𝜃,2𝜋𝛾,e^i𝜃∈Γ×𝑆¹,𝑥∈𝑃_𝜔. (32) Denote by 𝑆𝑃_𝜔 the subspace of 𝑃_𝜔 consisting of all 𝜔-periodic solutions of system (3) with 𝜏=𝜏_𝑠. Then 𝑆𝑃_𝜔=²²_𝑞=1𝑝=1𝑦_𝑝𝑞𝜖^𝑝𝑞+𝑧_𝑝𝑞𝜍^𝑝𝑞,𝑦_𝑝𝑞,𝑧_𝑝𝑞∈ℝ,𝑝,𝑞=1,2, (33) where 𝜖^𝑝𝑞(𝑡)=Ree^i𝛽𝑡𝑣_𝑝𝑞=cos𝛽𝑡Re𝑣_𝑝𝑞−sin𝛽𝑡Im𝑣_𝑝𝑞,𝜍^𝑝𝑞(𝑡)=Ime^i𝛽𝑡𝑣_𝑝𝑞=sin𝛽𝑡Re𝑣_𝑝𝑞+cos𝛽𝑡Im𝑣_𝑝𝑞. (34) Therefore, for 𝜖^𝑝𝑞(𝑡),𝜍^𝑝𝑞(𝑡)(𝑝,𝑞=1,2), we have the following properties.

Lemma 6. 𝜌,1⋅𝜖^𝑝𝑞=𝜖^𝑝𝑞cos2𝑝𝜋3−𝜍^𝑝𝑞sin2𝑝𝜋3,𝜌,1⋅𝜍^𝑝𝑞=𝜍^𝑝𝑞cos2𝑝𝜋3+𝜖^𝑝𝑞sin2𝑝𝜋3,1,𝜌⋅𝜖^𝑝𝑞=𝜖^𝑝𝑞cos2𝑞𝜋3−𝜍^𝑝𝑞sin2𝑞𝜋3,1,𝜌⋅𝜍^𝑝𝑞=𝜍^𝑝𝑞cos2𝑞𝜋3+𝜖^𝑝𝑞sin2𝑞𝜋3,(𝜅,1)⋅𝜖^𝑝𝑞=𝜖^{3−𝑝,𝑞},(𝜅,1)⋅𝜍^𝑝𝑞=𝜍^{3−𝑝,𝑞},(1,𝜅)⋅𝜖^𝑝𝑞=𝜖^{𝑝,3−𝑞},(1,𝜅)⋅𝜍^𝑝𝑞=𝜍^{𝑝,3−𝑞}. (35)

Proof. For 𝑖,𝑗(mod3) and 𝑡∈ℝ, note that 𝜖^𝑝𝑞_𝑖,𝑗(𝑡)=cos𝛽𝑡+2𝑖𝑝𝜋3+2𝑗𝑞𝜋3,𝜍^𝑝𝑞_𝑖,𝑗(𝑡)=sin𝛽𝑡+2𝑖𝑝𝜋3+2𝑗𝑞𝜋3, (36) we have 𝜌,1⋅𝜖^𝑝𝑞(𝑡)_𝑖,𝑗=𝜖^𝑝𝑞_𝑖+1,𝑗(𝑡)2()=cos𝛽𝑡+𝑖+1𝑝𝜋3+2𝑗𝑞𝜋3=cos𝛽𝑡+2𝑖𝑝𝜋3+2𝑗𝑞𝜋3+2𝑝𝜋3=𝜖^𝑝𝑞_𝑖,𝑗cos2𝑝𝜋3−𝜍^𝑝𝑞_𝑖,𝑗sin2𝑝𝜋3,(𝜌,1)⋅𝜍^𝑝𝑞(𝑡)_𝑖,𝑗=𝜍^𝑝𝑞_𝑖+1,𝑗(𝑡)2()=sin𝛽𝑡+𝑖+1𝑝𝜋3+2𝑗𝑞𝜋3=sin𝛽𝑡+2𝑖𝑝𝜋3+2𝑗𝑞𝜋3+2𝑝𝜋3=𝜍^𝑝𝑞_𝑖,𝑗cos2𝑝𝜋3+𝜖^𝑝𝑞_𝑖,𝑗sin2𝑝𝜋3,(𝜅,1)⋅𝜖^𝑝𝑞(𝑡)_𝑖,𝑗=𝜖^𝑝𝑞_{3−𝑖,𝑗}(𝑡)2()=cos𝛽𝑡+3−𝑖𝑝𝜋3+2𝑗𝑞𝜋3𝜋=cos𝛽𝑡+2𝑖3−𝑝3+2𝑗𝑞𝜋3=𝜖^{3−𝑝,𝑞}_𝑖,𝑗,(𝜅,1)⋅𝜍^𝑝𝑞(𝑡)_𝑖,𝑗=𝜖^𝑝𝑞_{3−𝑖,𝑗}(𝑡)2()=sin𝛽𝑡+3−𝑖𝑝𝜋3+2𝑗𝑞𝜋3𝜋=sin𝛽𝑡+2𝑖3−𝑝3+2𝑗𝑞𝜋3=𝜍^{3−𝑝,𝑞}_𝑖,𝑗. (37) Therefore 𝜌,1⋅𝜖^𝑝𝑞=𝜖^𝑝𝑞cos2𝑝𝜋3−𝜍^𝑝𝑞sin2𝑝𝜋3,𝜌,1⋅𝜍^𝑝𝑞=𝜍^𝑝𝑞cos2𝑝𝜋3+𝜖^𝑝𝑞sin2𝑝𝜋3,(𝜅,1)⋅𝜖^𝑝𝑞=𝜖^{3−𝑝,𝑞},(𝜅,1)⋅𝜍^𝑝𝑞=𝜍^{3−𝑝,𝑞}. (38) Similarly, we can prove that 1,𝜌⋅𝜖^𝑝𝑞=𝜖^𝑝𝑞cos2𝑞𝜋3−𝜍^𝑝𝑞sin2𝑞𝜋3,1,𝜌⋅𝜍^𝑝𝑞=𝜍^𝑝𝑞cos2𝑞𝜋3+𝜖^𝑝𝑞sin2𝑞𝜋3,(1,𝜅)⋅𝜖^𝑝𝑞=𝜖^{𝑝,3−𝑞},(1,𝜅)⋅𝜍^𝑝𝑞=𝜍^{𝑝,3−𝑞}. (39)

It is clear that if 𝑥 is a periodic solution of system (2), then so is (𝛾,e^i𝜃)𝑥 for every (𝛾,e^i𝜃)∈Γ×𝑆¹. The spatial-temporal symmetry of a bifurcation of periodic solutions 𝑥(𝑡) can be completely characterized by the isotropy group Σ_𝑥={(𝛾,e^i𝜃)∈Γ×𝑆¹∣(𝛾,e^i𝜃)𝑥=𝑥}≤Γ×𝑆¹, and it is easy to verify that the isotropy group of (𝛾,e^i𝜃)𝑥 is (𝛾,e^i𝜃)Σ_𝑥(𝛾,e^i𝜃)⁻¹, which is conjugate to Σ_𝑥. The maximal isotropy subgroups Σ^𝑚(𝑚=1,2,…,20) of Γ×𝑆¹ are listed in Table 1. For each subgroup Σ^𝑚, the Σ^𝑚-fixed-point set FixΣ^𝑚,𝑆𝑃_𝜔=𝑥∈𝑆𝑃_𝜔∣𝛾,e^i𝜃𝑥=𝑥,∀𝛾,e^i𝜃∈Σ^𝑚 (40) is a subspace of 𝑆𝑃_𝜔. According to Lemma 6, it is easy to verify that dimFix(Σ^𝑚,𝑆𝑃_𝜔)=2,𝑚=1,2,…,20.

Table 1 

The maximal isotropy subgroups of Γ×𝑆¹ and associated fixed point subspaces.

Together with Lemma 2, Lemmas 4–6 allow us to apply the equivariant Hopf bifurcation theorem for delay differential equations due to Wu [3] to obtain the following result on the spatio-temporal patterns of the bifurcated periodic solutions.

Theorem 7. Assume that (𝐇) is satisfied. Then for system (2), near each 𝜏_𝑠(𝑠=0,1,…) there exist 100 distinct branches of periodic solutions bifurcated from the equilibrium 𝑥=0. More precisely, for each isotropy group Σ^𝑚(𝑚=1,2,…,20) and a chosen basis {𝛿₁,𝛿₂} of Fix(Σ^𝑚,𝑆𝑃_𝜔), there exist 𝛼₀>0,𝜏^∗₀>0,𝜎₀>0, and a 𝐶¹-smooth mapping (𝜏^∗,𝜔^∗,𝑥^∗)∶ℝ²_𝛼0→ℝ×ℝ⁺×𝐶(ℝ,ℝ⁹), where ℝ²_𝛼0={𝛼∈ℝ²∣|𝛼|<𝛼₀}, such that for each 𝛼∈ℝ²_𝛼0, 𝑥^∗=𝑥^∗(𝑡;𝛼) is an 𝜔^∗(𝛼)-periodic solution of system (2) with 𝜏=𝜏_𝑠+𝜏^∗(𝛼), and 𝛾⋅𝑥^∗(𝑡)=𝑥^∗𝜔𝑡−^∗(𝛼)𝜃,2𝜋𝛾,e^i𝜃∈Σ^𝑚,𝜔^∗(0)=2𝜋𝛽,𝜏^∗(0)𝑥=0,_∗(𝑡;𝛼)=𝛿₁,𝛿₂𝛼+𝑜(|𝛼|),𝑎𝑠|𝛼|⟶0. (41) Furthermore, for |𝜏−𝜏_𝑠|<𝜏^∗₀, |𝜔−2𝜋/𝛽|<𝜎₀, every 𝜔-periodic solution 𝑥(𝑡) of system (2) with ‖𝑥_𝑡‖<𝜎₀, 𝛾⋅𝑥(𝑡)=𝑥(𝑡−(𝜔/2𝜋)𝜃) for (𝛾,e^i𝜃)∈Σ^𝑚 and 𝑡∈ℝ must be of the above type.

4. An Example and Numerical Simulations

As a simple example, we consider the following specific Hopfield model with 𝔻₃×𝔻₃-symmetry: 𝑥̇_𝑖,𝑗(𝑡)=−𝑥_𝑖,𝑗(𝑡)−0.6tanh𝑥_𝑖,𝑗(𝑥𝑡−𝜏)+0.3tanh_{𝑖,𝑗−1}()𝑥𝑡−𝜏+0.3tanh_𝑖,𝑗+1()𝑥𝑡−𝜏+0.4tanh_{𝑖−1,𝑗}()𝑥𝑡−𝜏+0.4tanh_𝑖+1,𝑗(𝑡−𝜏), (42) where 𝑖,𝑗=0,1,2(mod3).

Clearly, the origin is an equilibrium of system (42). It is easy to compute 𝑎=−0.6,𝑏=0.3,and𝑐=0.4 and verify that the hypothesis (𝐇) holds. According to (10), we obtain 𝜏_𝑠=2.9476+7.5641𝑠,𝑠=0,1,…. (43) From Theorem 3, the origin is asymptotically stable if 𝜏<𝜏₀=2.9476. It follows from Theorem 7 that the equivariant Hopf bifurcation occurs at 𝜏₀=2.9476 and there exist 100 distinct branches of periodic solutions bifurcated from the origin. To illustrate the analytical results found, we give some numerical simulations. Figure 2 shows that the origin of system (42) is stable when 𝜏=2.6. Figure 3 shows that four periodic orbits occur simultaneously when 𝜏=3. Unfortunately, we cannot verify the existence of all other bifurcated periodic orbits since they may be unstable.

(a)

(b)

(c)

Figure 2 

Trajectories 𝑥𝑖,𝑗(𝑖,𝑗=0,1,2) of system (42) when 𝜏=2.6.

Figure 3 

System (42) has multiple periodic solutions when 𝜏=3.

5. Discussion

In this paper, we have studied a coupled system of nine identical cells with delays and 𝔻₃×𝔻₃-symmetry. By choosing the time delay 𝜏 as a bifurcation parameter and analyzing the corresponding characteristic equation, we have shown that under some assumption, the equilibrium of the model loses its stability and periodic solutions via Hopf bifurcation occur when 𝜏 passes through a critical value. This implies that the time delay can be regarded as a source of instability and oscillatory response of the networks and is able to alter the dynamics of system (2) significantly. Moreover, the spatio-temporal patterns of bifurcating periodic solutions are explored clearly by employing the symmetric bifurcation theory of delay differential equations combined with representation theory of Lie groups. From Theorem 7, we have obtained the conclusion that the small-scale network with a special structure may have a large number of periodic oscillations. Therefore, it is natural that the large-scale network possesses complicated dynamics generally.

Further investigations such as the stability, direction and global existence of the periodic solutions bifurcating from the local Hopf bifurcations are essential in order to fully understand the periodic phenomenon of the system. We can compute the normal forms directly by using the method due to Faria and Magalháes [24, 25]. However, this is a complex and prolix task. In addition, we would like to point out that codimension two mode interactions may take place if the assumption (𝐇) is not satisfied. For example, if 𝑎+2𝑏+2𝑐=1, |𝑎+2𝑐−𝑏|<1,|𝑎+2𝑏−𝑐|<1, and 𝑎−𝑏−𝑐<−1, then system (2) undergoes a fold-Hopf bifurcation when 𝜏 passes through a critical value. Therefore, it is possible to study secondary bifurcations and more complex behaviours in this coupled network. We leave them for our future work.

Acknowledgments

This work is partially supported by National Natural Science Foundation of China (nos. 11101435 and 11101053), the Scientific Research Fund of Hunan Provincial Education Department of China (no. 12C0028), and the Open Fund Project of Key Research Institute of Philosophies and Social Sciences in Hunan Universities (no. 12FEFM07).