{'fraction_non_alphanumeric': 0.04295727767552962, 'fraction_numerical': 0.01883541527979011, 'mean_word_length': 4.43593018305662, 'pattern_counts': {'":': 0, '<': 4, '<?xml version=': 0, '>': 4, 'https://': 0, 'lorem ipsum': 0, 'www.': 0, 'xml': 0}, 'pii_count': 0, 'substrings_counts': 0, 'word_list_counts': {'cursed_substrings.json': 16, 'profanity_word_list.json': 0, 'sexual_word_list.json': 0, 'zh_pornsignals.json': 0}}
Double Hopf bifurcation in delayed reaction-diffusion systems

Double Hopf bifurcation analysis can be used to reveal some complicated dynamical behavior in a dynamical system, such as the existence or coexistence of periodic orbits, quasi-periodic orbits, or even chaos. In this paper, an algorithm for deriving the normal form near a codimension-two double Hopf bifurcation of a reaction-diffusion system with time delay and Neumann boundary condition is rigorously established, by employing the center manifold reduction technique and the normal form method. We find that the dynamical behavior near bifurcation points are proved to be governed by twelve distinct unfolding systems. Two examples are performed to illustrate our results: for a stage-structured epidemic model, we find that double Hopf bifurcation appears when varying the diffusion rate and time delay, and two stable spatially inhomogeneous periodic oscillations are proved to coexist near the bifurcation point; in a diffusive predator-prey system, we theoretically proved that quasi-periodic orbits exist on two- or three-torus near a double Hopf bifurcation point, which will break down after slight perturbation, leaving the system a strange attractor.


INTRODUCTION
A Hopf bifurcation refers to the phenomenon that steady states lose their stability and give rise to periodic solutions, as a parameter crosses critical values [2,29,33,42,51]. Since the condition of existence of Hopf bifurcation can be easily verified, and both the direction of Hopf bifurcation and the stability of the bifurcating periodic solution can be determined by formulas derived from center manifold reduction technique, Hopf bifurcation analysis has been an effective method of the study on the existence of periodic solutions of differential equations. In recent years, Hopf bifurcation has been intensely studied to investigate the dynamics of differential equations in various fields [6,32,43,47]. Normal forms theory provides powerful tools to bifurcation analysis, whose basic idea is to employ near-identity nonlinear transformations that lead the original system to a qualitatively equivalent differential equation with the simplest form.
In a reaction-diffusion system, especially reaction-diffusion system with delays, there may exist several Hopf bifurcation critical points with parameters' varying. Thus, like the case in ordinary differential equations [36], double Hopf bifurcations will be characterized by intersections of two Hopf bifurcation curves in a two-parameter plane. In fact, a double Hopf bifurcation occurs from a critical point at which the Jacobian evaluated there involves two conjugate pairs of pure imaginary eigenvalues. Guckenheimer and Holmes [26] gave the dynamics in the neighborhood of a codimension-two double Hopf point in an ordinary differential equation. When bifurcation parameters are closed to the double Hopf bifurcation point, the system may exhibit rich dynamics, such as periodic oscillations, quasi-periodic oscillations, coexisting of several oscillations, two-or three-dimensional invariant torus, and even chaos [26,36,60]. For delayed systems governed by delay differential equations, due to the fact that delay differential equations are of infinite dimensional, the center manifold method was adopted to reduce delay differential equations into finite-dimensional systems, and then normal forms can be calculated on the center manifold [8,10,17,20,21,30].
Correspondingly, the multiple timescale method has also been successfully applied to analyze codimension-two nonresonant double Hopf bifurcations [39,59]. Research on the double Hopf bifurcation and complex dynamics of delay differential equations is also carried out by these methods we mentioned above [7,9,40,58,62].
In recent years, the analysis of Hopf bifurcation in a reaction-diffusion systems has at-tracted much attention. Based on the method provided in [31,52], Yi et al. [57] carried out Hopf bifurcation analysis in a reaction-diffusion system without delays. By decomposing the equation with the orthonormal Fourier basis corresponding to the eigenvalue problem, an explicit expression is given for several key parameters in determining the properties of bifurcation according to the method given in [31]. Faria [18] investigated the effect of both time delay and diffusion on Hopf bifurcation, and improved the method given in [19,21] to study normal forms with perturbation parameters of Hopf bifurcation in delayed reaction-diffusion equations with Neumann boundary condition. The results in [18] have been further applied to equations with various practical backgrounds, such as the predator-prey model [13] and so on. Recently, Hopf bifurcation analysis of delayed reaction-diffusion equations with Dirichlet boundary condition has also made considerable progress, such as in [24,50,56], the existence of Hopf bifurcation and the stability and direction of bifurcating periodic solutions for some population models with Dirichlet boundary condition were considered; Furthermore, Guo [23], Chen and Yu [12], combining the non-local phenomenon, considered Hopf bifurcation of the delayed reaction-diffusion equations with Dirichlet boundary condition.
Comparing the research on Dirichlet boundary conditions with the research on Neumann boundary conditions, the main difficulty when investigating such models is the existence of positive steady-state solutions and the complexity of the corresponding eigenvalue problem, so practical analysis is usually carried out by combining the method of phase space analysis, topological degree theory or Liapunov-Schmidt reduction theory [35]. In this paper, only homogeneous Neumann boundary condition is considered, as we can assume that zero is an equilibrium of the system.  [14]. Later, Meixner et al. [41] conducted a preliminary analysis of the problem through codimension-two bifurcation analysis. Parallel research was also seen in literature [5] on the predator-prey system. Song et al. [48] investigated the Turing-Hopf bifurcation from the point of view of analysis of bifurcation and normal forms, in which they extended the method given in [18,57], and deduced the normal form with parameters of Turing-Hopf bifurcation. Xu and Wei [53] applied this method to study the Turing-Hopf bifurcation of a predator-prey model. Turing-Hopf bifurcation can be seen as a special case of Hopf-zero bifurcation of dynamic system from the point of view of classification of bifurcation. Results about Hopf-zero bifurcations and normal form analysis in delayed reaction-diffusion equations can be found in An & Jiang [1]. For double Hopf bifurcation, Lewis and Nagata [37] studied the transitions from axisymmetric steady solutions to nonaxisymmetric waves in a Navier-Stokes model of the differentially heated rotating annulus experiment, and an analytical-numerical center manifold reduction is used to analyze the double Hopf bifurcation points that occur at this transition. Besides [37], to our best knowledge, there have not been any results on double Hopf bifurcation analysis and the corresponding normal form calculation in the reaction-diffusion equations with or without delays.
Based on the general method of normal form simplification in ordinary differential equation or delay differential equation, the ratio of two pairs of imaginary roots appearing at the double Hopf bifurcation point is important to determine the final simplest normal form.
When the ratio is not a rational number, the normal form has a universal form as those provided in [7-9, 39, 40, 58, 59, 62]. However, if the ratio is rational, there must exist some additional terms which coudn't be eliminated because of the resonance. Precisely, if we are about to obtain the topological classification of the dynamical behavior near the double Hopf points by analyzing the third order normal form, this only require that the ratio is not m : n for m, n = 1, 2, 3, 4, namely, weak resonance case. For all weakly resonant double Hopf bifurcation, we can treat them in the same method as that used for nonresonant cases.
In this paper, we aim at the development of a method of computing normal forms on center manifold near an equilibrium of a system of a general delayed reaction-diffusion equations with nonresonant or weakly resonant double Hopf singularity. Consider a general reaction-diffusion model with homogeneous Neuman boundary conditions and u t i (θ)(x) = u i (t + θ, x). Define f i (µ, 0, · · · , 0) = 0, i = 1, 2, · · · , n, such that (1) has zero equilibrium, and assume that every f i is at least C 3 .
We will extend the normal form method given by Faria and Magalhães [21,22] to investigating non-resonant or weakly resonant double Hopf bifurcations in reaction-diffusion systems with delays. The first step of this method is to decompose the phase space, and then to decompose the equations we investigate. Then we can calculate the second and third order normal form on the center manifold, and get the normal form for double Hopf bifurcation, which is a four-dimensional system up to the third order with unfolding parameters restricted to the center manifold. By cylindrical coordinate transformation, the four-dimensional system can be reduced to a two-dimensional amplitude system, in which the unfolding parameters can be expressed by the parameters in the original system. Due to the unfolding analysis given in [26], we can determine the type of twelve distinct kinds of unfoldings by the coefficients we calculated in the normal form, and get corresponding bifurcation sets. By detailed analysis of bifurcation sets, all the dynamical behaviors near the double Hopf bifurcation point can be figured out. In this paper, a series of explicit calculation formulae of the normal form for non-resonant or weakly resonant double Hopf bifurcation is given, with three different cases of wave number: n 1 = 0, n 2 = 0, n 1 = 0, n 2 = 0, and n 1 = 0, n 2 = 0. The wave number is related to the spatial profile of bifurcating periodic solutions.
To illustrate the calculation process of the normal forms, we perform two examples. The first one a diffusive, age-structured epidemic model with two delays. Xiao and Chen [55] proposed an SIS epidemic model with stage structure and time delays in the following form where S(t) and I(t) are the population of susceptible and infected immature individuals, y(t) represents the population of mature individuals, τ is the maturity delay, and ω is freely moving delay. Hopf bifurcation results for such a model are given in Du et al. [15].
Adding random diffusion of individuals into (2), we have a diffusive epidemic model, which is proposed in section 4. We take the time delay ω and a diffusion coefficient as bifurcation parameters. We show that double Hopf bifurcation can appear with two different cases wave numbers: n 1 = 0, n 2 = 0, or n 1 = 0, n 2 = 0. Thus, the spatio-temporal dynamics turn out to be very complicated near the double Hopf bifurcation point. In some regions, there are two stable and spatially nonhomogeneous periodic solutions or a homogeneous and a nonhomogeneous periodic solution coexisting.
The second example we considered is a predator-prey system with time delay. Predatorprey systems have been widely investigated [27,28,54]. Here, we choose a simple system which was studied in [47] to illustrate the procedure of normal form calculation Incorporating diffusion terms into (3), the system becomes a reaction-diffusion system, which is proposed in section 5. Taking r 1 and τ as bifurcation parameters, we will show that such a simple system will exhibit complicated dynamical behavior near the double Hopf bifurcation point such as the existence of quasi-periodic solution on a 2-torus and quasi-periodic solution on a 3-torus. Generally, a vanishing 3-torus might accompany strange attractors, and lead a system into chaos, which is called the "Ruelle-Takens-Newhouse" scenario [4,16,46]. Thus a strange attractor is found near the double Hopf bifurcation point in this prey-predator system.
The paper is organized as follows. In section 2, we discuss eigenfunctions and decomposition of the phase space. Section 3 is devoted to calculating the normal forms of nonresonant or weakly resonant double Hopf bifurcation of (1), and explicit formula of normal form truncated to the third order are obtained in three cases of different wave numbers.
As applications of the method, a diffusive epidemic model with two delays and a diffusive predator-prey system with a delay are considered in sections 4 and 5, respectively, where the conditions for the existence of double Hopf bifurcation are obtained, the normal forms are calculated using the method and formula in section 3. Two-parameter unfoldings and bifurcation diagrams near the critical point are analyzed, and some numerical simulations about the rich dynamics near these bifurcation points are demonstrated. Finally, we close the paper with some conclusions.

EIGENFUNCTIONS AND DECOMPOSITION OF THE PHASE SPACE
Normal form analysis usually depends on a decomposition of the corresponding phase space [36,51]. In case of reaction-diffusion model with time delay, it also requires the decomposition with respect to spatial variables which can be accomplished by finding the eigenfunctions of Laplacian operator. Thus, we first rewrite system (1) into an abstract ordinary differential equation in an appropriate phase space.

Eigenfunctions
Using the general approach to put an partial differential equations into an abstract form introduced in [52], we define the real-valued Hilbert space Since we are about to deal with the double Hopf bifurcations with two pairs of purely imaginary eigenvalues, we usually extend the space X into the corresponding complexification space of X in a natural way by with the general complex-value L 2 inner product Now we can rewrite system (1) into an abstract form as where The linearized equation of (4) at (0,0,...,0) is It is well known that the eigenvalue problem Define B m the subspace of C, by For simplification of notations, in the following we write The characteristic equation associated with (5) is where e λ· (θ)y = e λθ y, for θ ∈ [−r, 0]. By using the Fourier expansion we find that the characteristic equation (6) is equivalent to a sequence of characteristic To consider the double Hopf bifurcation, we consider the case with two parameters, i.e., p = 2 and µ = (µ 1 , µ 2 ) T , and assume that for some µ 0 ∈ R 2 , the following conditions hold [36]: (H 1 ) There exists a neighborhood O of µ 0 such that for µ ∈ O, (5) has two pairs of complex simple conjugate eigenvalues α 1 (µ) ± iω 1 (µ) and α 2 (µ) ± iω 2 (µ), all continuously differential in µ with α 1 (µ 0 ) = 0 , ω 1 (µ 0 ) = ω 1 > 0, α 2 (µ 0 ) = 0 , ω 2 (µ 0 ) = ω 2 > 0, and all other eigenvalues of (5) have non-zero real part for µ ∈ O.

Decomposition of the phase space
In order to adapt the center manifold reduction technique, we have to operate on an equation about U t , which requires an enlarged phase space BC defined by Rewrite (8) as an abstract ordinary differential equation on BC [22]: where A is the infinitesimal generator of the C 0 -semigroup of solution maps of the linear equation (5), defined by with dom(A) = {ϕ ∈ C :φ ∈ C, ϕ(0) ∈ dom(∆)} and X 0 is given by Then on B m , the linear equation is equivalent to the retarded functional differential equation on C n : For the n k (k = 1, 2) mentioned in (H 3 ), define functions of bounded variation η k ∈ Let A k (k = 1, 2) denote the infinitesimal generator of the semigroup defined by (11) with µ = µ 0 , m = n 1 , n 2 , and A * k denote the formal adjoint of A k under the bilinear form be the basis of the generalized eigenspace of A k , A * k corresponding to the eigenvalues iω 1 and iω 2 , respectively, and satisfying , and Ψ(s) = (Ψ 1 (s), Ψ 2 (s)) T . Now, we can decompose BC into a center subspace and its orthocomplement, i.e., where π : BC → P is the projection defined by n k , · · · β (n) n k for c = (c 1 , c 2 , · · · , c n ) T ∈ C. Decompose U t ∈ C 1 0 according to (12) as where z k (t) = (Ψ k , U t , β n k ) k , and y t ∈ C 1 Kerπ := Q 1 for any t. Then in BC the system (9) is equivalent to the systeṁ where z = ( z 1 , z 2 ) T := (z 1 , z 2 , z 3 , z 4 ) T , B = diag(iω 1 , −iω 1 , iω 2 , −iω 2 ), and A 1 is the restric-

CENTER MANIFOLD REDUCTION AND NORMAL FORM
In the previous section, we have calculated the basis of the center space and defined a projection onto it. Applying the center manifold theory [38,52], we know there exists an invariant local center manifold of the equilibrium, which will be calculated approximately and be used to obtain the normal form.

Center manifold reduction
Consider the formal Taylor expansion where F j is the jth Fréchet derivation of F . Then (14) can be written aṡ where z = (z 1 , z 2 , z 3 , z 4 ) ∈ C 4 , y ∈ Q 1 , and f j = (f 1 j , f 2 j ), j ≥ 2, are defined by Let us introduce the following notations as those in [22]: for a normed space Y , denotes the space of homogeneous polynomials of degree j in N +p variables z = (z 1 , · · · , z N ), According to [22], the normal forms are obtained by a recursive transformations of variables of the form (z, y, α) = ( z, y, α) . This recursive process transforms (15) into the where g j = (g 1 j , g 2 j ), j ≥ 2, are the new terms of order j, given by with P is the projection operator, and f j = (f 1 j , f 2 j ) denote the terms of order j in (z, y) obtained after the computation of normal forms up to order j − 1.

Third order normal form
To find the third order normal form of the double Hopf singularity, neglecting the high order terms of the perturbation parameters and noticing the assumption (H 2 ), we first have According to [22], the normal form up to the third order is given bẏ given by (19). After the calculation of the second order normal forms, we obtain the third order f 1 3 (z, 0, 0) given by (27). Noticing that g 1 2 (z, 0, 0) = 0, we divide the computation of the third term g 1 3 (z, 0, 0) = Proj Im(M 1 3 ) c f 1 3 (z, 0, 0) into the following three parts which will be calculated in the following part.

Remark 4.
In section 4, case Ib appears, thus we draw bifurcation set and phase portraits for the unfolding of case Ib in Figure 1 [26]. For case Ib, near the double Hopf bifurcation point, the (α 1 , α 2 ) plane is divided into six regions D1-D6. In region D1, the equilibrium is

APPLICATION TO A DIFFUSIVE EPIDEMIC MODEL WITH DELAY AND STAGE STRUCTURE
In this section, a diffusive epidemic model with delay and stage structure is considered.
Taking the time delay ω and a diffusive coefficient as bifurcation parameters, we show that double Hopf bifurcation can undergo with the two wave numbers in different cases, such as n 1 = 0, n 2 = 0, and n 1 = 0, n 2 = 0. Following the steps of Remark 3, the normal form can be calculated, and the unfolding system can be get. Simulations illustrate that the spatiotemporal dynamics turn out to be very complicated near the double Hopf bifurcation point.
In some regions, there are two stable nonhomogeneous periodic solutions or a homogeneous and a nonhomogeneous periodic solution coexisting.

Model and the existence of double Hopf bifurcation
The stage-structured epidemic model with the maturation delay and freely-moving delay is given by where S(x, t) and I(x, t) represent the population of susceptible and infected immature individuals, and y(x, t) represents the population of mature individuals. τ is the time taken for the immature individuals to maturity, ω is the time taken for the immature individuals from birth to moving freely. α is the natural birth rate, d is the natural death rate of the immature stage, β is the death rate of the mature individuals of logistic nature, µ is the disease transmission rate, and γ is the recovery rate.
By (66), we can draw the curves of Hopf bifurcation values when d 2 varies.
For the double Hopf bifurcation points HH1 and HH3, we can calculate the parameters in Table 2. Obviously, the two points are also of type Ib. We notice that during the calculation of step 3 in third order normal form of the three double Hopf bifurcation points, the cases n 1 = 0, n 2 = 0 (e.g. HH3), and n 1 = 0, n 2 = 0 (e.g. HH1, HH2) both occur.

APPLICATION TO A DIFFUSIVE PREDATOR-PREY MODEL WITH DELAY
In this section, a diffusive predator-prey model with delay is considered. Taking r 1 and τ as bifurcation parameters, we can find the double Hopf bifurcation point. By the steps in Remark 3, we can obtain the normal form and get the unfolding system. We will show that the system will exhibit complex dynamical behavior near the double Hopf bifurcation point: the existence of quasi-periodic solution on a 2-torus, quasi-periodic solution on a 3-torus, and even chaos.

Model and the existence of double-Hopf bifurcation
Let us consider a predator-prey system x ∈ [0, lπ] (72) equipping with homogeneous Neumann boundary condition, where X(x, t) and Y (x, t) represent the densities of prey and predator populations at time t and location x, respectively, d 1 , d 2 > 0 denote the diffusion coefficients of prey and predator, τ is the generation time of the prey species, and r i , a ij (i, j = 1, 2) are positive constants.
By the tanslation u = X − X * , v = Y − Y * , (72) can be written as The linearization of (73) at the origin is , whose characteristic equation is with n ∈ N 0 , Clearly, λ = 0 is not a root of (74), which excludes the existence of Turing bifurcation.
When τ = 0, Eq. (74) becomes the following sequence of quadratic polynomial equations where A n + C = d 1 l 2 )a 11 X * > 0. We know that all roots of Eq. (75) have negative real part.
We would like to seek critical values of τ such that there exists a pair of imaginary Then, we obtian Suppose that holds, Eq. (76) has two positive roots given by By calculation and the results in [47], we can verify the transversality condition Theorem 3. Suppose that (H 5 ) holds, system (72) undergoes a Hopf bifurcation at the origin when τ = τ j− n or τ j+ n .
We fix the other parameters, and choose r 1 and τ as double Hopf bifurcation parameters.
Theorem 4. Suppose that (H 5 ) holds, and there exists r * 1 such that τ j 1 − n 1 = τ j 2 + n 2 . Then system (72) undergoes a double Hopf bifurcation at the origin when

Simulations
In this section, we fix r 2 = 1, a 11 = 1, a 12 = 1.2, a 21 = 2.8, a 22 = 1, d 1 = 0.1, d 2 = 0.2, l = 3, and take r 1 and τ as bifurcation parameters. We can find that when 0.316 < r 1 < 1.47, (H 5 ) holds for all n = 0, and when 0.316 < r 1 < 1.19, (H 5 ) holds for all n = 1. By (77), we can draw the curves of Hopf bifurcation values when r 1 varies, which is shown in Figure 5 Figure 2. System (72) has a quasi-periodic solution on a 2-torus and a quasi-periodic solution on a 3-torus, which will be eliminated by a saddle connection orbit. According to the "Ruelle-Takens-Newhouse" scenario to chaos [4,16,46], a vanishing 3-torus might accompany strange attractors, and lead a system into chaos. Thus, we know near point HH, there may exist a strange attractor, which will be shown numerically in Figure 6.
In Figure 6, we show the existence of quasi-periodic solutions with multiple frequencies, and show the results through the Poincaré map on a Poincaré section. Since all periodic or quasi-periodic solutions are spatially homogeneous, we choose the solution curve of (X(0, t), Y (0, t)), i.e., x = 0. Moreover, because the periodic solutions oscillate in an infinite dimensional phase space [52], we give simulations near the double Hopf bifurcation point on the Poincaré section X(0, t − τ ) = X * . We can see that the system exhibits rich dynamical behavior near the bifurcation point HH in Figure 6. System (72) has a quasi-periodic solution on a 2-torus which becomes a quasi-periodic solution on a 3-torus and then vanishes through a heteroclinic orbit, which is shown in Figure 6 (a) and (b) respectively. After the vanishing of 3-torus, system (72) has a strange attractor, and exhibits chaotic behavior, which is shown in 6 (c).

CONCLUSION REMARKS
In this paper, we have extended the center manifold reduction and normal form method to analyze the dynamical behavior near the double Hopf bifurcation point in a delayed reaction-diffusion system. The method appears to be quite complicated as proceeded in this paper, but it is still an explicit algorithm and is not difficult to be implanted into a computer program. What we should care more about in this method is to calculate the double Hopf point and all the rest part can be exactly obtained by following the procedure we have given.
After normal form derivation with respect to a delayed reaction-diffusion system, the dynamics near the bifurcation point, like the case in ordinary differential equations, is also governed by twelve distinct kinds of unfolding systems, and the bifurcation set about ordinary differential equations for each of the twelve types of unfoldings still applies to our model. To show that our method is a powerful method for analyzing local behavior around a double Hopf bifurcation point, we give two examples: in the stage-structured epidemic model, we show that two spatially inhomogeneous periodic oscillations coexist near the singular point; in the predator-prey system, we show the existence of quasi-periodic solution on a 2-torus, quasi-periodic solution on a 3-torus, and strange attractor appears near the bifurcation point.