Multipulse Homoclinic Orbits and Chaotic Dynamics of a Reinforced Composite Plate with Carbon Nanotubes

Faculty of Mathematics and Physics, Huaiyin Institute of Technology, Huaian 223003, China Department of Mathematics, Nanjing University of Aeronautics and Astronautics, Nanjing 210016, China College of Mathematics and Systems Science, Shandong University of Science and Technology, Qingdao 266590, China School of Mathematics and Physics, Yancheng Institute of Technology, Yancheng 224051, China School of Mathematics and Physics, Anhui University of Technology, Ma’anshan 243002, China


Introduction
Carbon nanotubes (CNTs), as a new type of advanced materials, have attracted a lot of attention of researchers.
is is because CNTs possess high strength and stiffness with high aspect ratio and low density. Due to these properties of CNTs, a number of studies on the nonlinear vibrations and dynamic responses of CNT-reinforced composite plates have been carried out by many researchers in recent years. By employing an equivalent continuum model, Formica et al. [1] investigated the vibration behaviors of CNTreinforced composite plates. Zhu et al. [2] presented the bending and free vibration analysis of CNT-reinforced composite plates using the finite element method based on the first-order shear deformation plate theory. In their work, the authors showed the effects of the volume fractions of CNTs and the edge-to-thickness ratios on the bending responses, natural frequencies, and mode shapes of the plates. In addition, Wang and Shen [3] examined the nonlinear dynamic response of CNT-reinforced composite plates resting on elastic foundations in thermal environments. e motion equations were derived based on a higher-order shear deformation theory with a von Kármántype of kinematic nonlinearity. In [4], a mixed Navierlayerwise differential quadrature method was employed by Malekzadeh and Heydarpour for the free vibration and static response analysis of functionally graded carbon nanotube-(FG-CNT-) reinforced composite laminated plates. Based on the conventional Ritz method accompanied with the Lagrangian multiplier technique, Kiani [5] analyzed the free vibration characteristics of FG-CNTreinforced composite plates located on point supports. Rafiee et al. [6] applied Galerkin's method to deal with the nonlinear dynamic stability of initially imperfect piezoelectric FG-CNT-reinforced composite plates under combined thermal and electrical loadings. Using the Fourier series expansion and state-space technique, Alibeigloo and Liew [7] investigated the bending behavior of FG-CNT-reinforced composite plates with simply supported edges subjected to thermomechanical loads. Subsequently, this work was extended by Alibeigloo and Emtehani [8] for various boundary conditions by employing the differential quadrature method. Recently, Zhang et al. [9] carried out the analysis of geometrically nonlinear large deformation of triangular FG-CNT-reinforced composite plates using the element-free improved moving least-squares Ritz (IMLS-Ritz) method.
In addition, with the increasing applications of functionally graded materials (FGM) in modern technology, the buckling and vibration analysis of FGM structures began to attract a widespread attention. us, many studies on the buckling and vibration analysis of FGM shells have been published in the literature. Zhang and Li [10] discussed the dynamic buckling of FGM truncated conical shells subjected to normal impact loads. Bagherizadeh et al. [11] investigated the mechanical buckling of simply-supported FGM cylindrical shells using third-order shear deformation shell theory. ey found that system parameters have great influence on the buckling characteristics of FGM shells. Buckling behavior analysis based on Reddy's high-order shear deformation theory has also been performed by Sun et al. [12] for FGM cylindrical shells subjected to an axial compression in the thermal environment. In recent years, there have been some numerical methods for analyzing the plates and shells. e multiquadric radial basis function (MQ) method was developed and applied by Ferreira [13] to discuss the effects of system parameters on the laminated composite plates. e method of discrete singular convolution (DSC) gives a fast and accurate solution of the mathematical physics and engineering problems. en, Civalek and his coworkers made a number of remarkable studies using the DSC method; such investigations involved the analysis of composite conical shells and panels [14][15][16], plates on elastic foundations [17,18], and so on.
However, in the course of our study, we found that there are only few studies on the global bifurcations and multipulse chaotic dynamics for the CNT-reinforced composite plate. In order to eliminate or suppress large nonlinear vibrations and chaotic motions of the CNTreinforced composite plate, we should deepen and complete the theoretical analysis on the CNT-reinforced composite plate model, discuss the complex dynamic behaviors, explore the existence conditions of the multipulse Shilnikov-type orbits, and analyze the impact of parameters on the system, so as to ensure the stability and controllability of the CNT-reinforced composite plate. e present work is therefore motivated to examine the global bifurcations and multipulse chaotic dynamics of the CNTreinforced composite plate considered by Guo and Zhang [19]. e energy-phase method developed by Haller and Wiggins [20] is applied to study the multipulse homoclinic bifurcations and chaotic dynamics for the CNT-reinforced composite plate. In recent decades, two theories that are usually employed to detect the global bifurcations and multipulse Shilnikov-type chaotic dynamics in high-dimensional systems are the energy-phase method [20] and the extended Melnikov method [21]. Also, some researchers have utilized the developed theories to engineering applications, see, for example, [22][23][24][25][26][27][28][29].
Based on the research by Guo and Zhang [19], the averaged equations are obtained for the case of 1 : 1 internal resonance, principal parametric resonance, and 1/2 subharmonic resonance of the CNT-reinforced composite plate. From the averaged equations, the explicit expressions of the 3-order normal form are derived employing the normal form theory. e energy-phase method is utilized to detect the existence of the multipulse Shilnikov-type homoclinic orbits in the CNT-reinforced composite plate. Homoclinic trees are constructed to illustrate the repeated bifurcations of multipulse solutions. In addition, numerical simulations are given to verify the analytical predictions. e organization of this paper is as follows. In Section 2, we derive the averaged equations and normal form for the CNT-reinforced composite plate. According to the normal form derived in Section 2, the dynamics of the unperturbed system is discussed in Section 3. In Section 4, a detailed dynamic analysis of the perturbed system is presented, and the existence of multipulse homoclinic orbits is confirmed. Numerical simulations are given in Section 5 to show that the chaotic motions can occur in the CNT-reinforced composite plate. A short conclusion is reported in Section 6.

Equations of Motion and Perturbation Analysis.
is paper focuses on studies of the multipulse homoclinic orbits and chaotic dynamics of a four-edge simply-supported CNT-reinforced composite rectangular plate subjected to the in-plane and transverse excitations. e model is shown in Figure 1 [19]. e edge length and width of the plate in the x and y directions are, respectively, a and b, and the thickness is h. A Cartesian coordinate system is located in the middle surface of the CNT-reinforced composite rectangular plate. e displacements of an arbitrary point within the plate are u, v, and w in the x, y, and z directions, respectively. e in-plane excitation along the y direction at x � a is given by P � P 0 − P 1 cos Ω 2 t. e transverse excitation subjecting to the CNT-reinforced composite rectangular plate is of the form F � F 0 − F 1 cos Ω 1 t.
According to the research by Guo and Zhang [19], the axial and transverse displacement fields at any point for the CNT-reinforced composite plate are given as v(x, y, t) � v 0 (x, y, t) + zφ y (x, y, t) − z 3 4 3h 2 φ y + zw 0 zy , where u 0 , v 0 , and w 0 represent the displacements at the midplane of the plate in the x, y, and z directions. φ x and φ y denote the rotations of the transverse normal at the midplane about the x and y axes. Employing the von Kármántype plate theory and in terms of the displacements, the strains ε i (i � xx, yy) and the curvatures c i (i � xy, yz, zx) in the midplane can be expressed as (2) e force and moment resultants associated with the strains and curvatures of the plate constitutive equations are where the coupling matrix B ij and E ij is equal to zero in the stiffness matrix. en, the nonlinear governing equations of motion for the CNT-reinforced composite plate are obtained as follows [19]: Figure 1: e model of a CNT-reinforced composite rectangular plate.
Mathematical Problems in Engineering 3 where μ is the damping coefficient, and all the other coefficients can be found in [19]. e associated boundary conditions for the simply-supported CNT composite plate can be written as Also, the following nondimensional parameters are introduced: Mathematical Problems in Engineering (i � 1, 2, 3, 4, 5, 6; j � 1, 2, 3, 4, 5, 6).

(7)
In this paper, our research is focused on the multipulse global bifurcations and chaotic dynamics of the CNT composite plate in its first two modes. Hence, we express w in the following form: where w 1 and w 2 are the amplitudes of two modes. e other variables and transverse excitation are given as en, substituting these expressions into equations (5a)-(5e) and applying the Galerkin integration procedure, we can obtain the displacements u 0 , v 0 , φ x , and φ y with respect to w 0 . Consequently, the nonlinear ordinary differential equations of this system in terms of transverse displacements are derived. According to the results of convergence studies given by Hao et al. [30], the dimensionless equations of motion for the CNT-reinforced composite plate in the first two modes are derived as [19] € where all the coefficients can be found in [19]. We consider the case of 1 : 1 internal resonance, principal parametric resonance, and 1/2 subharmonic resonance of the CNTreinforced composite plate, for which the resonant relations are listed as follows: where σ 1 and σ 2 are two detuning parameters. For convenience, we let Ω � 2 in the following analysis. e uniform solutions of equations (10a) and (10b) take the following form: where T 0 � t and T 1 � εt. We employ the method of multiple scales to equations (10a) and (10b) and obtain the averaged equations as follows: where f 0 � (1/4)β 1 .

Dynamics of the Unperturbed System
Setting ε � 0 in equations (26a)-(26d) results in the completely integrable equations defined as the unperturbed system. Hence, we will now study the nonlinear dynamics of the following unperturbed system: Consider the first two equations of systems (28a)-(28d): Note that homoclinic bifurcations occur in systems (29a) and (29b) for η < 0. Also, when μ − α 5 I 2 > 0, we can obtain that systems (29a) and (29b) have only the trivial zero solution (z 1 , z 2 ) � (0, 0) being a center. On a curve expressed by μ � α 5 I 2 , that is, three solutions bifurcate from the trivial zero solution through a pitchfork bifurcation denoted by q 0 (I) � (0, 0) and q ± (I) � (B, 0), respectively, where For all I ∈ [I 1 , +∞), systems (29a) and (29b) possess one hyperbolic saddle point q 0 (I) connected to itself by a pair of homoclinic orbits, z h ± (T 1 , I), that is, lim T 1 ⟶ ± ∞ z h ± (T 1 , I) � q 0 (I). erefore, the set defined by is a two-dimensional invariant manifold in the full fourdimensional phase space. Based on the analysis given by Haller and Wiggins [20], we know that the two-dimensional invariant manifold M is normally hyperbolic and contains three-dimensional stable and unstable manifolds. Let W s (M) and W u (M), respectively, denote the stable and unstable manifolds. e existence of the homoclinic orbit in systems (29a) and (29b) to q 0 (I) � (0, 0) means that W s (M) and W u (M) intersect nontransversally along a three-dimensional homoclinic manifold Γ described as follows: where D I H(z, I) � − (zH/zI). Letting ξ 1 � − μ + α 5 I 2 , δ 1 � − η, systems (29a) and (29b) can be rewritten as e Hamiltonian for Hamilton systems (34a) and (34b) is When H � 0, there is a homoclinic loop Γ 0 formed by one saddle point q 0 (I) and a pair of homoclinic orbits z ± (T) 1 . Using equations (34a) and (34b) and (35), the analytical expressions for the homoclinic orbits of (34a) and (34b) are then obtained as We now consider the dynamics of the unperturbed system of (26a) and (26d) restricted to M given by _ I _ ϕ � σ 2 I + α 5 Iq 2 0 (I) + us, the resonant value I r is derived as For I � I r , the phase shift Δϕ can be defined as

Mathematical Problems in Engineering
which will be used in subsequent studies to determine the condition under which the multipulse Shilnikov-type homoclinic orbits may exist. Substituting the first equations of (36a) and (36b) into the fourth equations of (28a)-(28d) leads to After integration, we obtain Accordingly, the phase shift is also obtained as (43)

Dynamics of the Perturbed System
After obtaining the detailed properties of the subspace (z 1 , z 2 ) for systems (28a)-(28d), the next step is to examine the effects of small perturbation terms (0 < ε ≪ 1) on systems (28a)-(28d). e energy-phase method developed by Haller and Wiggins [20] is utilized to determine the existence of multipulse orbits and chaotic dynamics in the CNT-reinforced composite rectangular plate.

e Case of Dissipative Perturbations.
We start by considering the influence of such small perturbations on manifold M. Based on the research by Haller and Wiggins [20], we know that the manifold M is invariant for small nonzero ε. erefore, the perturbed annulus M ε is taken to be the same as M. at is, In order to study the dynamic behaviors of the perturbed vector field restricted to M ε near the resonance I � I r , we now introduce the following scale transformations: e last two equations of (26a)-(26d) then become where the prime represents the differentiation with respect to τ. When ε � 0, systems (46a) and (46b) become

Systems (47a) and (47b) are a Hamiltonian system with the Hamiltonian
A simple calculation indicates that systems (47a) and (47b) have two singular points expressed by It is known that if 6β 3 f 0 I 2 r cos 2ϕ c < 0, the singular point P 0 is a center point. And if 6β 3 f 0 I 2 r cos 2ϕ s > 0, the singular point Q 0 is a saddle connected to itself by a homoclinic orbit. Figure 2(a) shows the phase portrait of systems (47a) and (47b). It is observed that, for sufficiently small parameter ε, the singular point Q 0 still remains a hyperbolic saddle point Q ε , while for small perturbations, the singular point P 0 changes to a hyperbolic sink P ε , see Figure 2(b) for the phase portrait of perturbed systems (46a) and (46b).
On the basis of equation (48), at h � 0, the basin of the attractor for ϕ min can be estimated as follows: Substituting ϕ s in (49a) and (49b) into (50) yields e neighborhood of I r is defined as where the constant C is sufficiently large so that the unperturbed homoclinic orbit is enclosed within the annulus.

e Existence of Multipulse Homoclinic
Orbits. Based on the results obtained by Haller and Wiggins [20], the dissipative energy-difference function is of the form where and A represents the region contained between a pair of homoclinic orbits in the (z 1 , z 2 ) plane, zA l is the boundary of A, and Δϕ denotes the phase difference given in (43). e third term in (53) can be calculated as follows: us, the dissipative energy-difference function becomes where ξ 1 � ξ 1 + (3/4)α 5 I r (I r − 1). A dissipation factor is defined as d � (μ/f 0 ). Equation (56) can then be rewritten in the following form: Δ n H D (ϕ) � f 0 I r sin(2ϕ + nΔϕ)sin(nΔϕ) + 2nξ 1 d 3α 5 Δϕ .
(57) e zeros of Δ n H D (ϕ) can be achieved by solving the following equation: e upper bound of the dissipation factor d is derived from equation (58): According to (59), it is found that, for small dissipative effects d < 1, the following upper bound on the pulse numbers is obtained: which means that the homoclinic tree is finite even for arbitrarily small dissipation. e result in equation (60) also indicates that the upper bound n max is inversely proportional to the dissipation factor d. We now identify zeros of the energy-difference function Δ n H D (ϕ) in the interval [− (π/2), (π/2)]. From (57), we find that, for any n satisfying the zeros of Δ n H D (ϕ) in [− (π/2), (π/2)] are where φ � − (1/2)arcsin(2nξ 1 dΔϕ/3α 5 I r sin(nΔϕ)). Define a set that contains all transverse zeros of Δ n H D (ϕ) as It is easy to know that both zeros are transverse under condition (61). We now need to introduce the two additional angles and utilize them to construct the two sets of transverse zeros of Δ n H D (ϕ) as Consider a domain S 0 ∈ A ε enclosed inside the homoclinic orbit of (47a) and (47b) located in the interval ϕ ∈ [− (π/2), (π/2)]. e periodic orbits in S 0 are classified according to their pulse numbers, see Figure 3. Note that, under condition (61), all internal orbits outside S 0 intersect Z 1 − transversally. us, for any periodic or homoclinic orbit c outside S 0 , the pulse number is N(c) � 1. Define the energy sequence such that h n provides the energy level related to an orbit closer to the center. Define the open set of internal orbits in S 0 as e pulse sequence is then defined as Since the energy of the periodic orbits in S 0 decreases monotonically with the orbits shrinking to the center, we necessarily have e layer sequence is defined as where Int(·) denotes the interior of a set. And the construction of the layer sequence is shown in Figure 3. We can easily see that all these sequences defined above are finite by (60). For any periodic orbit c ∈ S 0 , the pulse number is N(c) � N k . As illustrated in Figure 3, the layer radii are defined as A recursive algorithm for the calculation of the distribution of pulse numbers and the layer radii is implemented on a computer. For different values of the dissipation factor d, we obtain the corresponding pulse diagrams and layer radius diagrams as a function of the phase shift Δϕ for N k < 20, illustrated in Figures 4 and 5. As can be seen in Figure 4, the horizontal line segments at each level N k identify that an infinity of N k -pulse orbits exist for all values of the phase shift in the interval below that line. Besides, the diagrams in Figure 5 exhibit a gradual breakup of the homoclinic tree as d is increased.
is observation means that the system parameters have a vital effect on the distribution of the pulse numbers and the layer radii.
We now detect the existence of the multipulse Shilnikovtype orbits based on eorem 4.5 in [20]. Note that the sinks of (46a) and (46b) become the centers of (47a) and (47b). erefore, systems (47a) and (47b) have a nondegenerate equilibrium e values of (d, Δϕ, N) for which the center P 0 falls in the zero set Z n − given in (63) must be determined, that is, which yields

12
Mathematical Problems in Engineering

Conclusions
In this paper, the multipulse homoclinic orbits and chaotic dynamics of a CNT-reinforced composite plate under combined in-plane and transverse excitations are studied. e method of multiple scales is adopted to acquire the averaged equations in the case of 1 : 1 internal resonance, principal parametric resonance, and 1/2 subharmonic resonance. On the basis of the averaged equations obtained, the normal form theory is employed to derive the expressions of normal form associated with a double zero and a pair of pure imaginary eigenvalues. e energy-phase method proposed by Haller and Wiggins [20] is utilized to detect the presence of the multipulse homoclinic orbits and chaotic dynamics for the resonant case. As the trajectory of motion comes close to the sink point P ε , every Shilnikov-type orbit takes off again and repeats this similar motion in the full four-dimensional phase space and eventually leads to the multipulse Shilnikov-type orbits. Our analysis in Section 4 also indicates that the multipulse Shilnikov-type orbits depend on the system parameters and dissipative perturbations. It is known that the existence of multipulse Shilnikov-type orbits implies the existence of chaos in the sense of Smale horseshoes. Homoclinic trees are presented to describe the repeated bifurcations of multipulse solutions. From the diagrams, we can see a gradual breakup of the homoclinic tree with the increase of the dissipation factor. is observation denotes that the damping coefficient and system parameters affect the distribution of the pulse numbers and the layer radii. e analytical results obtained here are extensions of those appearing in [19].
To confirm the theoretical results, numerical simulations are applied to examine the chaotic dynamics of the CNTreinforced composite plate. e numerical results demonstrate that there exist multipulse homoclinic orbits and chaotic dynamics in systems (13a)-(13d) when the conditions are satisfied. As we all know, under certain conditions, the multipulse chaotic motions in the averaged equations can result in the multipulse amplitude-modulated chaotic vibrations in the original system, which indicates the existence of multipulse chaotic dynamics in the CNT-reinforced composite plate. e research in this paper provides a detailed explanation of the multipulse jumping behaviors for the CNT-reinforced composite plate under combined in-plane and transverse excitations. It is established simultaneously that the damping coefficient and system parameters have a significant influence on the nonlinear dynamics of the CNT-reinforced composite plate. erefore, as mentioned above, the nonlinear dynamic behaviors of the CNT-reinforced composite plate can be controlled by varying the system parameters.

Data Availability
e data used to support the findings of this study are available from the corresponding author upon request.

Conflicts of Interest
e authors declare no conflicts of interest.