Free Vibration Analysis of a Graphene-Platelet-Reinforced, Porous, Two-Cylindrical-Panel System

In this study, a novel, dynamic model of a graphene-platelet-reinforced, porous (GPLRP) double-cylindrical-panel system is proposed. The material properties of a graphene-platelet-reinforced, porous, double-cylindrical-panel system were determined by the Halpin–Tsai micromechanics model and the typical mechanical properties of open-cell metal foams. Different types of porosity distribution and graphene platelet (GPL) distribution patterns were considered. Love’s shell theory was utilized to derive the theoretical formulation, and the Rayleigh–Ritz method was used to calculate the natural frequencies of the system. The proposed model was validated by several comparison studies with the natural frequencies in the existing literature. Finally, the effects of stiffness of Winkler springs, boundary condition, porosity coefficient, porosity distribution, GPL distribution pattern, and GPL weight fraction on the free vibration characteristics of the system were evaluated.


Introduction
The two-cylindrical-panel system, an important component of modern lightweight structures, is widely used in mechanical, vehicle, and marine structures [1,2]. A twocylindrical-panel system may possess ideal structural properties, such as better vibration attenuation and lighter weight, than a single cylindrical panel [3,4]. Therefore, the in-depth study of such systems has great theoretical and practical significance.
Recently, research on the vibration characteristics of multi-beam systems and multiplate systems has been of great interest [5][6][7][8][9][10][11][12]. Many experts have conducted a large number of studies on the vibration characteristics of multi-beam systems [13][14][15][16][17][18]. Kim et al. [19,20] established a double-beam system and analyzed the free vibration characteristics of a system. Deng et al. [21] studied the vibration characteristics of a double-functionallygraded (FG) beam system on an elastic foundation. Hao et al. [22] utilized the modified Fourier-Ritz method to analyze the free vibration characteristics of a double-beam system. Rahman and Lee [23] developed a novel harmonic balance method to analyze the nonlinear vibration characteristics of a double-beam system. Kelly and Srinivas [24] established a system of multiple beams connected by elastic layers and investigated the free vibration of the system. Han et al. [25] discussed the vibration characteristics of a double-beam system by using the improved Wittrick-Williams algorithm.
Many researchers have also conducted numerous studies on the vibration characteristics of a multi-plate system. Oniszczuk [26] investigated the free vibration of two plates connected by a Winkler elastic layer. Based on the Rayleigh-Ritz method, Jeong and Kang [27] developed a novel theoretical method to study the vibration characteristics of multiple rectangular plates coupled by a liquid. Hedrih [28,29] investigated the free vibration characteristics of double-plate systems. Stojanović et al. [30] revealed the vibration and stability characteristics of a system of multiple rectangular plates connected with elastic layers by using an analytical method. The above studies on multi-plate structures are limited to rectangular plates, though many scholars have also carried out detailed studies on circular plates. Hedrih and Simonović [31][32][33] presented a dynamic model of the double circular plate system and investigated the non-linear characteristics of the system. Noga [34] established a dynamic model of double annular and circular membranes connected with the Winkler elastic layer and studied the free vibration of the system by using analytical and numerical methods.
Metal foam materials are widely employed in energy-absorbing systems due to their energy-absorbing capability [35] and other particular characteristics [36]. Keleshteri and Jelovica [37] investigated the free vibration characteristics of the FG porous cylindrical panels. Wang and Wu [38] analyzed the influences of the porous coefficient on the free vibration of an FG porous cylindrical shell.
Due to their pores, porous metal structures are light, but they also can reduce the stiffness of such a structure, so filling materials are needed to enhance the structural stiffness. GPLs possess good mechanical properties, which is why they are employed as reinforcement materials and widely used in various composite structures [39,40]. Using the GPLRP material, Yang et al. [41] analyzed the buckling and free vibration characteristics of a plate, and Wang et al. [42] studied the nonlinear free vibration characteristics of circular cylindrical shells. Chai and Wang [43] studied the traveling wave vibration characteristics of GPLRP-joined conical-cylindrical shells in a spinning motion. Ye and Wang [44] analyzed the nonlinear forced vibration of GPLRP cylindrical shells. Xu et al. [45] investigated the free vibration of a rotating GPLRP beam by using the differential transformation method. Teng and Wang [46] studied the primary, superharmonic, and subharmonic resonances of GPLRP plates. Zhou et al. [47] established an accurate nonlinear buckling analysis of the GPLRP composite cylindrical shells under an axial compressive load. Twinkle and Jeyaraj [48] analyzed the buckling and vibration characteristics of a GPLRP, cylindrical panel. The above studies indicate that a GPLRP structure has good vibration absorption performance and high strength, so it has good prospects for engineering structures.
The above review shows that most studies mainly focus on the vibration characteristics of plane systems, and there are almost no studies on curved systems. To fill the gap in the studies of curved systems, a novel dynamic model of a graphene-platelet-reinforced, porous (GPLRP) two-cylindrical-panel system based on Love's shell theory is established.
To validate the proposed model, several comparisons of the present results with those from open literature are made. Subsequently, the effects of stiffness of the Winkler springs, boundary condition, porosity coefficient, porosity distribution, GPL distribution pattern, and GPL weight fraction on the free vibration characteristics of the system are investigated using the Rayleigh-Ritz method.

Description of the Model
The schematic diagram of the GPLRP, two-cylindrical-panel system with the radius of curvature R, thickness h, subtended angle θ 0 , axial length L, and arc length S (S = Rθ 0 ) coupled by a polymer matrix is shown in Figure 1. As shown in Figure 1, the local cylindrical coordinate system (x i , θ i , z i ) of the ith (i = 1, 2) GPLRP, cylindrical panel is established on the midplane of the GPLRP, cylindrical panel, the origin o i of which is located at the midpoint on the left side of the midplane. The displacements of the arbitrary point of the midplane are denoted by u i , v i , and w i along with the directions x i , θ i , and z i , respectively. The polymer matrix between the GPLRP double cylindrical panels can be considered equivalent to continuously distributed Winkler springs with stiffness K. Both of the GPLRP, cylindrical panels have identical material parameters, geometric parameters, and boundary conditions. of the GPLRP, cylindrical panels have identical material parameters, geometric parameters, and boundary conditions.  Figure 2, three types of porosity distribution along the thickness direction, denoted by porosity-Ⅰ, porosity-Ⅱ, and porosity-Ⅲ , are considered in this study. For different porosity distributions, the Young's modulus, mass density, and Poisson's ratio along the thickness direction of GPLRP, cylindrical panel can be given as follows [41]: in which E*, ρ*, and μ* represent Young's modulus, mass density, and Poisson's ratio of the GPL-reinforced cylindrical panel without pores, respectively. e1, e2, and e3 are the porosity coefficients; em1, em2, and em3 are the mass density coefficients.  As shown in Figure 2, three types of porosity distribution along the thickness direction, denoted by porosity-I, porosity-II, and porosity-III, are considered in this study. For different porosity distributions, the Young's modulus, mass density, and Poisson's ratio along the thickness direction of GPLRP, cylindrical panel can be given as follows [41]: in which E*, ρ*, and µ* represent Young's modulus, mass density, and Poisson's ratio of the GPL-reinforced cylindrical panel without pores, respectively. e 1 , e 2 , and e 3 are the porosity coefficients; e m1 , e m2 , and e m3 are the mass density coefficients. As shown in Figure 2, three types of porosity distribution along the thickness direction, denoted by porosity-Ⅰ, porosity-Ⅱ, and porosity-Ⅲ, are considered in this study. For different porosity distributions, the Young's modulus, mass density, and Poisson's ratio along the thickness direction of GPLRP, cylindrical panel can be given as follows [41]: in which E*, ρ*, and μ* represent Young's modulus, mass density, and Poisson's ratio of the GPL-reinforced cylindrical panel without pores, respectively. e1, e2, and e3 are the porosity coefficients; em1, em2, and em3 are the mass density coefficients.  According to the Halpin-Tsai micromechanics model, the effective Young's modulus of the GPL-reinforced cylindrical panel without porosities can be given as follows [41]: in which where E GPL and E m represent Young's moduli of the GPLs and metal matrix; l GPL , w GPL , and t GPL represent the average length, width, and thickness of the GPLs, respectively; V GPL (z) is the volume fraction of the GPLs. The mass density ρ* and Poisson's ratio µ* of the GPL reinforced cylindrical panel without porosities can be defined by [41] in which ρ GPL and µ GPL represent the mass density and Poisson's ratio of the GPLs, respectively; ρ m and µ m denote the mass density and Poisson's ratio of the metal matrix, respectively. The relationship between the modulus of elasticity and mass density of the open-cell metal foams is as follows [41]: The relationship between the mass density coefficients and the porosity coefficients can be obtained by substituting Equations (1) and (2) into Equation (9): To obtain a valuable comparison of different combinations, the total mass of the GPLRP, cylindrical panel was set to be the same for various types of porosity and GPL distribution, resulting in in which e 2 and e 3 can be determined by e 1 .
To reinforce the porous cylindrical panel, the GPLs are applied as fillers inside the porous materials. The three GPL distribution patterns in the thickness direction of the cylindrical panel are shown in Figure 3. The volume fraction V GPL (z) of the GPLs can be written as follows [41]: in which s j1 , s j2 , and s j3 are the coefficients of the GPL volume fraction for different porosity distributions and GPL distribution patterns, and j = 1, 2, 3 correspond to the Porosity-I, -II, and -III distributions, respectively. The s j1 , s j2 , and s j3 can be given as follows [41]: (13) where V T GPL is the total volume fraction of the GPLs and given by [41]  where W GPL is the weight fraction obtained by the ratio of the mass of the GPL to the mass of the cylindrical panel.
is the total volume fraction of the GPLs and given by [41] GPL GPL GPL GPL m GPL where WGPL is the weight fraction obtained by the ratio of the mass of the GPL to the mass of the cylindrical panel.

Governing Equations and Solution
According to Love's shell theory [49], the strains of the arbitrary point at the distance z from the midplane can be written as where,

Governing Equations and Solution
According to Love's shell theory [49], the strains of the arbitrary point at the distance z from the midplane can be written as where, ε The stress-strain relationship of the ith GPLRP, cylindrical panel can be defined as [50] The force and moment resultants of the ith GPLRP, cylindrical panel can be obtained by integrating the stresses along with the z-direction of the GPLRP, cylindrical panel, which are expressed as  By substituting Equations (15) and (18) into Equations (20) and (21), we get: where ε T (i) and N T (i) are given by and S (i) is a stiffness matrix defined by pq (p, q = 1, 2, 6) are stretching, coupling, and bending stiffness coefficients and can be expressed as The strain energy and kinetic energy of the GPLRP, two-cylindrical-panel system can be written as where The potential energy stored in the polymer matrix between the two-cylindrical-panel system can be written as where K is the stiffness of Winkler springs. For generality and convenience, the dimensionless parameters x i and θ i are introduced Substituting Equation (30) into Equations (16) and (17) yields By separating variables, the admissible displacement functions of the ith GPLRP, cylindrical panel can be given by where t is the time variable; the ω is the natural angular frequency; The displacement amplitude functions of the ith GPLRP, cylindrical panel expressed by the Chebyshev polynomials can be given as [41] where U where P m (x i ) and P n (θ i ) are a set of orthogonal polynomials in the interval [−1,1]. In addition, in Equation (34), functions f α (x i , θ i ), α = u, v, w are the auxiliary functions related to boundary conditions.
in which p, q, r, and s are equal to zero, one, or two, for different boundary conditions. For instance, in the GPLRP, two-cylindrical-panel system with clamped boundaries (C) on x = 0 and x = L, free (F) at θ = −θ 0 /2, and simply supported (S) at θ = θ 0 /2, auxiliary functions can be written as The strain energy of the GPLRP, two-cylindrical-panel system can be obtained by substituting Equations (24) The kinetic energy of the GPLRP, two-cylindrical-panel system can be obtained by substituting Equation (33) into Equation (27): The potential energy of the polymer matrix can be acquired by substituting Equation (33) into Equation (29): To determine the solutions, the Rayleigh-Ritz method is employed. The Lagrangian energy function Π of the GPLRP, two-cylindrical-panel system can be expressed as (41) in which U max and T max are the maximum strain energy and maximum kinetic energy of the GPLRP, two-cylindrical-panel system; U Smax is the maximum potential energy of the polymer matrix. Then, the Lagrangian energy function Π is minimized with respect to the unknown coefficient, yielding The eigenvalue equations derived from Equation (42) can be summarized in the matrix form as K − ω 2 M X = 0 (43) in which M and K denote the mass and stiffness matrices, respectively. The explicit forms of the modal vector X can be written as follows: For convenience, the dimensionless natural frequency of the GPLRP, two-cylindricalpanel system is defined by

Model Validation and Convergence Analysis
To prove the validity of the present formulations, the dimensionless fundamental frequency ω = ωL ρ m (1 − µ 2 m )/E m of a single GPLRP plate under simply supported boundary conditions is compared with the results in the existing literature [41]. The material parameters used are consistent with the literature [41]. The geometric parameters of a single GPLRP cylindrical panel are as follows: h = 0.1 m, L/S = 1, L/h = 20, W GPL = 1.0%, e 1 = 0.5, i = 1, K = 0, and R = ∞. As shown in Table 1, it can be clearly seen that the results of this paper are well consistent with the literature [41]. To further clarify the correctness of the derivation of this study, Table 2 shows the comparison of the natural angular frequencies ω (rad/s) in this article with existing literature for the double-plate system under the simply supported boundary condition [26]. The geometric and material parameters are as follows: L = 1 m, S = 2 m, E = 1 × 10 10 N/m 2 , h = 1 × 10 −2 m, i = 1, 2, K = 6 × 10 4 N/m 3 , µ = 0.3, ρ = 5 × 10 3 kg/m 3 , and R = ∞. As can be observed in Table 2, very good agreement was obtained in this comparison study. The above comparisons with the open literature are sufficient to prove the accuracy of the derivation in this paper. In this study, the Rayleigh-Ritz method was applied to solve the natural frequencies of the GPLRP, two-cylindrical-panel system. The accuracy of this method depends on the truncation numbers (M, N) of the Chebyshev polynomials. Therefore, it is necessary to discuss the influence of truncation orders M and N on the calculation results. The simulation parameters are given in Section 4. As shown in Table 3, when the truncation numbers M = N = 11, the dimensionless natural frequencies ω = ωR ρ m /E m obtained good convergence. Therefore, the truncation numbers M and N were set to 11 in the following study. Table 3. Convergence of the first six dimensionless natural frequencies of the GPLRP, two-cylindricalpanel system (Porosity-I, GPL A, e 1 = 0.5).

Results and Discussion
In this section, the material and geometric parameters of the GPLRP, two-cylindricalpanel system in this paper used are as follows: E m = 130 GPa, ρ m = 8960 kg/m 3 , µ m = 0.34, E GPL = 1.01 Tpa, ρ GPL = 1062.5 kg/m 3 , µ GPL = 0.186, t GPL = 1.5 nm, l GPL = 2.5 µm, w GPL = 1.5 µm, W GPL = 1.0%, e 1 = 0.5, R = 1 m, h = 0.01R, L = R, θ 0 = π/6, and K = 1 × 10 14 N/m 3 . Figure 4 shows the variation un dimensionless natural frequencies of the system with the spring stiffness K, considering different boundary conditions. An interesting phenomenon is that the variation of spring stiffness has no influence on the first-order natural frequency of the system. In addition, there is an elastic domain for other natural frequencies except for the first-order natural frequency, in which the natural frequency increases rapidly with the increase in the spring stiffness K. This finding is further confirmed by the results in Table 4. The reason for this phenomenon is that the natural frequency of the GPLRP, two-cylindrical-panel system can be divided into two parts, one of which is the same as the natural frequency of a single GPLRP, cylindrical panel, and the other is identical to the natural frequency of a single GPLRP, cylindrical panel on an elastic foundation with the stiffness 2K [26]. Therefore, except for the first-order natural frequency, other natural frequencies can be adjusted by adjusting the value of the spring stiffness to obtain the desired natural frequency. This means that the stiffness of the system can be changed by adjusting the spring stiffness. For the sake of brevity, only the fundamental frequency of the system is discussed below. frequency, other natural frequencies can be adjusted by adjusting the value of the spring stiffness to obtain the desired natural frequency. This means that the stiffness of the system can be changed by adjusting the spring stiffness. For the sake of brevity, only the fundamental frequency of the system is discussed below.     6 demonstrate the effects of the porosity coefficient on dimensionless fundamental frequencies of the system for specific porosity distributions and boundary conditions. It can be clearly seen that the natural frequency of the system decreases significantly with the increase in the porosity coefficient, except for Porosity-Ⅰ under the CCCC boundary condition. This indicates that in most cases, the greater the porosity coefficient e1, the weaker the stiffness of the system. Then, comparing Figure 5a and Figure  6a, it can also be seen that the natural frequency of the system under the SSSS boundary condition decreases rapidly with the increase in the porosity coefficient. However, the natural frequency of the system under the CCCC boundary condition changes slightly with the increase in the porosity coefficient, which indicates that the combined influence of the porosity coefficient and GPL distribution on the natural frequency is related to the boundary conditions. It is also notable that when e1 ≤ 0.5, for all porosity distributions, GPL A possesses the largest natural frequency and GPL B possesses the smallest natural  Figures 5 and 6 demonstrate the effects of the porosity coefficient on dimensionless fundamental frequencies of the system for specific porosity distributions and boundary conditions. It can be clearly seen that the natural frequency of the system decreases significantly with the increase in the porosity coefficient, except for Porosity-I under the CCCC boundary condition. This indicates that in most cases, the greater the porosity coefficient e 1 , the weaker the stiffness of the system. Then, comparing Figures 5a and 6a, it can also be seen that the natural frequency of the system under the SSSS boundary condition decreases rapidly with the increase in the porosity coefficient. However, the natural frequency of the system under the CCCC boundary condition changes slightly with the increase in the porosity coefficient, which indicates that the combined influence of the porosity coefficient and GPL distribution on the natural frequency is related to the boundary conditions. It is also notable that when e 1 ≤ 0.5, for all porosity distributions, GPL A possesses the largest natural frequency and GPL B possesses the smallest natural frequency. However, when e 1 > 0.5, GPL A does not always possess the largest natural frequency; for example, see Porosity-II under the SSSS boundary condition. It can be concluded that the effect of the porosity coefficient on the natural frequency of the system is associated with the porosity distributions, GPL distribution patterns, and boundary conditions. frequency. However, when e1 > 0.5, GPL A does not always possess the largest natural frequency; for example, see Porosity-ⅠⅠ under the SSSS boundary condition. It can be concluded that the effect of the porosity coefficient on the natural frequency of the system is associated with the porosity distributions, GPL distribution patterns, and boundary conditions. To explain the effect of porosity distribution more clearly, Figures 7 and 8 indicate the dimensionless fundamental frequencies versus the porosity coefficient for different porosity distributions with specific GPL distribution patterns and boundary conditions. As can be seen in Figures 7 and 8, for all kinds of porosity distributions, Porosity-Ⅰ possesses the largest natural frequency. This indicates that Porosity-Ⅰ can significantly To explain the effect of porosity distribution more clearly, Figures 7 and 8 indicate the dimensionless fundamental frequencies versus the porosity coefficient for different porosity distributions with specific GPL distribution patterns and boundary conditions. As can be seen in Figures 7 and 8, for all kinds of porosity distributions, Porosity-I possesses the largest natural frequency. This indicates that Porosity-I can significantly enhance the stiffness of the system. It is worth noting that the natural frequency of the system decreases significantly with the increase in the porosity coefficient, except for the Porosity-I under the CCCC boundary condition. This further confirms that the effect of the porosity distributions on the natural frequency of the system is associated with the boundary conditions. From the above results, it can be confirmed that the porosity coefficient, porosity distribution, GPL distribution pattern, and the boundary condition interact with each other and have coupled effect on the free vibration characteristics of the GPLRP, two-cylindrical-panel system.   Figure 9 reveals the variations in dimensionless fundamental frequencies with GPL weight fraction WGPL for different porosity distributions and GPL distribution patterns. It is clearly shown that the natural frequency increases significantly with the increase in the weight fraction WGPL for different porosity distributions and GPL distribution patterns. This means that the GPL filling material can significantly improve the effective stiffness of the system. As shown in Figure 9b, it was interesting to find that all GPL distribution patterns in Porosity-ⅠⅠ have almost the same increase in efficiency with the increase in the GPL weight fraction WGPL. This indicates that the porosity distributions play a dominant role, and it also shows that choosing an appropriate porosity distribution is of great significance to improving the stiffness of the system. This means that the strength of the system can be improved by choosing the appropriate porosity distribution.  It is clearly shown that the natural frequency increases significantly with the increase in the weight fraction W GPL for different porosity distributions and GPL distribution patterns. This means that the GPL filling material can significantly improve the effective stiffness of the system. As shown in Figure 9b, it was interesting to find that all GPL distribution patterns in Porosity-II have almost the same increase in efficiency with the increase in the GPL weight fraction W GPL . This indicates that the porosity distributions play a dominant role, and it also shows that choosing an appropriate porosity distribution is of great significance to improving the stiffness of the system. This means that the strength of the system can be improved by choosing the appropriate porosity distribution.

Conclusions
This study proposed a dynamic model of the GPLRP, two-cylindrical-panel system with general boundary conditions based on Love's shell theory to study the free vibration characteristics of the system. The Rayleigh-Ritz method is used to calculate the solutions. From the results, some valuable conclusions can be drawn as follows: (1) The variation in the stiffness of the Winkler springs in the GPLRP, two-cylindricalpanel system has no effect on the first-order natural frequency, but other natural frequencies can be adjusted by controlling the value of the spring stiffness to obtain the desired natural frequency. (2) The porosity coefficient, porosity distribution, GPL distribution pattern, and boundary condition interact with each other and have a coupled influence on the vibration characteristics of the GPLRP, two-cylindrical-panel system. (3) Increasing the proportion of the GPL filling material can significantly improve the stiffness of the GPLRP, two-cylindrical-panel system, thereby increasing the natural frequency of the system.
Funding: This research was supported by the National Natural Science Foundation of China (grant number 12002080).

Conflicts of Interest:
The authors declare no conflict of interest.

Conclusions
This study proposed a dynamic model of the GPLRP, two-cylindrical-panel system with general boundary conditions based on Love's shell theory to study the free vibration characteristics of the system. The Rayleigh-Ritz method is used to calculate the solutions. From the results, some valuable conclusions can be drawn as follows: (1) The variation in the stiffness of the Winkler springs in the GPLRP, two-cylindricalpanel system has no effect on the first-order natural frequency, but other natural frequencies can be adjusted by controlling the value of the spring stiffness to obtain the desired natural frequency. (2) The porosity coefficient, porosity distribution, GPL distribution pattern, and boundary condition interact with each other and have a coupled influence on the vibration characteristics of the GPLRP, two-cylindrical-panel system. (3) Increasing the proportion of the GPL filling material can significantly improve the stiffness of the GPLRP, two-cylindrical-panel system, thereby increasing the natural frequency of the system.
Funding: This research was supported by the National Natural Science Foundation of China (grant number 12002080).