Electromechanical free vibration analysis of piezoelectric-integrated cylindrical panels and shells using the two-dimensional GDQ method

Abstract

This study investigates the free vibration behavior of cylindrical panels and shells of revolution integrated with piezoelectric patches. Based on the first-order shear deformation theory and Hamilton’s principle, the mechanical equations of the shell and panel are formulated. The governing electromechanical equations are derived using the constitutive relations of piezoelectric materials and Maxwell’s equations. To model the host system with piezoelectric material as both a complete shell and a panel, and to represent the piezoelectric material as a patch, the two-dimensional Generalized Differential Quadrature (2D–GDQ) method is employed. Utilizing of the 2D–GDQ approach is essential for accurately capturing variations in both the axial and circumferential directions, which cannot be achieved with a one-dimensional formulation. The natural frequencies of the shell and panel with piezoelectric patches are computed under various electrical and mechanical boundary conditions. The accuracy of proposed model and numerical approach is validated through comparisons with published results and finite element simulations conducted in Abaqus software. Furthermore, the study examines the influence of key geometric parameters—including the radius and length of the host panel, as well as the thicknesses of the piezoelectric patch and the host shell—on the natural frequencies.

Keywords

free vibration piezoelectric patch smart structures two-dimensional generalized differential quadratic first-order shear deformation theory cylindrical shell/panel

1. Introduction

Shells of revolution, such as cylindrical, conical, and spherical shells, are extensively used in various engineering applications, including aerospace, mechanical, and civil industries. Due to their high strength-to-weight ratio, these structures are highly regarded for their ability to sustain significant loads while maintaining a low weight. Among these, cylindrical shells are particularly favored in industrial applications owing to their unique geometric characteristics. In recent decades, the integration of piezoelectric materials has gained considerable traction, especially in fields like energy harvesting and active vibration control. The combination of piezoelectric materials with structural elements such as shells offers promising opportunities to develop multifunctional systems capable of both structural performance and energy conversion.

The equations of motion of cylindrical shells have been derived based on various shell theories, each characterized by distinct assumptions. The foundational contributions to shell theory trace back to the works of Love (1888) and Flügge (2013). Over time, numerous variants of shell theories have been proposed, some simplifying Flügge’s original equations by neglecting certain terms. Extensive research has been conducted over the years to refine these models and enhance their applicability (Sanders, 1960; Holmes, 1989; Carrera, 1991, 1999; Bhimaraddi, 1984). Furthermore, various solution methods have been developed to solve the governing equations of shells, ranging from energy-based approaches such as the Galerkin and Rayleigh-Ritz methods (Li et al., 2022; Pang et al., 2023; Zhou et al., 2024; Weiwei et al., 2024; Xu et al., 2025), to exact analytical solutions explored in studies of Li et al. (2018); Forsberg (1969); Khalifa (2011); Goldman (1974).

In recent years, the increasing global demand for sustainable and renewable energy sources has propelled research efforts in energy harvesting technologies. Among the various techniques, piezoelectric energy harvesting has emerged as a highly effective solution due to its direct conversion of ambient mechanical vibrations into usable electrical energy. This feature is particularly valuable in environments where traditional power sources are impractical or where long-term, maintenance-free operation is desired. Structural components such as beams, plates, and shells frequently encounter dynamic loads and environmental vibrations in diverse engineering applications, including aerospace structures, bridges, and mechanical systems. By embedding piezoelectric materials into these structural elements, it becomes feasible to capture mechanical vibrations and transform them into electrical power, thereby transforming passive structures into active energy-harvesting devices. Numerous studies have investigated the implementation of piezoelectric materials on beams, plates, and shells to optimize energy extraction efficiency. These studies typically focus on developing coupled electromechanical models to accurately describe the interaction between mechanical deformation and piezoelectric response. Optimization of factors such as material properties, boundary conditions, and excitation frequencies plays a crucial role in enhancing energy harvesting performance.

In 1997, for the first time, Umeda et al. (1997) explored the fundamental concepts of a generator that converted mechanical energy into electrical energy using a piezoelectric vibrator and a steel ball. They examined the effects of various piezoelectric vibrator characteristics and proposed an equivalent electrical circuit to simulate the power generation mechanism. Their findings revealed that maximum frequencies corresponded to the resonant vibration modes, and the calculated output voltage waveforms closely matched experimental measurements. The study concluded that increasing the electromechanical coupling coefficient and reducing the dielectric constant led to improved efficiency. Smits and Choi (1991) derived the fundamental equations for a unimorph piezoelectric actuator mounted on a cantilever beam. They analyzed the electromechanical behavior of the beam under various boundary conditions, considering both applied moment at the free end and distributed force along the beam. In a subsequent study Smits et al. (1991) extended their work to a beam with two bonded piezoelectric layers under similar boundary and loading conditions. Qing et al. (2006) introduced a semi-analytical method for static and dynamic analysis of a plate with multiple piezoelectric patches. Their model incorporated shear stresses and inertial moments within the plate equations. In the realm of advanced composite plates, research has been conducted to investigate the free vibration and vibro-acoustic responses of two types of piezoelectric functionally graded carbon nanotube-reinforced composite (FG-CNTRC) plates. The dynamic equations were derived using third-order shear deformation theory and Hamilton’s principle. Employing Navier’s method, the vibration frequencies and structural-acoustic responses of simply supported plates under concentrated harmonic excitation were obtained, with theoretical results validated against finite element simulations. The study compared the two configurations and discussed the influence of parameter variations on the dynamic and acoustic behavior (Li et al., 2025). Furthermore, the free vibration and dynamic response of bidirectional functionally graded (2D-FGM) plates integrated with piezoelectric layers and subjected to a moving mass have been analyzed. A finite element model based on a four-variable refined plate theory (HSDT4) was developed, and the Newmark method was employed to compute transient responses. The effects of parameters such as moving mass velocity and weight, along with electrical and mechanical boundary conditions, were comprehensively examined (Tran et al., 2023). In a complementary methodological approach, a Spectral-Tchebychev (S-T) dynamic model was proposed to analyze the free and random vibration characteristics of plates with bonded piezoelectric patches under random excitation. Utilizing a multi-partition strategy and first-order shear deformation theory (FSDT), the governing equations were derived for sub-plates, with boundary and coupling conditions handled via an artificial spring technique. Random excitation was introduced using the pseudo-excitation method (PEM). The model’s accuracy was verified against finite element results and existing data, and the influence of key parameters—such as patch position, size, and type—on the vibrational response was systematically investigated Gao et al. (2023).

Sheng and Wang (2009), employing the first-order shear deformation theory (FSDT) and Maxwell’s equations, formulated the electromechanical equations for a cylindrical shell subjected to external moving excitations. They utilized a hybrid approach combining modal analysis and integral methods to accurately capture the relationship between displacement and electric potential. Alibeigloo and Kani (2010) applied three-dimensional elasticity theory and combined the differential quadrature method (DQM) with state-space techniques to study the free vibrations of a multilayered cylindrical shell with a piezoelectric layer. They developed a general solution using Fourier series expansions for simply supported boundary conditions, validated their results against established references, and examined the influence of various parameters such as boundary conditions, radius-to-thickness ratio, length-to-radius ratio, and piezoelectric layer thickness on the natural frequencies. Another work provided a rigorous three-dimensional elasticity solution for the static and free vibration behavior of functionally graded carbon nanotube-reinforced cylindrical panels embedded in piezoelectric layers, focusing on the influence of the material gradient index and CNT distribution (Asadinia et al., 2025). Further, Bodaghi and Shakeri (2012) investigated the free and excited vibration of a simply supported functionally graded material (FGM) cylindrical panel covered with a piezoelectric layer. They subjected the panel to time-dependent impact forces and derived the dynamic response using FSDT and Maxwell’s equations. This study also considered open-circuit and short-circuit electrical boundary conditions and examined their influence on system behavior, employing Fourier series and Laplace transform techniques for the solution. Sayyaadi and Farsangi (2014) derived the dynamic equations of a thick doubly curved shell with a piezoelectric layer using third-order shear deformation theory. By applying simply supported boundary conditions and eliminating spatial derivatives, they achieved an analytical solution for the panel’s dynamic response, evaluating the effect of the piezoelectric layer. Building upon these foundations, recent investigations have further explored piezoelectric-integrated structures with advanced geometries and materials. For doubly curved shells, one study analyzed the free vibration of composite micro-panels with piezoelectric layers using a modified Sander’s shell theory and the Chebyshev–Ritz method, examining the effects of graphene platelet distribution and applied voltage (Wang et al., 2025). Arefi et al. (2016) conducted vibrational analysis of a cylindrical FGM shell integrated with two piezoelectric layers acting as sensor and actuator. Utilizing FSDT, von Kármán assumptions, and Hamilton’s principle, they formulated the nonlinear governing equations, discretized them via Navier’s method, and performed numerical analysis to explore the effects of geometrical and material parameters. The Generalized Differential Quadrature (GDQ) method is an efficient numerical technique for solving partial differential equations, first introduced by Shu Shu (2012) in 1997. Its high accuracy and ability to handle complex boundary conditions and geometries have made it a popular choice for analyzing advanced structures, including composite plates and shells and functionally graded beams with defects or reinforcements. Recent studies have effectively demonstrated its capabilities in this domain. For instance, the GDQ method has been applied to analyze composite plates with through-width delaminations using semi-layerwise modeling, successfully capturing stress fields and energy release rates at the delamination front (Szekrenyes, 2021, 2022). The method has also been extended to study the free vibration of bi-directional functionally graded (2D-FGM) beams, where it accurately determines natural frequencies and coupled mode shapes by transforming the governing differential equations into an algebraic eigenvalue problem (Huang et al., 2024). Furthermore, GDQ has been employed to predict the buckling loads and free vibration frequencies of cenosphere-reinforced syntactic foam beams under various boundary conditions and loading types, showing excellent agreement with both theoretical and experimental results (Duryodhana et al., 2023). These studies collectively validate GDQ as a robust and versatile tool for the static, buckling, and dynamic analysis of complex composite and functionally graded structures. Regarding the application of this method to the vibration of piezoelectric shells, several studies have also been conducted. Azaripour and Baghani (2019) applied the 1D-GDQ method to analyze the free vibrations of a circular FGM plate with a piezoelectric layer. Their formulation integrated FSDT and Maxwell’s equations, considering various mechanical and electrical boundary conditions. They systematically studied the influence of boundary conditions, plate and piezoelectric layer thicknesses, and plate radius on the system’s natural frequencies. Khorshidi et al. (2020) examined the free vibrations of a doubly curved panel with two piezoelectric layers on its upper and lower surfaces. They employed FSDT, incorporated FGMs and carbon nanotubes, and subjected the panel to thermal environments. Using the 1D-GDQ method, they solved the governing differential equations and demonstrated that natural frequencies decrease under elevated temperatures and short-circuit conditions. In the realm of aeroelasticity, research has addressed the supersonic flutter of cylindrical sandwich shells with reinforced honeycomb cores and piezoelectric face sheets, employing a semi-analytical approach to reveal that core stiffening does not always enhance flutter stability (Givi et al., 2025). Furthermore, studies have demonstrated the active control of structural dynamics through piezoelectric elements. One approach numerically and experimentally showed how piezoelectric-induced prestress can tune the natural frequencies of plates and shallow shells (Kamenskikh et al., 2024), while another investigated using shunted piezoelectric circuits attached to cylindrical shells for passive flutter suppression and enhanced vibration damping (Bochkarev and Lekomtsev, 2024).

While the majority of the aforementioned studies have focused on beams and plates, the complex nature of shell structures and the associated challenges in solving their electromechanical equations have resulted in fewer investigations on piezoelectric cylindrical shells. Moreover, most studies have modeled complete shells, whereas in practical industrial applications, cylindrical panels are frequently used. Additionally, piezoelectric materials are often modeled as continuous layers on the host structures, which is not economically viable. Instead, using piezoelectric patches is more cost-effective and practical. According to these research gaps, the present study focuses on the free vibration analysis of cylindrical shells and panels with bonded piezoelectric patches under both open-circuit and short-circuit electrical conditions. To achieve higher accuracy and account for moderately thick shells and panels, FSDT is adopted, and Hamilton’s principle is employed to derive the governing equations. By coupling the constitutive relations of piezoelectric materials and Maxwell’s equations with mechanical governing equations of the cylindrical shell, the electromechanical equations of cylindrical shells with piezoelectric patches are formulated. Subsequently, the 2–Dimensional Generalized Differential Quadrature method (2D–GDQ) is utilized to discrete equations of motion. Then, the effective parameters on natural frequencies of propose structures will be studied.

2. Mathematical modeling

2.1. Governing equations of motion

Cylindrical panels embedded with piezoelectric patch are defined by vital parameters in Figure 1. Cylindrical and piezoelectric patch are defined by cylindrical coordinate system (x, θ, z). Thickness and length of cylindrical panel and piezoelectric patch are h, L and h_p, L_p, respectively. In this study, a cylindrical panel radius of R_c is considered for sector angle θ ranging from zero to θ_c and patch angle spanning from θ₁ to θ₂. Considering θ_c = 2π, the cylindrical panel will be changed into the cylindrical shell.

Figure 1.

Schematic of a cylindrical panel embedded with piezoelectric patch respective coordinate system.

As the piezoelectric patch is always connected to the panel, the relations of the strain and displacement function of the panel also hold for the patch. Based on FSDT, the displacement function field of panel and patch are defined as follows (Reddy, 2006),

\begin{matrix} u_{1} (x, θ, z, t) & = & u (x, θ, t) + z ϕ_{1} (x, θ, t) \\ u_{2} (x, θ, z, t) & = & v (x, θ, t) + z ϕ_{2} (x, θ, t) \\ u_{3} (x, θ, z, t) & = & w (x, θ, t) \end{matrix}

(1)

The strain relations for an arbitrary point on panel and patch are also formulated as,

\{\begin{matrix} ε_{x x} \\ ε_{θ θ} \\ γ_{x θ} \\ γ_{x z} \\ γ_{θ z} \end{matrix}\} = \{\begin{matrix} ε_{1}^{0} + z κ_{1} \\ ε_{2}^{0} + z κ_{2} \\ \{w_{1}^{0} + w_{2}^{0}\} + z \{w_{1}^{1} + w_{2}^{1}\} \\ ε_{5}^{0} \\ ε_{6}^{0} \end{matrix}\} = \{\begin{matrix} ε_{1}^{0} + z κ_{1} \\ ε_{2}^{0} + z κ_{2} \\ ε_{12}^{0} + z κ_{12} \\ ε_{5}^{0} \\ ε_{6}^{0} \end{matrix}\}

(2)

where the 0 superscript represents the strain of mid-plane surface and κ is the changing of curvature, which may be defined as follows,

\begin{matrix} \{\begin{matrix} ε_{1}^{0} \\ ε_{2}^{0} \\ ε_{12}^{0} \\ ε_{5}^{0} \\ ε_{6}^{0} \end{matrix}\} & = & \{\begin{matrix} u_{, x} \\ \frac{v_{, θ}}{R} + \frac{w}{R} \\ v_{, x} + \frac{u_{, θ}}{R} \\ w_{, x} + ϕ_{1} \\ \frac{w_{, θ}}{R} - \frac{v}{R} + ϕ_{2} \end{matrix}\} \end{matrix}

(3)

\begin{matrix} \{\begin{matrix} κ_{1} \\ κ_{2} \\ κ_{12} \end{matrix}\} & = & \{\begin{matrix} ϕ_{1, x} \\ \frac{ϕ_{2, θ}}{R} \\ ϕ_{2, x} + \frac{ϕ_{1, θ}}{R} \end{matrix}\} \end{matrix}

(4)

In this study, the cylindrical panel is considered as an isotropic panel. Therefore, the constitutive relations between stress and strain can be expressed as follows,

\{\begin{matrix} σ_{x x} \\ σ_{θ θ} \\ σ_{x θ} \\ σ_{x z} \\ σ_{θ z} \end{matrix}\} = [\begin{matrix} C_{11} & C_{12} & 0 & 0 & 0 \\ C_{21} & C_{22} & 0 & 0 & 0 \\ 0 & 0 & C_{44} & 0 & 0 \\ 0 & 0 & 0 & C_{55} & 0 \\ 0 & 0 & 0 & 0 & C_{66} \end{matrix}] \{\begin{matrix} ε_{x x} \\ ε_{θ θ} \\ γ_{x θ} \\ γ_{x z} \\ γ_{θ z} \end{matrix}\}

(5)

where the components of stiffness matrix C_ij will be defined as,

\begin{matrix} C_{11} = C_{22} = \frac{E}{1 - υ^{2}}; C_{12} = C_{21} = \frac{E υ}{1 - υ^{2}}; \\ C_{44} = C_{55} = C_{66} = \frac{E υ}{2 (1 + υ)} \end{matrix}

(6)

where E and ν represent the elastic modulus and Poisson’s ratio of the panel, respectively. Moreover, the constitutive equation of piezoelectric patch is formulated as Alibeigloo and Kani (2010),

\begin{matrix} \{σ\} & = & [Q^{P}] \{ε\} - {[e^{p}]}^{T} \{E_{i}\} \end{matrix}

(7)

\begin{matrix} \{D\} & = & [e^{p}] \{ε\} + [Ξ^{p}] \{E_{i}\} \end{matrix}

(8)

where Ξ^p, e^p, and Q^p represent piezoelectric constant, dielectric constant, and the elastic constant, respectively. The ɛ and E_i denote the strain and electric field of the piezoelectric patch. Consequently, the stress and electric displacement in a cylindrical piezoelectric patch can be defined as follows

\begin{matrix} \{\begin{matrix} σ_{x x}^{P} \\ σ_{θ θ}^{P} \\ σ_{x θ}^{P} \\ σ_{x z}^{P} \\ σ_{θ z}^{P} \end{matrix}\} & = & [\begin{matrix} Q_{11}^{P} & Q_{12}^{P} & 0 & 0 & 0 \\ Q_{12}^{P} & Q_{22}^{P} & 0 & 0 & 0 \\ 0 & 0 & Q_{44}^{P} & 0 & 0 \\ 0 & 0 & 0 & Q_{55}^{P} & 0 \\ 0 & 0 & 0 & 0 & Q_{66}^{P} \end{matrix}] \{\begin{matrix} ε_{x x} \\ ε_{θ θ} \\ γ_{x θ} \\ γ_{x z} \\ γ_{θ z} \end{matrix}\} \\ - & [\begin{matrix} 0 & 0 & e_{31} \\ 0 & 0 & e_{32} \\ \begin{matrix} 0 \\ e_{15} \\ 0 \end{matrix} & \begin{matrix} 0 \\ 0 \\ e_{24} \end{matrix} & \begin{matrix} 0 \\ 0 \\ 0 \end{matrix} \end{matrix}] \{\begin{matrix} E_{x} \\ E_{θ} \\ E_{z} \end{matrix}\} \end{matrix}

(9)

\begin{align} \{\begin{matrix} D_{x} \\ D_{θ} \\ D_{z} \end{matrix}\} & = [\begin{matrix} 0 & 0 & 0 & e_{15} & 0 \\ 0 & 0 & 0 & 0 & e_{24} \\ e_{31} & e_{32} & 0 & 0 & 0 \end{matrix}] \{\begin{matrix} ε_{x x} \\ ε_{θ θ} \\ γ_{x θ} \\ γ_{x z} \\ γ_{θ z} \end{matrix}\} \\ - [\begin{matrix} Ξ_{11} & 0 & 0 \\ 0 & Ξ_{22} & 0 \\ 0 & 0 & Ξ_{33} \end{matrix}] \{\begin{matrix} E_{x} \\ E_{θ} \\ E_{z} \end{matrix}\} \end{align}

(10)

where the components of equation (9) and equation (10) are obtained by considering polarization in the thickness direction as follows,

\begin{matrix} Q_{11}^{P} = C_{11}^{P} - \frac{{C_{13}^{P}}^{2}}{C_{33}^{P}}; Q_{12}^{P} = C_{12}^{P} - \frac{C_{13}^{P} C_{23}^{P}}{C_{33}^{P}}; \end{matrix}

\begin{matrix} Q_{44}^{P} = \frac{Q_{11}^{P} - Q_{12}^{P}}{2}; Q_{22}^{P} = C_{22}^{P} - \frac{{C_{23}^{P}}^{2}}{C_{33}^{P}}; \end{matrix}

\begin{matrix} Q_{55}^{P} = C_{55}^{P}; Q_{66}^{P} = C_{66}^{P} \end{matrix}

(11)

\begin{matrix} e_{31} = e_{31}^{P} - \frac{C_{13}^{P}}{C_{33}^{P}} e_{33}^{P}; e_{24} = e_{24}^{P}; e_{15} = e_{15}^{P}; \end{matrix}

\begin{matrix} e_{32} = e_{32}^{P} - \frac{C_{23}^{P}}{C_{33}^{P}} e_{33}^{P} \end{matrix}

(12)

\begin{matrix} Ξ_{11} = Ξ_{11}^{P}; Ξ_{22} = Ξ_{22}^{P}; Ξ_{33} = Ξ_{33}^{P} + \frac{{e_{33}^{P}}^{2}}{C_{33}^{P}} \end{matrix}

(13)

According to the piezoelectric constitutive equations, the electric field generated in the piezoelectric patch directly affects on the stresses. Furthermore, by knowing the electric potential generated within the piezoelectric material, the field in the piezoelectric patch can be calculated utilizing equations equations (14) and (15). In these equations, the function Ψ(x, θ, z) represents the electric potential, which can be determined based on the piezoelectric boundary conditions.

E = - ▿ \nabla Ψ (x, θ, z)

(14)

where

\begin{matrix} E_{x} = - \frac{\partial Ψ (x, θ, z)}{\partial x}; \\ E_{θ} = - \frac{1}{R (1 + z / R)} \frac{\partial Ψ (x, θ, z)}{\partial θ}; \\ E_{z} = - \frac{\partial Ψ (x, θ, z)}{\partial z} \end{matrix}

(15)

The piezoelectric patch is considered as a sensor, where the electric potential is generated due to elastic deformations caused by vibrations. The shell has no electric potential, also the bottom layer of the piezoelectric patch is connected to the panel/shell. Therefore, the bottom surface of piezoelectric patch must follow the same condition based on the compatibility requirement. Thus, the bottom layer of the piezoelectric patch must have zero electric potential, which known as grounding.

Furthermore, the top layer of piezoelectric condition are divided into two phases, namely, zero and non-zero potential. If both of top and bottom layer’s conditions are assumed ground condition, it will be represented as short-circuit condition. If the top layer will be fully isolated, the electrical displacement and voltage have zero and non-zero value, respectively. This condition will be represented open-circuit. The relations of short and open-circuit are defined as the following, for Short-circuit:

Ψ (x, θ, \frac{h}{2}, t) = 0 a n d Ψ (x, θ, \frac{h}{2} + h_{p}) = 0

(16)

for Open-circuit:

Ψ (x, θ, \frac{h}{2}, t) = 0 a n d D_{z} (x, θ, \frac{h}{2} + h_{p}) = 0

(17)

Assuming the electric potential function as follows Bodaghi and Shakeri (2012); Khorshidi et al. (2020),

\begin{matrix} Ψ (x, θ, z, t) & = & A (x, θ, t) + z B (x, θ, t) \\ + & (z_{p}^{2} - {(\frac{h_{p}}{2})}^{2}) ψ (x, θ, t); \end{matrix}

(18)

\begin{matrix} z_{p} & = & z - \frac{h + h_{p}}{2} \end{matrix}

(19)

The electric potential function of the patch can be calculated by applying the boundary conditions from equation (17) into equation (18) as follows, for Short-circuit:

Ψ (x, θ, z, t) = (z_{p}^{2} - {(\frac{h_{p}}{2})}^{2}) ψ (x, θ, t)

(20)

for Open-circuit:

\begin{matrix} Ψ (x, θ, z, t) & = & (z_{p}^{2} - {(\frac{h_{p}}{2})}^{2} + h z_{p} - 2 z z_{p}) ψ (x, θ, t) \\ + & (z - \frac{h}{2}) \hat{C} \end{matrix}

(21)

where

\hat{C} = \frac{e_{33}}{Ξ_{33}} (ε_{1}^{0} + (\frac{h}{2} + h_{p}) κ_{1}) + \frac{e_{33}}{Ξ_{33}} (ε_{2}^{0} + (\frac{h}{2} + h_{p}) κ_{2})

(22)

The resultant forces and moments in the cylindrical panel and shell are defined as follows,

\{\begin{matrix} N_{x x} \\ N_{θ θ} \\ N_{x θ} \\ Q_{x z} \\ Q_{θ x} \end{matrix}\} = \int_{- \frac{h}{2}}^{\frac{h}{2}} \{\begin{matrix} σ_{x x} \\ σ_{θ θ} \\ τ_{x θ} \\ K_{s} τ_{x z} \\ K_{s} τ_{θ x} \end{matrix}\} d z

(23)

\{\begin{matrix} M_{x x} \\ M_{θ θ} \\ M_{x θ} \end{matrix}\} = \int_{- \frac{h}{2}}^{\frac{h}{2}} \{\begin{matrix} σ_{x x} \\ σ_{θ θ} \\ τ_{x θ} \end{matrix}\} z d z

(24)

In the above equation, κ_s is the shear correction factor, which is taken as 5/6. By substituting the stress values into the resultant forces and moments equations, they will be defined as follows,

\begin{matrix} \begin{matrix} \{\begin{matrix} N_{x x} \\ N_{θ θ} \\ N_{x θ} \\ M_{x x} \\ M_{θ θ} \\ M_{x θ} \end{matrix}\} = [\begin{matrix} A_{11} & A_{12} & 0 & B_{11} & B_{12} & 0 \\ A_{12} & A_{22} & 0 & B_{12} & B_{22} & 0 \\ 0 & 0 & A_{44} & 0 & 0 & B_{44} \\ B_{11} & B_{12} & 0 & D_{11} & D_{12} & 0 \\ B_{12} & B_{22} & 0 & D_{12} & D_{22} & 0 \\ 0 & 0 & B_{44} & 0 & 0 & D_{44} \end{matrix}] \{\begin{matrix} ε_{1}^{0} \\ ε_{2}^{0} \\ ε_{12}^{0} \\ κ_{1} \\ κ_{2} \\ κ_{12} \end{matrix}\} \\ \{\begin{matrix} Q_{x z} \\ Q_{θ z} \end{matrix}\} = K_{s} [\begin{matrix} A_{55} & 0 \\ 0 & A_{66} \end{matrix}] \{\begin{matrix} ε_{5}^{0} \\ ε_{6}^{0} \end{matrix}\} \end{matrix} \end{matrix}

(25)

By substituting Equation Equation (9) into Equations (23) and (24), the resultant forces and moments in the piezoelectric patch will be as follows,

\begin{matrix} \begin{matrix} \{\begin{matrix} N_{x x}^{P} \\ N_{θ θ}^{P} \\ N_{x θ}^{P} \\ M_{x x}^{P} \\ M_{θ θ}^{P} \\ M_{x θ}^{P} \end{matrix}\} = [\begin{matrix} A_{11}^{P} & A_{12}^{P} & 0 & B_{11}^{P} & B_{12}^{P} & 0 \\ A_{12}^{P} & A_{22}^{P} & 0 & B_{12}^{P} & B_{22}^{P} & 0 \\ 0 & 0 & A_{44}^{P} & 0 & 0 & B_{44}^{P} \\ B_{11}^{P} & B_{12}^{P} & 0 & D_{11}^{P} & D_{12}^{P} & 0 \\ B_{12}^{P} & B_{22}^{P} & 0 & D_{12}^{P} & D_{22}^{P} & 0 \\ 0 & 0 & B_{44}^{P} & 0 & 0 & D_{44}^{P} \end{matrix}] \{\begin{matrix} ε_{1}^{0} \\ ε_{2}^{0} \\ ε_{12}^{0} \\ κ_{1}^{0} \\ κ_{2}^{0} \\ κ_{12}^{0} \end{matrix}\} - \{\begin{matrix} H_{x x}^{P} \\ H_{θ θ}^{P} \\ 0 \\ R_{x x}^{P} \\ R_{θ θ}^{P} \\ 0 \end{matrix}\} \\ \{\begin{matrix} Q_{x z}^{P} \\ Q_{θ z}^{P} \end{matrix}\} = K_{s} [\begin{matrix} A_{55}^{P} & 0 \\ 0 & A_{66}^{P} \end{matrix}] \{\begin{matrix} ε_{5}^{0} \\ ε_{6}^{0} \end{matrix}\} - \{\begin{matrix} H_{x z}^{P} \\ H_{θ z}^{P} \end{matrix}\} \end{matrix} \end{matrix}

(26)

where the

H_{i j}^{P}

and

R_{i j}^{P}

of the equation (26) will be defined as follows,

\begin{matrix} (H_{x x}^{P}, R_{x x}^{P}) = \int_{\frac{h}{2}}^{\frac{h}{2} + h_{p}} (1, z) e_{31} E_{z} d z \end{matrix}

(27)

\begin{matrix} (H_{θ θ}^{P}, R_{θ θ}^{P}) = \int_{\frac{h}{2}}^{\frac{h}{2} + h_{p}} (1, z) e_{32} E_{z} d z \end{matrix}

(28)

\begin{matrix} H_{x z}^{P} = \int_{\frac{h}{2}}^{\frac{h}{2} + h_{p}} (1, z) e_{15} E_{x} d z \end{matrix}

(29)

\begin{matrix} H_{θ z}^{P} = \int_{\frac{h}{2}}^{\frac{h}{2} + h_{p}} (1, z) e_{24} E_{θ} d z \end{matrix}

(30)

The matrix components of Equations (25) and (26) are obtained by integrating the stiffness matrices over the thickness direction as follows,

\begin{matrix} (A_{i j}, B_{i j}, D_{i j}) = \int_{- \frac{h}{2}}^{\frac{h}{2} + h_{p}} (1, z, z^{2}) C_{i j} d z \end{matrix}

(31)

\begin{matrix} (A_{i j}^{P}, B_{i j}^{P}, D_{i j}^{P}) = \int_{- \frac{h}{2}}^{\frac{h}{2} + h_{p}} (1, z, z^{2}) Q_{i j}^{P} d z \end{matrix}

(32)

Moreover, the governing equations are derived using Hamilton’s principle, which is formulated as follows,

δ H = \int_{t_{1}}^{t_{2}} δ L d t = 0; L = T - U + W

(33)

where T, U, and W denote the kinetic energy, potential energy, and the work of external forces, respectively. Therefore, the calculus of variations and Hamilton’s principle are utilized to derive the equations of motion of the cylindrical shell embedded with Piezoelectric patch as follows,

\begin{matrix} \frac{\partial}{\partial x} ({\bar{N}}_{x x}) + \frac{1}{R} \frac{\partial}{\partial θ} ({\bar{N}}_{x θ}) + q_{x} = {\bar{I}}_{0} \ddot{u} + {\bar{I}}_{1} {\ddot{ϕ}}_{1} \\ \frac{\partial}{\partial x} ({\bar{N}}_{x θ}) + \frac{1}{R} \frac{\partial}{\partial θ} ({\bar{N}}_{θ θ}) + \frac{1}{R} ({\bar{Q}}_{θ z}) + q_{θ} = {\bar{I}}_{0} \ddot{v} + {\bar{I}}_{1} {\ddot{ϕ}}_{2} \\ \frac{\partial}{\partial x} ({\bar{Q}}_{x z}) - \frac{1}{R} ({\bar{N}}_{θ θ}) + \frac{1}{R} \frac{\partial}{\partial θ} ({\bar{Q}}_{θ z}) + q_{z} = {\bar{I}}_{0} \ddot{w} \\ \frac{\partial}{\partial x} ({\bar{M}}_{x x}) + \frac{1}{R} \frac{\partial}{\partial θ} ({\bar{M}}_{x θ}) - {\bar{Q}}_{x z} = {\bar{I}}_{1} \ddot{u} + {\bar{I}}_{2} {\ddot{ϕ}}_{1} \\ \frac{\partial}{\partial x} ({\bar{M}}_{x θ}) + \frac{1}{R} \frac{\partial}{\partial θ} ({\bar{M}}_{θ θ}) - {\bar{Q}}_{θ z} = {\bar{I}}_{1} \ddot{v} + {\bar{I}}_{2} {\ddot{ϕ}}_{2} \end{matrix}

(34)

where q_x, q_θ, and q_z represent the external forces applied to the shell in the axial, circumferential, and thickness directions, respectively. The resultant force and moments are also defined as follows,

\begin{matrix} {\bar{N}}_{x x} = N_{x x} + N_{x x}^{P} (P (x, θ)) \\ {\bar{N}}_{θ θ} = N_{θ θ} + N_{θ θ}^{P} (P (x, θ)) \\ {\bar{M}}_{x x} = M_{x x} + M_{x x}^{P} (P (x, θ)) \\ {\bar{M}}_{θ θ} = M_{θ θ} + M_{θ θ}^{P} (P (x, θ)) \\ {\bar{N}}_{x θ} = N_{x θ} + N_{x θ}^{P} (P (x, θ)) \\ {\bar{M}}_{x θ} = M_{x θ} + M_{x θ}^{P} (P (x, θ)) \\ {\bar{Q}}_{x z} = Q_{x z} + Q_{x z}^{P} (P (x, θ)) \\ {\bar{Q}}_{θ z} = Q_{θ z} + Q_{θ z}^{P} (P (x, θ)) \end{matrix}

(35)

The ${\bar{I}}_{0,1,2}$ known as the inertia terms, which are defined as follows,

\begin{matrix} {\bar{I}}_{0,1,2} = I_{0,1,2} + I_{0,1,2}^{P} (P (x, θ)) \\ I_{0,1,2} = \int_{- \frac{h}{2}}^{\frac{h}{2}} ρ (1, z, z^{2}) d z \\ I_{0,1,2}^{P} = \int_{\frac{h}{2}}^{\frac{h}{2} + h_{P}} ρ^{P} (1, z, z^{2}) d z \end{matrix}

(36)

The superscript p in Equations (35) and (36) represents the piezoelectric portion. Additionally, the function P (x, θ) is defined as,

P (x, θ) = (H (x - x_{1}) - H (x - x_{2})) (H (θ - θ_{1}) - H (θ - θ_{2}))

(37)

where H is represent the Heaviside function. To establish the relationship between the mechanical displacement and electric potential (voltage), Maxwell’s equations in cylindrical coordinates are employed, leading to the derivation of the sixth equation, commonly referred to as the electromechanical coupling equation as follows,

P (x, θ) \int_{h / 2}^{h / 2 + h_{P}} (\partial / \partial x D_{x} + 1 / R \partial / \partial θ D_{θ} + \partial / \partial z D_{z} + D_{z} / R) d z = 0

(38)

Therefore, in regions without piezoelectric material, the system consists of five coupled mechanical equations, while in regions with piezoelectric patches, it includes six coupled electromechanical equations, which must be solved simultaneously. It is important to note that the derivative of the Heaviside function with respect to the x and θ directions appears in the above equations, which corresponds to the Dirac delta function δ, defined as follows,

\begin{matrix} \frac{\partial}{\partial x} P (x, θ) = (δ (x - x_{1}) - δ (x - x_{2})) (H (θ - θ_{1}) - H (θ - θ_{2})) \\ \frac{\partial}{\partial θ} P (x, θ) = (H (x - x_{1}) - H (x - x_{2})) (δ (θ - θ_{1}) - δ (θ - θ_{2})) \end{matrix}

(39)

If the cylindrical panel is defined in the θ direction from 0 to 2π radians and the edges at 0 and 2π are joined, the cylindrical panel becomes a cylindrical shell. To ensure these two edges align, a compatibility condition must be applied at this boundary. The compatibility conditions are also applied to the piezoelectric patch. Specifically, if the piezoelectric patch is considered as a ring with edges at 0 and 2π, in addition to the mechanical compatibility conditions, electrical compatibility conditions must also be enforced between these two edges. Section. 2.2 will be represented the equations of boundary and compatibility condition.

2.2. Boundary and compatibility condition

If the host structure with the piezoelectric patch is considered as a panel, it will have four edges, which has shown in Figure 2. Each of these four edges can be subjected to one of three mechanical boundary conditions: clamped (C), simply supported (S), or free (F). The mathematical expressions for these boundary conditions at each edge are defined as follows,

1. Clamp boundary condition

\begin{matrix} θ = θ_{0} o r θ = θ_{1}; 0 < x \leq L \\ u = v = w = ϕ_{1} = ϕ_{2} = 0 \\ x = 0 o r x = L; θ_{0} < θ \leq θ_{1} \\ u = v = w = ϕ_{1} = ϕ_{2} = 0 \end{matrix}

(40)

2. Simply support boundary condition

\begin{matrix} θ = θ_{0} o r θ = θ_{1}; 0 < x \leq L \\ u = v = w = ϕ_{1} = 0, M_{θ θ} = 0 \\ x = 0 o r x = L; θ_{0} < θ \leq θ_{1} \\ u = v = w = ϕ_{2} = 0, M_{x x} = 0 \end{matrix}

(41)

3. Free boundary condition

\begin{matrix} θ = θ_{0} o r θ = θ_{1}; 0 < x \leq L \\ N_{θ θ} = N_{x θ} = M_{θ θ} = M_{x θ} = Q_{θ z} = 0 \\ x = 0 o r x = L; θ_{0} < θ \leq θ_{1} \\ N_{x x} = N_{x θ} = M_{x x} = M_{x θ} = Q_{x z} = 0 \end{matrix}

(42)

Figure 2.

Boundary conditions of cylindrical panels.

For the piezoelectric patch attached to the panel, the four edges can be subjected to electrical boundary conditions of either grounding or insulation. The mathematical expressions for these conditions are expressed as follows,

1. Ground condition

\begin{matrix} x = x_{1} o r x_{2} \Rightarrow Ψ = 0 \\ θ = θ_{1} o r θ_{2} \Rightarrow Ψ = 0 \end{matrix}

(43)

2. Insulation boundary condition

\begin{matrix} x = x_{1} o r x_{2} \Rightarrow \int_{\frac{h}{2}}^{\frac{h}{2} + h_{P}} D_{x} d z = 0 \\ θ = θ_{1} o r θ_{2} \Rightarrow \int_{\frac{h}{2}}^{\frac{h}{2} + h_{P}} D_{θ} d z = 0 \end{matrix}

(44)

If the edges at θ = 0 and θ = 2π of the cylindrical panel are joined, transforming the panel into a complete cylindrical shell, compatibility conditions must be applied at these two edges instead of boundary conditions. These compatibility conditions are divided into two categories: displacement-based mechanical conditions and force-based mechanical conditions, which are defined as follows,

\begin{matrix} u (x, 0, t) = u (x, 2 π, t) \\ v (x, 0, t) = v (x, 2 π, t) \\ w (x, 0, t) = w (x, 2 π, t) \\ ϕ_{1} (x, 0, t) = ϕ_{1} (x, 2 π, t) \\ ϕ_{2} (x, 0, t) = ϕ_{2} (x, 2 π, t) \end{matrix}

(45)

\begin{matrix} N_{θ θ} (x, 0, t) = N_{θ θ} (x, 2 π, t) \\ N_{x θ} (x, 0, t) = N_{x θ} (x, 2 π, t) \\ M_{θ θ} (x, 0, t) = M_{θ θ} (x, 2 π, t) \\ M_{x θ} (x, 0, t) = M_{x θ} (x, 2 π, t) \\ Q_{θ z} (x, 0, t) = Q_{θ z} (x, 2 π, t) \end{matrix}

(46)

Additionally, if the patch is considered as a ring, after applying the mechanical compatibility conditions, the electrical compatibility conditions must also be imposed on the sixth equation of the system (Maxwell’s equation). Due to the alignment of the edges at θ = 0 and θ = 2π, the electric potential and electric displacements at this boundary must also be equal, just as the forces and displacements were equal at this interface.

Ψ (x, 0, z, t) = Ψ (x, 2 π, z, t)

(47)

\int_{\frac{h}{2}}^{\frac{h}{2} + h_{P}} D_{θ} (x, 0, z, t) d z = \int_{\frac{h}{2}}^{\frac{h}{2} + h_{P}} D_{θ} (x, 2 π, z, t) d z

(48)

3. Solution approach

The partial differential equations and the boundary conditions as well as compatibility conditions will be discretized using the 2D–GDQ method. Therefore, the PDEs are transformed into a linear algebraic system of equations. Shu Shu (1991) developed the GDQ method to improve the Differential Quadrature (DQ) technique for calculating the weighting coefficients. DQ method is a technique which transforms the partial differential system of equations into the linear algebraic system of equations. According to this method, the partial derivative of a smooth function with respect to a variable is approximated by a weighted sum of function values at all discrete points in that direction. The weighting coefficients depend on the grid points, which is written with arbitrary grid distribution, and the derivative order. Shu Shu (2012) improved the DQ method based on analysis of linear vector space and high-order polynomial approximation. It is well known that a smooth function Φ(x) in the closed interval [a, b] can be accurately approximated by a high-order polynomial. Thus, if Φ(x) is the solution of a PDE, then it can be approximated by a polynomial which can be written as Motazedian et al. (2024),

Φ (x) ≅ \sum_{j = 1}^{N} p_{j} (x) Φ (x_{j})

(49)

where N and x_j are the number of grid points in the whole domain and the coordinates of grid points in the closed interval [a, b], respectively. Φ(x_j) and p_j(x) are the function values at the grid points and the Lagrange interpolated polynomials, respectively. The functions can be defined as follows,

p_{j} (x) = \frac{L (x)}{(x - x_{j}) L^{(1)} (x_{j})}, j = 1,2, \dots, N

(50)

where

L (x) = \prod_{i = 1}^{N} (x - x_{i}), L^{(1)} (x_{j}) = \prod_{i = 1, i \neq j}^{N} (x_{j} - x_{i})

(51)

derivative of Φ(x) with respect to x at a certain point of the function domain is possible to obtain as,

\begin{matrix} Φ^{(1)} & ≅ & \sum_{j = 1}^{N} p_{j}^{(1)} (x_{i}) Φ (x_{j}) \\ = & \sum_{j = 1}^{N} ζ_{i j}^{(1)} Φ (x_{j}); i = 1,2, \dots, N \end{matrix}

(52)

where,

p_{j}^{(1)} (x_{i}) = ζ_{i j}^{(1)} = \frac{L^{(1)} (x_{i})}{(x_{i} - x_{j}) L^{(1)} (x_{j})}

(53)

where

ζ_{i j}^{(1)}

is the GDQ weighting coefficients of the first-order derivative. As can be observed from equation (52), the GDQ weighting coefficients are calculated using the coordinates of the grid points in the closed interval [a, b]. The different grid point distributions have difference effects on the convergence and accuracy of the GDQ method. One of the best distributions with fast convergence rate is the Chebyshev–Gauss–Lobatto distribution, which is assumed as follow,

x_{i} = \frac{1 - \cos (\frac{i - 1}{N - 1}) π}{2}; i = 1, 2, \dots, N

(54)

The GDQ method can also be employed for bi variate functions such as Φ(x, y). Approximating the function of Φ(x, y) with respect to x and y to a polynomial allows us to calculate the n-order derivatives with respect to x and the m-order derivatives with respect to y of the function Φ(x, y) as follows,

\begin{matrix} \frac{\partial^{n}}{\partial x^{n}} Φ (x_{i}, y_{j}) = \sum_{k = 1}^{N_{x}} A_{i k}^{(n)} Φ_{k j}; \begin{matrix} i = 1,2,3, \dots, N_{x} \\ j = 1,2,3, \dots, N_{y} \end{matrix} \end{matrix}

(55)

\begin{matrix} \frac{\partial^{m}}{\partial y^{m}} Φ (x_{i}, y_{j}) = \sum_{l = 1}^{N_{y}} B_{j l}^{(m)} Φ_{i l}; \begin{matrix} i = 1,2,3, \dots, N_{x} \\ j = 1,2,3, \dots, N_{y} \end{matrix} \end{matrix}

(56)

The partial derivative with respect to both x and y is defined by equation (57).

\frac{\partial^{m + n}}{\partial y^{m} \partial x^{n}} Φ (x_{i}, y_{j}) = \sum_{k = 1}^{N_{x}} A_{i k}^{(n)} \sum_{l = 1}^{N_{y}} B_{j l}^{(m)} Φ_{k l}; \begin{matrix} i = 1,2,3, \dots, N_{x} \\ j = 1,2,3, \dots, N_{y} \end{matrix}

(57)

where

A_{i j}^{(n)}

and

B_{i j}^{(n)}

are weighting coefficients of the GDQ method. To calculate these weighting coefficients, as mentioned earlier, it is necessary to discretize the x and y domains. Weighting coefficients can be defined as follows,

\begin{matrix} A_{i j}^{(1)} = \frac{L^{(1)} (x_{i})}{(x_{i} - x_{j}) L^{(1)} (x_{j})}; i, j = 1,2, \dots, N_{x} & i \neq j . \\ B_{i j}^{(1)} = \frac{G^{(1)} (y_{i})}{(y_{i} - y_{j}) G^{(1)} (y_{j})}; i, j = 1,2, \dots, N_{y} & i \neq j . \\ G^{(1)} (y_{i}) = \prod_{i = 1, i \neq j}^{N_{y}} (y_{j} - y_{i}) \end{matrix}

(58)

where x_i and y_i are defined as,

\begin{matrix} x_{i} = \frac{1 - \cos (\frac{i - 1}{N_{x} - 1}) π}{2}; i = 1,2, \dots, N_{x} . \\ y_{j} = \frac{1 - \cos (\frac{j - 1}{N_{y} - 1}) π}{2}; j = 1,2, \dots, N_{y} . \end{matrix}

(59)

N_x and N_y are the number of grid points in the x and y domains, respectively. When i equals j the weighting coefficients matrix can be obtained as follows,

\begin{matrix} \sum_{j = 1}^{N_{x}} A_{i j}^{(1)} = 0 \Rightarrow A_{i i}^{(1)} = - \sum_{j = 1, j \neq i}^{N_{x}} A_{i j}^{(1)} . \\ \sum_{j = 1}^{N_{y}} B_{i j}^{(1)} = 0 \Rightarrow B_{i i}^{(1)} = - \sum_{j = 1, j \neq i}^{N_{y}} B_{i j}^{(1)} . \end{matrix}

(60)

In addition, weighting coefficient matrix of the n^th-order and m^th-order derivative of the function $Φ (x, y)$ is defined as the following,

\begin{matrix} [A^{(n)}] = [A^{(1)}] [A^{(n - 1)}] \end{matrix}

(61)

\begin{matrix} [B^{(m)}] = [B^{(1)}] [B^{(m - 1)}] \end{matrix}

(62)

4. Numerical solution

After deriving the Navier equations of motion for the panels or shells, the 2D–GDQ method is employed to discretize the partial differential equations. During this process, internal and boundary points on the grid are distinguished, and the appropriate boundary and compatibility conditions are applied. Consequently, the simplified equations can be expressed as follows,

[M] \{\ddot{Δ}\} + [\bar{K_{1}}] \{Δ\} + [\bar{L_{1}}] \{ψ\} = 0

(63)

[\bar{K_{2}}] \{Δ\} + [\bar{L_{2}}] \{ψ\} = 0

(64)

where

{[Δ]}^{T} = \{u, v, w, ϕ_{1}, ϕ_{2}\}

(65)

The ψ is the electric potential (voltage) of the piezoelectric patch. The equation (65) represents Maxwell’s equation. Separating the boundary and internal points of the grid Equations (63) and (65) is reformulated as follows,

\begin{matrix} \begin{matrix} [\begin{matrix} 0 & 0 \\ 0 & M_{d d} \end{matrix}] {\begin{matrix} {\ddot{Δ}}_{b} \\ {\ddot{Δ}}_{d} \end{matrix}} + [\begin{matrix} K_{1bb} & K_{1bd} \\ K_{1db} & K_{1dd} \end{matrix}] \{\begin{matrix} Δ_{b} \\ Δ_{d} \end{matrix}\} \\ + [\begin{matrix} L_{1bb} & L_{1bd} \\ L_{1db} & L_{1dd} \end{matrix}] \{\begin{matrix} ψ_{b} \\ ψ_{d} \end{matrix}\} = \{\begin{matrix} 0 \\ 0 \end{matrix}\} \end{matrix} \end{matrix}

(66)

\begin{matrix} [\begin{matrix} K_{2bb} & K_{2bd} \\ K_{2db} & K_{2dd} \end{matrix}] \{\begin{matrix} Δ_{b} \\ Δ_{d} \end{matrix}\} + [\begin{matrix} L_{2bb} & L_{2bd} \\ L_{2db} & L_{2dd} \end{matrix}] \{\begin{matrix} ψ_{b} \\ ψ_{d} \end{matrix}\} = \{\begin{matrix} 0 \\ 0 \end{matrix}\} \end{matrix}

(67)

If the piezoelectric patch is located at the boundary, it interacts with the boundary nodes, causing the L_1bb, L_1bd, and L_1db terms in equation (66) and the L_2bb, L_2bd, L_2db, K_2bb, K_2bd, and K_2db terms in equation (67) to be nonzero. Otherwise, if the piezoelectric patch isn’t located at the boundary, values of terms will be zero. In the following, the equation (67) can be rewritten as follows,

\begin{matrix} \begin{matrix} [\begin{matrix} L_{2bb} & L_{2bd} \\ L_{2db} & L_{2dd} \end{matrix}] \{\begin{matrix} ψ_{b} \\ ψ_{d} \end{matrix}\} = - [\begin{matrix} K_{2bb} & K_{2bd} \\ K_{2db} & K_{2dd} \end{matrix}] \{\begin{matrix} Δ_{b} \\ Δ_{d} \end{matrix}\} \\ \Rightarrow \{\begin{matrix} ψ_{b} \\ ψ_{d} \end{matrix}\} = - {[\begin{matrix} L_{2bb} & L_{2bd} \\ L_{2db} & L_{2dd} \end{matrix}]}^{- 1} [\begin{matrix} K_{2bb} & K_{2bd} \\ K_{2db} & K_{2dd} \end{matrix}] \{\begin{matrix} Δ_{b} \\ Δ_{d} \end{matrix}\} \end{matrix} \end{matrix}

(68)

Therefore, equation (68) substitutes into equation (66) as follows,

\begin{matrix} \begin{matrix} [\begin{matrix} 0 & 0 \\ 0 & M_{d d} \end{matrix}] {\begin{matrix} {\ddot{Δ}}_{b} \\ {\ddot{Δ}}_{d} \end{matrix}} + [\begin{matrix} K_{1bb} & K_{1bd} \\ K_{1db} & K_{1dd} \end{matrix}] \{\begin{matrix} Δ_{b} \\ Δ_{d} \end{matrix}\} \\ - [\begin{matrix} L_{1bb} & L_{1bd} \\ L_{1db} & L_{1dd} \end{matrix}] {[\begin{matrix} L_{2bb} & L_{2bd} \\ L_{2db} & L_{2dd} \end{matrix}]}^{- 1} [\begin{matrix} K_{2bb} & K_{2bd} \\ K_{2db} & K_{2dd} \end{matrix}] \{\begin{matrix} Δ_{b} \\ Δ_{d} \end{matrix}\} = \{\begin{matrix} 0 \\ 0 \end{matrix}\} \end{matrix} \end{matrix}

(69)

The equation (69) may be defined as follows,

\begin{matrix} \begin{matrix} ([\begin{matrix} K_{1bb} & K_{1bd} \\ K_{1db} & K_{1dd} \end{matrix}] - [\begin{matrix} L_{1bb} & L_{1bd} \\ L_{1db} & L_{1dd} \end{matrix}] {[\begin{matrix} L_{2bb} & L_{2bd} \\ L_{2db} & L_{2dd} \end{matrix}]}^{- 1} [\begin{matrix} K_{2bb} & K_{2bd} \\ K_{2db} & K_{2dd} \end{matrix}]) \{\begin{matrix} Δ_{b} \\ Δ_{d} \end{matrix}\} \\ + [\begin{matrix} 0 & 0 \\ 0 & M_{d d} \end{matrix}] \{\begin{matrix} {\ddot{Δ}}_{b} \\ {\ddot{Δ}}_{d} \end{matrix}\} = \{\begin{matrix} 0 \\ 0 \end{matrix}\} \end{matrix} \end{matrix}

(70)

\Rightarrow [\begin{matrix} 0 & 0 \\ 0 & M_{d d} \end{matrix}] \{\begin{matrix} {\ddot{Δ}}_{b} \\ {\ddot{Δ}}_{d} \end{matrix}\} + [\begin{matrix} K_{b b} & K_{b d} \\ K_{d b} & K_{d d} \end{matrix}] \{\begin{matrix} Δ_{b} \\ Δ_{d} \end{matrix}\} = \{\begin{matrix} 0 \\ 0 \end{matrix}\}

(71)

The equation (71) may be simplified as follows,

[K_{b b}] Δ_{b} + [K_{b d}] Δ_{d} = \{0\} \Rightarrow Δ_{b} = - {[K_{b b}]}^{- 1} [K_{b d}] Δ_{d}

(72)

[M_{d d}] {\ddot{Δ}}_{d} + [K_{b d}] Δ_{b} + [K_{d d}] Δ_{d} = \{0\}

(73)

Substituting Δ_b from equation (72) into equation (73) as,

[M_{d d}] {\ddot{Δ}}_{d} - [K_{b d}] {[K_{b b}]}^{- 1} [K_{b d}] Δ_{d} + [K_{d d}] Δ_{d} = \{0\}

(74)

\Rightarrow [M_{d d}] {\ddot{Δ}}_{d} + [K_{total}] Δ_{d} = \{0\};

(75)

where

[K_{total}] = [K_{d d}] - [K_{b d}] {[K_{b b}]}^{- 1} [K_{b d}]

(76)

In this study, the free vibration of shells will be investigated. Therefore, the displacement functions may be defined as follows,

\begin{matrix} u (x, θ, t) = U (x, θ) \sin (ω_{n} t) \\ v (x, θ, t) = V (x, θ) \sin (ω_{n} t) \\ w (x, θ, t) = W (x, θ) \sin (ω_{n} t) \\ ϕ_{1} (x, θ, t) = Φ_{1} (x, θ) \sin (ω_{n} t) \\ ϕ_{2} (x, θ, t) = Φ_{2} (x, θ) \sin (ω_{n} t) \end{matrix}

(77)

Substituting equation (77) into equation (72), the natural frequencies will be obtained as follows,

(K_{total} - ω_{n}^{2} M_{d d}) Δ_{d} = 0 \Rightarrow |K_{total} - ω_{n}^{2} M_{d d}| = 0

(78)

Due to the Chebyshev–Gauss–Lobatto point distribution and the non-uniform spacing between nodes, MATLAB software is unable to accurately detect the location of the patch using the Dirac delta function, rendering the implementation of the 2D–GDQ method impractical. Consequently, an approximation of the Dirac delta function must be used Eftekhari (2016). Various functions have been proposed for this purpose, with the approximation adopted in this study given by Equation (79). In practice, wherever the Dirac delta function appears in the equations, it is replaced by the approximation in Equation (79).

\begin{matrix} δ (x - x_{1}) = \frac{1}{2 \sqrt{π n}} e^{- \frac{{(x - x_{1})}^{2}}{4 n}}, n \to 0 \\ δ (θ - θ_{1}) = \frac{1}{2 \sqrt{π n}} e^{- \frac{{(θ - θ_{1})}^{2}}{4 n}}, n \to 0 \end{matrix}

(79)

5. Results and discussions

In the following, the results of dynamic behavior of cylindrical panels/shells with embedded piezoelectric patch are provided in three phases, namely, Convergence, Verification, and Parametric Studies. To this aim, the mechanical properties of cylindrical panels/shells and piezoelectric patch are tabulated in Table 1 and Table 2, respectively.

Table 1.

The mechanical properties of panels and shells.

Property	unit	Value
Density (ρ)	kg/m³	7800
Young modulus (E)	GPa	210
Poisson’s ratio (ν)	—	0.3

Table 2.

The mechanical properties of piezoelectric of PZT–4.

Elastic properties (GPa)
C ₁₁	C ₁₂	C ₁₃	C ₂₂	C ₂₃	C ₃₃	C ₄₄	C ₅₅	C ₆₆
139	78	74	139	74	115	25.6	25.6	30.5

Electrical constant					Dielectric constant			Density
(c/m²)					(n/Fm²)			(kg/m³)
e ₁	e ₂	e ₃	e ₄	e ₅	Ξ₁₁	Ξ₂₂	Ξ₃₃	ρ
−5.2	−5.2	15.1	12.7	12.7	6.5	6.5	5.8	7500

In order to verify the results, the finite element software Abaqus is also utilized. The piezoelectric patch is modeled as a solid part, with the piezoelectric element type specified as C3D20RE, a 20-node element applied to the piezoelectric patch. The model was analyzed using a 3D stress formulation with C3D20 R elements. The term “3D stress element” refers to a three-dimensional continuum element. Also, The C3D20 R is a 20-node quadratic brick solid element with reduced integration. Tie constraint is employed to model the connection between the panel and the piezoelectric patch. Finally, the natural frequencies of the system are extracted using the frequency solver from the Linear Perturbation module. Mesh convergence studies were conducted for all dimensions and geometries analyzed in Abaqus, and the reported frequencies are independent of the mesh size and number.

5.1. Convergence

The geometry of cylindrical shell and panel are effective on convergence of grid points of 2D–GDQ method. Also, schematic piezoelectric such as being patch or layer are effective on convergence of grid points. Therefore, the convergence of first ten non–dimensional frequencies of cylindrical shell with embedded piezoelectric layer are tabulated in Table 3. The natural frequencies are provided for fully clamped boundary condition and have been non-dimensionalized by equation (80).

ω_{N D}^{i} = ω R \sqrt{ρ (1 - ν^{2}) / E}

(80)

where i indicates the number of natural frequency. Based on results, the frequencies are converged by considering N_x = N_y = N = 23.

Table 3.

Convergence study of the first ten non-dimensional frequencies (ω_ND) of a cylindrical shell embedded with a piezoelectric layer (N_x = N_y = N).

N	$ω_{N D}^{1}$	$ω_{N D}^{2}$	$ω_{N D}^{3}$	$ω_{N D}^{4}$	$ω_{N D}^{5}$	$ω_{N D}^{6}$	$ω_{N D}^{7}$	$ω_{N D}^{8}$	$ω_{N D}^{9}$	$ω_{N D}^{10}$
5	0.6988	0.7335	0.8937	1.008	1.4567	1.4730	1.5200	1.5291	1.5523	1.7075
6	0.6310	0.6310	0.8956	1.0056	1.1277	1.1277	1.1844	1.4582	1.4923	1.4923
7	0.6287	0.6292	0.8281	0.8957	1.0025	1.0057	1.1353	1.1355	1.1807	1.1840
8	0.5517	0.5789	0.6291	0.6304	0.8957	1.0025	1.0714	1.0725	1.121	1.1210
9	0.5346	0.5358	0.6291	0.6291	0.8957	1.0025	1.0588	1.0589	1.0641	1.1212
10	0.5358	0.5382	0.6290	0.6290	0.7868	0.8031	0.8957	1.0025	1.0598	1.0598
11	0.5335	0.5343	0.6290	0.6290	0.7020	0.7021	0.8957	1.0025	1.0584	1.0584
12	0.5322	0.5324	0.6290	0.6290	0.6839	0.6843	0.8957	1.0025	1.0576	1.0576
13	0.5325	0.5326	0.6290	0.6290	0.6597	0.6598	0.8957	1.0025	1.0577	1.0577
14	0.5326	0.5327	0.6290	0.6290	0.6490	0.6490	0.8957	1.0025	1.0177	1.0578
15	0.5326	0.5326	0.6290	0.6290	0.6498	0.6499	0.8957	0.9560	0.9751	1.0025
16	0.5326	0.5326	0.6290	0.6290	0.6507	0.6507	0.8957	0.9311	0.9338	1.0025
17	0.5326	0.5326	0.6290	0.6290	0.6502	0.6503	0.8957	0.9309	0.9311	1.0025
18	0.5326	0.5326	0.6290	0.6290	0.6501	0.6501	0.8957	0.9315	0.9325	1.0025
19	0.5326	0.5326	0.6290	0.6290	0.6502	0.6502	0.8957	0.9297	0.9301	1.0025
20	0.5326	0.5326	0.6290	0.6290	0.6502	0.6502	0.8957	0.9291	0.9293	1.0025
21	0.5326	0.5326	0.6290	0.6290	0.6502	0.6502	0.8957	0.9295	0.9295	1.0025
22	0.5326	0.5326	0.6290	0.6290	0.6502	0.6502	0.8957	0.9296	0.9296	1.0025
23	0.5326	0.5326	0.6290	0.6290	0.6502	0.6502	0.8957	0.9295	0.9295	1.0025
24	0.5326	0.5326	0.6290	0.6290	0.6502	0.6502	0.8957	0.9295	0.9295	1.0025
25	0.5326	0.5326	0.6290	0.6290	0.6502	0.6502	0.8957	0.9295	0.9295	1.0025
26	0.5326	0.5326	0.6290	0.6290	0.6502	0.6502	0.8957	0.9295	0.9295	1.0025
27	0.5326	0.5326	0.6290	0.6290	0.6502	0.6502	0.8957	0.9295	0.9295	1.0025
28	0.5326	0.5326	0.6290	0.6290	0.6502	0.6502	0.8957	0.9295	0.9295	1.0025
29	0.5326	0.5326	0.6290	0.6290	0.6502	0.6502	0.8957	0.9295	0.9295	1.0025

Table 4 indicates the effect of a cylindrical panel embedded with a piezoelectric layer geometry on number of convergence. According to the results, only with N_x = N_y = N = 15 grid points all of first 10 frequencies will be converged.

Table 4.

Convergence study of the first ten non-dimensional frequencies (ω_ND) of a cylindrical panel embedded with a piezoelectric layer (N_x = N_y = N).

N	$ω_{N D}^{1}$	$ω_{N D}^{2}$	$ω_{N D}^{3}$	$ω_{N D}^{4}$	$ω_{N D}^{5}$	$ω_{N D}^{6}$	$ω_{N D}^{7}$	$ω_{N D}^{8}$	$ω_{N D}^{9}$	$ω_{N D}^{10}$
5	0.8837	1.2418	1.6470	1.7844	1.9924	2.0254	2.2799	2.3024	2.3913	2.6050
6	0.6479	0.8432	0.9476	0.9706	1.7012	1.7835	2.0280	2.2637	2.4153	2.4555
7	0.5813	0.7338	0.9210	0.9361	1.1033	1.2251	1.2501	1.2808	1.5129	1.7226
8	0.5615	0.7031	0.8958	0.9209	1.0654	1.1899	1.2746	1.3048	1.3803	1.5333
9	0.5590	0.6859	0.8942	0.9168	1.0033	1.1370	1.2200	1.2598	1.3492	1.4580
10	0.5588	0.6844	0.8938	0.9166	1.0000	1.1349	1.2173	1.2214	1.2629	1.4592
11	0.5588	0.6842	0.8938	0.9166	0.9969	1.1334	1.2135	1.2205	1.2623	1.4580
12	0.5588	0.6842	0.8938	0.9166	0.9967	1.1333	1.2079	1.2205	1.2622	1.4533
13	0.5588	0.6842	0.8938	0.9166	0.9966	1.1333	1.2077	1.2204	1.2622	1.4531
14	0.5588	0.6842	0.8938	0.9166	0.9966	1.1333	1.2075	1.2204	1.2622	1.4528
15	0.5588	0.6842	0.8938	0.9166	0.9966	1.1333	1.2075	1.2204	1.2622	1.4528
16	0.5588	0.6842	0.8938	0.9166	0.9966	1.1333	1.2075	1.2204	1.2622	1.4528
17	0.5588	0.6842	0.8938	0.9166	0.9966	1.1333	1.2075	1.2204	1.2622	1.4528
18	0.5588	0.6842	0.8938	0.9166	0.9966	1.1333	1.2075	1.2204	1.2622	1.4528
19	0.5588	0.6842	0.8938	0.9166	0.9966	1.1333	1.2075	1.2204	1.2622	1.4528
20	0.5588	0.6842	0.8938	0.9166	0.9966	1.1333	1.2075	1.2204	1.2622	1.4528
21	0.5588	0.6842	0.8938	0.9166	0.9966	1.1333	1.2075	1.2204	1.2622	1.4528
22	0.5588	0.6842	0.8938	0.9166	0.9966	1.1333	1.2075	1.2204	1.2622	1.4528
23	0.5588	0.6842	0.8938	0.9166	0.9966	1.1333	1.2075	1.2204	1.2622	1.4528
24	0.5588	0.6842	0.8938	0.9166	0.9966	1.1333	1.2075	1.2204	1.2622	1.4528
25	0.5588	0.6842	0.8938	0.9166	0.9966	1.1333	1.2075	1.2204	1.2622	1.4528
26	0.5588	0.6842	0.8938	0.9166	0.9966	1.1333	1.2075	1.2204	1.2622	1.4528
27	0.5588	0.6842	0.8938	0.9166	0.9966	1.1333	1.2075	1.2204	1.2622	1.4528
28	0.5588	0.6842	0.8938	0.9166	0.9966	1.1333	1.2075	1.2204	1.2622	1.4528
29	0.5588	0.6842	0.8938	0.9166	0.9966	1.1333	1.2075	1.2204	1.2622	1.4528

In addition, the effect of piezoelectric patch will be studied. In order, the first 10 frequencies of a cylindrical shell embedded with a piezoelectric patch are provided in Table 5. When the piezoelectric material is modeled as a patch, the presence of Dirac delta terms in the system equations, which are a source of instability in the solution, causes the solution method to converge more slowly. Therefore, the frequencies will be converged with N_x = N_y = N = 27. The influence of various boundary conditions on the number of grid points required for convergence has been investigated by Motazedian et al. (2024).

Table 5.

Convergence study of the first ten non-dimensional frequencies (ω_ND) of a cylindrical shell embedded with a piezoelectric patch (N_x = N_y = N).

N	$ω_{N D}^{1}$	$ω_{N D}^{2}$	$ω_{N D}^{3}$	$ω_{N D}^{4}$	$ω_{N D}^{5}$	$ω_{N D}^{6}$	$ω_{N D}^{7}$	$ω_{N D}^{8}$	$ω_{N D}^{9}$	$ω_{N D}^{10}$
5	0.6884	0.7138	0.9272	0.9590	0.9769	1.0055	1.0107	1.2095	1.2448	1.2491
6	0.6073	0.6073	0.9279	0.9409	0.9409	0.9799	1.0146	1.2071	1.2071	1.2441
7	0.6031	0.6042	0.7638	0.9278	0.9398	0.9801	0.9854	1.0156	1.2137	1.2139
8	0.4474	0.4913	0.6038	0.6048	0.8187	0.8282	0.9278	0.9397	0.9397	0.9799
9	0.4341	0.4343	0.6037	0.6037	0.8074	0.8076	0.8842	0.9278	0.9397	0.9397
10	0.4345	0.4348	0.5069	0.5894	0.6034	0.6037	0.8052	0.8053	0.8356	0.8484
11	0.4350	0.4351	0.4432	0.4566	0.6037	0.6037	0.7826	0.7850	0.8101	0.8104
12	0.4363	0.4442	0.6038	0.6038	0.7787	0.7808	0.8063	0.8093	0.8093	0.8754
13	0.4248	0.4295	0.4354	0.4356	0.6038	0.6038	0.6696	0.6701	0.7723	0.7732
14	0.4189	0.4190	0.4357	0.4357	0.6038	0.6038	0.6207	0.6230	0.7683	0.7683
15	0.4189	0.4192	0.4357	0.4357	0.5701	0.5728	0.6038	0.6038	0.7683	0.7683
16	0.4198	0.4200	0.4357	0.4357	0.5428	0.5433	0.6038	0.6038	0.7687	0.7687
17	0.4196	0.4197	0.4357	0.4357	0.5402	0.5402	0.6038	0.6038	0.7686	0.7686
18	0.4196	0.4196	0.4358	0.4358	0.5417	0.5417	0.6038	0.6038	0.7685	0.7685
19	0.4196	0.4196	0.4358	0.4358	0.5402	0.5402	0.6038	0.6038	0.7605	0.7620
20	0.4197	0.4197	0.4358	0.4358	0.5395	0.5395	0.6038	0.6038	0.7588	0.7610
21	0.4197	0.4197	0.4358	0.4358	0.5398	0.5398	0.6038	0.6038	0.7538	0.7550
22	0.4198	0.4198	0.4358	0.4358	0.5399	0.5399	0.6038	0.6038	0.7517	0.7519
23	0.4198	0.4198	0.4359	0.4359	0.5399	0.5399	0.6038	0.6038	0.7523	0.7525
24	0.4198	0.4198	0.4359	0.4359	0.5400	0.5400	0.6038	0.6038	0.7527	0.7528
25	0.4199	0.4199	0.4359	0.4359	0.5400	0.5400	0.6038	0.6038	0.7526	0.7527
26	0.4199	0.4199	0.4359	0.4359	0.5400	0.5400	0.6038	0.6038	0.7526	0.7526
27	0.4199	0.4199	0.4359	0.4359	0.5401	0.5401	0.6038	0.6038	0.7527	0.7527
28	0.4199	0.4199	0.4359	0.4359	0.5401	0.5401	0.6038	0.6038	0.7527	0.7527
29	0.4199	0.4199	0.4359	0.4359	0.5401	0.5401	0.6038	0.6038	0.7527	0.7527

5.2. Verification

In the following, the accuracy of mathematical model needs to examine. To this aim, the results are compared with other researchers and Finite element software. First, only the natural frequency of cylindrical shell without any piezoelectric patch/layer will be studied. Therefore, the first 10 frequencies of cylindrical panel for various boundary conditions are provided in Table 6. The results are compared with Tornabene et al. (2009), Motazedian and Jafari (2023), and Abaqus software, which indicates the mathematical model has good accuracy for various boundary conditions. Furthermore, the accuracy of natural frequencies of cylindrical shell without any piezoelectric material is investigated. Table 7 indicates that the provided results for various boundary conditions have good accuracy with Tornabene et al. (2009) and Abaqus software.

Table 6.

Validation of the natural frequency (Hz) for various boundary conditions of a cylindrical panel.

ω ⁱ	C-C-C-C			C-F-F-F			S-S-F-F
ω ⁱ	Abaqus	2D–GDQ	Motazedian and Jafari (2023)	Abaqus	2D–GDQ	Tornabene et al. (2009)	Abaqus	2D–GDQ	Tornabene et al. (2009)
$ω_{N D}^{1}$	486.40	485.61	485.61	58.92	58.64	58.32	77.03	76.23	76.18
$ω_{N D}^{2}$	595.18	593.18	593.18	91.82	90.98	90.62	187.40	188.06	188.14
$ω_{N D}^{3}$	776.89	777.06	777.06	144.59	146.01	146.35	232.79	232.31	232.26
$ω_{N D}^{4}$	799.31	799.79	799.79	232.46	231.05	230.72	284.73	285.97	285.37
$ω_{N D}^{5}$	875.44	872.97	872.97	266.06	264.00	263.63	423.22	423.91	428.84
$ω_{N D}^{6}$	986.57	984.45	984.45	278.88	278.84	278.56	466.84	467.51	467.63
$ω_{N D}^{7}$	1047.22	1039.62	1039.62	338.80	339.56	339.43	535.99	537.42	537.52
$ω_{N D}^{8}$	1061.42	1062.94	1062.94	427.43	430.52	430.81	572.25	573.61	573.73
$ω_{N D}^{9}$	1094.25	1094.11	1094.11	488.07	489.11	489.26	668.43	673.12	673.55
$ω_{N D}^{10}$	1260.58	1252.80	1252.80	512.94	512.13	511.30	724.67	730.11	731.76

Table 7.

Validation of the natural frequency (Hz) for various boundary conditions of a cylindrical shell.

ω ⁱ	C-C			C-F			C-S
ω ⁱ	Abaqus	2D–GDQ	Tornabene et al. (2009)	Abaqus	2D–GDQ	Tornabene et al. (2009)	Abaqus	2D–GDQ	Tornabene et al. (2009)
ω¹	359.50	360.16	360.36	145.18	146.11	146.17	343.18	344.58	344.78
ω²	359.50	360.16	360.36	145.18	146.11	146.17	343.18	344.58	344.78
ω³	375.19	375.70	375.86	209.52	210.19	210.32	360.82	361.37	361.52
ω⁴	375.19	375.70	375.86	209.52	210.19	210.32	360.82	361.37	361.52
ω⁵	463.29	462.99	463.29	242.81	242.56	242.55	449.50	450.89	451.18
ω⁶	463.29	462.99	463.29	242.81	242.56	242.55	449.50	450.89	451.18
ω⁷	523.64	523.50	523.55	367.02	366.39	366.58	516.31	515.47	515.53
ω⁸	523.64	523.50	523.55	367.02	366.39	366.58	516.31	515.47	515.53
ω⁹	648.01	646.12	646.56	402.24	402.24	401.91	626.18	628.46	628.74
ω¹⁰	648.01	646.12	646.56	409.05	412.23	412.36	626.18	628.46	628.74

In addition, the non-dimensional frequency for a cylindrical shell embedded with piezoelectric layer with dimension of L/R = 2, h/R = 0.05, h_p = 0.2h, and θ_c/2π = 1/3 is studied. Therefore, the results of a fully clamped cylindrical shell/panel embedded with piezoelectric layer for short-circuit boundary condition are compared with Abaqus software. Table 8 indicates the natural frequencies have good accuracy. Also, the effect of open-circuit boundary condition on natural frequencies is also tabulated in Table 9. Furthermore, the validation of a cylindrical shell embedded with piezoelectric patch with length of L_p = 0.25 L for two boundary conditions of short-circuit and open-circuit boundary conditions are provided in Table 10 and Table 11, respectively. Based on the results, the mathematical model as well as 2D–GDQ method have good accuracy for modeling the dynamic behavior of cylindrical shell embedded with piezoelectric material.

Table 8.

Verification of the non-dimensional frequency for a fully clamped cylindrical shell/panel embedded with piezoelectric layer (short circuit).

$ω_{N D}^{i}$	Shell			Panel
$ω_{N D}^{i}$	2D–GDQ	Abaqus	Error (%)	2D–GDQ	Abaqus	Error (%)
$ω_{N D}^{1}$	0.9145	0.9173	0.30	1.0268	1.0240	0.27
$ω_{N D}^{2}$	0.9145	0.9173	0.30	1.0590	1.0562	0.26
$ω_{N D}^{3}$	0.9386	0.9399	0.15	1.2744	1.2728	0.13
$ω_{N D}^{4}$	0.9386	0.9401	0.16	1.5634	1.5655	0.13
$ω_{N D}^{5}$	0.9896	0.9934	0.38	1.8795	1.8783	0.06
$ω_{N D}^{6}$	0.9896	0.9934	0.38	1.9548	1.9552	0.02
$ω_{N D}^{7}$	1.0612	1.0606	0.06	2.0556	2.0624	0.33
$ω_{N D}^{8}$	1.0612	1.0610	0.03	2.1320	2.1347	0.13
$ω_{N D}^{9}$	1.1523	1.1572	0.42	2.4305	2.4369	0.27
$ω_{N D}^{10}$	1.1523	1.1572	0.42	2.4583	2.4567	0.06

Table 9.

Verification of the non-dimensional frequency for a fully clamped cylindrical shell/panel embedded with piezoelectric layer (open circuit).

$ω_{N D}^{i}$	Shell			Panel
$ω_{N D}^{i}$	2D–GDQ	Abaqus	Error (%)	2D–GDQ	Abaqus	Error (%)
$ω_{N D}^{1}$	0.9228	0.9243	0.16	1.0299	1.0312	0.12
$ω_{N D}^{2}$	0.9228	0.9243	0.16	1.0607	1.0593	0.13
$ω_{N D}^{3}$	0.9452	0.9459	0.07	1.2811	1.2808	0.03
$ω_{N D}^{4}$	0.9452	0.9459	0.07	1.5847	1.5856	0.06
$ω_{N D}^{5}$	0.9982	1.0019	0.37	1.8985	1.8964	0.11
$ω_{N D}^{6}$	0.9982	1.0019	0.37	1.9830	1.9796	0.17
$ω_{N D}^{7}$	1.0631	1.0695	0.60	2.0819	2.0858	0.19
$ω_{N D}^{8}$	1.0631	1.0698	0.62	2.1628	2.1629	0.01
$ω_{N D}^{9}$	1.1654	1.1670	0.14	2.4110	2.4062	0.20
$ω_{N D}^{10}$	1.1654	1.1670	0.14	2.4647	2.4673	0.11

Table 10.

Verification of the non-dimensional frequency for a fully clamped cylindrical shell/panel embedded with piezoelectric patch (short circuit).

$ω_{N D}^{i}$	Shell			Panel
$ω_{N D}^{i}$	2D–GDQ	Abaqus	Error (%)	2D–GDQ	Abaqus	Error (%)
$ω_{N D}^{1}$	0.9234	0.9158	0.84	1.0384	1.0349	0.34
$ω_{N D}^{2}$	0.9234	0.9158	0.84	1.0712	1.0663	0.45
$ω_{N D}^{3}$	0.9537	0.9460	0.81	1.2858	1.2834	0.19
$ω_{N D}^{4}$	0.9537	0.9460	0.81	1.5816	1.5790	0.16
$ω_{N D}^{5}$	0.9935	0.9871	0.66	1.8911	1.8859	0.27
$ω_{N D}^{6}$	0.9935	0.9871	0.66	1.9670	1.9625	0.23
$ω_{N D}^{7}$	1.0817	1.0766	0.47	2.0864	2.0778	0.42
$ω_{N D}^{8}$	1.0817	1.0766	0.47	2.1544	2.1422	0.57
$ω_{N D}^{9}$	1.1527	1.1487	0.34	2.4526	2.4464	0.25
$ω_{N D}^{10}$	1.1527	1.1487	0.34	2.4805	2.4718	0.35

Table 11.

Verification of the non-dimensional frequency for a fully clamped cylindrical shell/panel embedded with piezoelectric patch (open circuit).

$ω_{N D}^{i}$	Shell			Panel
$ω_{N D}^{i}$	2D–GDQ	Abaqus	Error (%)	2D–GDQ	Abaqus	Error (%)
$ω_{N D}^{1}$	0.9240	0.9164	0.83	1.0393	1.0359	0.33
$ω_{N D}^{2}$	0.9240	0.9164	0.83	1.0722	1.0671	0.47
$ω_{N D}^{3}$	0.9543	0.9466	0.81	1.2864	1.2842	0.17
$ω_{N D}^{4}$	0.9543	0.9466	0.81	1.5819	1.5800	0.12
$ω_{N D}^{5}$	0.9943	0.9877	0.67	1.8917	1.8883	0.18
$ω_{N D}^{6}$	0.9943	0.9877	0.67	1.9681	1.9655	0.13
$ω_{N D}^{7}$	1.0854	1.0771	0.77	2.0872	2.0788	0.40
$ω_{N D}^{8}$	1.0854	1.0771	0.77	2.1558	2.1460	0.46
$ω_{N D}^{9}$	1.1534	1.1492	0.36	2.4540	2.4494	0.19
$ω_{N D}^{10}$	1.1534	1.1492	0.36	2.4817	2.4731	0.35

5.3. Parametric studies

In this section, the effect of various vital parameter on dynamic behavior of cylindrical shell/panels embedded with piezoelectric patch are studied. To this aim, the natural frequencies of a cylindrical panel embedded with piezoelectric patch with L/R = 2, h/R = 0.05, h_p = 0.2h, L_p = 0.25 L, and θ_c/2π = 1/3 are provided in Table 12. These results are obtained for various boundary conditions. Based on the results, the non-dimensional natural frequency in the open-circuit condition is slightly higher than in the short-circuit condition. As evident, the non-dimensional natural frequency in the open-circuit condition is slightly higher than in the short-circuit condition. This observation highlights the influence of the electrical boundary condition alongside the mechanical boundary conditions on the natural frequencies. It is noteworthy, nevertheless, that the mechanical boundary conditions have a more significant impact on the natural frequencies compared to the electrical boundary conditions.

Table 12.

The effect of various boundary conditions on non-dimensional frequencies for a cylindrical panel embedded with piezoelectric patch.

$ω_{N D}^{i}$	C-C-C-C		C-C-S-S		C-F-F-F
$ω_{N D}^{i}$	Short circuit	Open circuit	Short circuit	Open circuit	Short circuit	Open circuit
$ω_{N D}^{1}$	0.4406	0.4409	0.4317	0.4320	0.0663	0.0664
$ω_{N D}^{2}$	0.4476	0.4479	0.4384	0.4390	0.1729	0.1737
$ω_{N D}^{3}$	0.6655	0.6666	0.6540	0.6532	0.2008	0.2013
$ω_{N D}^{4}$	0.6873	0.6884	0.6766	0.6786	0.2740	0.2742
$ω_{N D}^{5}$	0.7170	0.7173	0.7011	0.7020	0.3761	0.3768
$ω_{N D}^{6}$	0.8717	0.8723	0.8576	0.8562	0.3989	0.3994
$ω_{N D}^{7}$	0.9117	0.9132	0.8915	0.8890	0.4050	0.4055
$ω_{N D}^{8}$	0.9228	0.9249	0.9088	0.9081	0.5027	0.5031
$ω_{N D}^{9}$	0.9260	0.9281	0.9132	0.9135	0.6137	0.6163
$ω_{N D}^{10}$	0.9398	0.9414	0.9301	0.9293	0.6439	0.6443

Table 13 indicates the effect of various mechanical and electrical boundary condition on non-dimensional frequencies where θ = 0 to 2π. According to the results, the fully boundary condition as well as open-circuit condition has the highest value of natural frequencies.

Table 13.

The effect of various boundary conditions on non-dimensional frequencies for a cylindrical shell embedded with piezoelectric patch.

$ω_{N D}^{i}$	C-C		C-S		C-F
$ω_{N D}^{i}$	Short circuit	Open circuit	Short circuit	Open circuit	Short circuit	Open circuit
$ω_{N D}^{1}$	0.3268	0.3272	0.3193	0.3202	0.1473	0.1477
$ω_{N D}^{2}$	0.3268	0.3272	0.3193	0.3202	0.1473	0.1477
$ω_{N D}^{3}$	0.3360	0.3363	0.3312	0.3316	0.1567	0.1571
$ω_{N D}^{4}$	0.3360	0.3363	0.3312	0.3316	0.1567	0.1571
$ω_{N D}^{5}$	0.4141	0.4145	0.4051	0.4054	0.2310	0.2314
$ω_{N D}^{6}$	0.4141	0.4145	0.4051	0.4054	0.2310	0.2314
$ω_{N D}^{7}$	0.4260	0.4267	0.4237	0.4242	0.2935	0.2937
$ω_{N D}^{8}$	0.4260	0.4267	0.4237	0.4242	0.2935	0.2937
$ω_{N D}^{9}$	0.5700	0.5704	0.5632	0.5646	0.3585	0.3586
$ω_{N D}^{10}$	0.5700	0.5704	0.5632	0.5646	0.3585	0.3586

Moreover, in addition to the mechanical and electrical boundary conditions, another factor influencing the system’s natural frequency is the number of piezoelectric patches. Table 14 presents the natural frequencies of the cylindrical panel without a patch, as well as with one and two piezoelectric patches. The piezoelectric patches are located at the two opposite edges of the cylindrical panel, where a fully clamped boundary condition is applied. As can be observed, increasing the number of piezoelectric patches from one to two leads to an increase in the natural frequency. Furthermore, it is possible to simultaneously vary the boundary conditions along with the number of patches, the results of which are shown in Table 15.

Table 14.

The effect of the number of piezoelectric patches on non-dimensional frequency.

$ω_{N D}^{i}$	Without patch	One patch		Two patches
$ω_{N D}^{i}$	Without patch	Open circuit	Short circuit	Open circuit	Short circuit
$ω_{N D}^{1}$	0.4220	0.4240	0.4237	0.4265	0.4258
$ω_{N D}^{2}$	0.4290	0.4298	0.4295	0.4319	0.4311
$ω_{N D}^{3}$	0.6427	0.6408	0.6392	0.6424	0.6393
$ω_{N D}^{4}$	0.6649	0.6621	0.6609	0.6627	0.6597
$ω_{N D}^{5}$	0.6875	0.6885	0.6880	0.6906	0.6898
$ω_{N D}^{6}$	0.8394	0.8379	0.8375	0.8375	0.8367
$ω_{N D}^{7}$	0.8810	0.8777	0.8762	0.8810	0.8775
$ω_{N D}^{8}$	0.8930	0.8900	0.8870	0.8934	0.8873
$ω_{N D}^{9}$	0.8944	0.8932	0.8902	0.8939	0.8876
$ω_{N D}^{10}$	0.9100	0.9057	0.9039	0.9005	0.8978

Table 15.

The effect of boundary conditions and two piezoelectric patches on non-dimensional frequency.

$ω_{N D}^{i}$	CCCC		CCFF		CFFF		SSSF
$ω_{N D}^{i}$	Open circuit	Short circuit	Open circuit	Short circuit	Open circuit	Short circuit	Open circuit	Short circuit
$ω_{N D}^{1}$	0.4265	0.4258	0.0637	0.0635	0.0480	0.0479	0.1642	0.1642
$ω_{N D}^{2}$	0.4319	0.4311	0.1630	0.1627	0.0762	0.0759	0.3371	0.3366
$ω_{N D}^{3}$	0.6424	0.6393	0.1940	0.1938	0.1207	0.1204	0.3395	0.3394
$ω_{N D}^{4}$	0.6627	0.6597	0.2565	0.2559	0.1850	0.1847	0.4404	0.4396
$ω_{N D}^{5}$	0.6906	0.6898	0.3594	0.3592	0.1864	0.1861	0.5641	0.5639
$ω_{N D}^{6}$	0.8375	0.8367	0.3862	0.3848	0.2175	0.2162	0.5803	0.5797
$ω_{N D}^{7}$	0.8810	0.8775	0.3899	0.3891	0.3108	0.3106	0.5864	0.5850
$ω_{N D}^{8}$	0.8934	0.8873	0.4810	0.4797	0.3440	0.3432	0.6544	0.6511
$ω_{N D}^{9}$	0.8939	0.8876	0.5926	0.5923	0.3629	0.3626	0.7026	0.7018
$ω_{N D}^{10}$	0.9005	0.8978	0.6202	0.6190	0.3832	0.3827	0.7857	0.7824

In addition, effective geometry parameters of piezoelectric patch on natural frequencies of cylindrical shell are investigated. Figure 3 indicates the variations thickness of piezoelectric patch. An increase in the thickness of the piezoelectric patch generally results in a decrease in the system’s natural frequency. This behavior suggests that the thickness predominantly influences the mass matrix, therefore increasing the value of inertia and leading to a reduction in natural frequency.

Figure 3.

The effect of h_p/h on non-dimensional frequencies of cylindrical panel embedded with piezoelectric patch.

Furthermore, the variation of the natural frequency with an increase in the length of the piezoelectric patch is shown in Figure 4. Based on the results, the value of natural frequencies initially increasing and subsequently decreasing. The determination of the system’s natural frequency involves the interplay between the mass and stiffness matrices. The relative influence of the patch length on these matrices dictates the observed behavior, where the frequency trend can shift from increasing to decreasing depending on the dominant effect.

Figure 4.

The effect of L_p/L on non-dimensional frequencies of cylindrical panel embedded with piezoelectric patch.

The effect of θ of the piezoelectric patch can also be observed in Figure 5, where an increase in θ results in a decrease in the natural frequency. This suggests that an increase in θ raises the system’s inertia, which consequently leads to a reduction in the natural frequency.

Figure 5.

The effect of θ_p/360° on non–dimensional frequencies of cylindrical panel embedded with piezoelectric patch.

6. Conclusion

In this study, the free vibration of cylindrical shell/panel embedded with piezoelectric patch is investigated. The equations of motion of cylindrical shell are derived by utilizing the FSDT and Hamilton’s principle. Also, Maxwell equations are considered to model the piezoelectric patch. The 2D–GDQ method is employed to discretize the PDEs into ODEs. In addition, the Harmonic response is considered to solve eigenvalue problem to obtain the natural frequencies. The effect of various parameters on dynamic behavior of cylindrical shell/panel embedded with piezoelectric patch are studied, which vital results are provided as follows,

• According to the results, changing the geometry of piezoelectric material from patch to layer can be affected on the number of grid points for converging the natural frequencies. The geometry of cylindrical shell or panel is also effective on convergence of grid points.

• The boundary conditions of cylindrical shell, which embedded with piezoelectric material, can be effected on value of natural frequencies. The fully clamped boundary condition and open-circuit condition have the highest value of natural frequencies.

• An increase in the value of θ in the piezoelectric patch generally leads to an increase in the system’s inertia, which consequently results in a decrease in the natural frequency.

• Increasing the thickness of the piezoelectric patch can also generally lead to an increase in the system’s inertia, which consequently results in a decrease in the system’s natural frequency.

• The natural frequency of the system initially increases and then decreases with an increase in the length of the piezoelectric patch.

• An increasing the number of piezoelectric patches from one to two leads to increase the value of natural frequency.

Footnotes

ORCID iD

Nadia Motazedian

Funding

The authors received no financial support for the research, authorship, and/or publication of this article.

Declaration of conflicting interests

The authors declared no potential conflicts of interest with respect to the research, authorship, and/or publication of this article.

References

Alibeigloo

Kani

(2010) 3d free vibration analysis of laminated cylindrical shell integrated piezoelectric layers using the differential quadrature method. Applied Mathematical Modelling 34(12): 4123–4137. https://doi.org/10.1016/j.apm.2010.04.010

Arefi

Karroubi

Irani-Rahaghi

(2016) Free vibration analysis of functionally graded laminated sandwich cylindrical shells integrated with piezoelectric layer. Applied Mathematics and Mechanics 37: 821–834. https://doi.org/10.1007/s10483-016-2098-9

Asadinia

Alibeigloo

Taheri Maslak

(2025) Vibration and static analysis of sandwich cylindrical panel with carbon nanotube reinforced functionally graded core embedded in piezoelectric sensor and actuator layers. Mechanics Based Design of Structures and Machines 54: 1–28. https://doi.org/10.1080/15397734.2025.2588311

Azaripour

Baghani

(2019) Vibration analysis of fg annular sector in moderately thick plates with two piezoelectric layers. Applied Mathematics and Mechanics 40: 783–804. https://doi.org/10.1007/s10483-019-2468-8

Bhimaraddi

(1984) A higher order theory for free vibration analysis of circular cylindrical shells. International Journal of Solids and Structures 20(7): 623–630. https://doi.org/10.1016/0020-7683(84)90019-2

Bochkarev

Lekomtsev

(2024) Natural vibration and aeroelastic stability of shallow shells with passive electric circuit. International Journal of Mechanical Sciences 269: 109064. https://doi.org/10.1016/j.ijmecsci.2024.109064

Bodaghi

Shakeri

(2012) An analytical approach for free vibration and transient response of functionally graded piezoelectric cylindrical panels subjected to impulsive loads. Composite Structures 94(5): 1721–1735. https://doi.org/10.1016/j.compstruct.2012.01.009

Carrera

(1991) The effects of shear deformation and curvature on buckling and vibrations of cross-ply laminated composite shells. Journal of Sound and Vibration 150(3): 405–433. https://doi.org/10.1016/0022-460x(91)90895-q

Carrera

(1999) A study of transverse normal stress effect on vibration of multilayered plates and shells. Journal of Sound and Vibration 225(5): 803–829. https://doi.org/10.1006/jsvi.1999.2271

10.

Duryodhana

Waddar

Bonthu

, et al. (2023) Buckling and free vibrations behaviour through differential quadrature method for foamed composites. Results in Engineering 17: 100894. https://doi.org/10.1016/j.rineng.2023.100894

11.

Eftekhari

(2016) A differential quadrature procedure with direct projection of the heaviside function for numerical solution of moving load problem. Latin American Journal of Solids and Structures 13(9): 1763–1781. https://doi.org/10.1590/1679-78252251

12.

Flügge

(2013) Stresses in Shells. Springer Science & Business Media.

13.

Forsberg

(1969) Axisymmetric and beam-type vibrations of thin cylindrical shells. AIAA Journal 7(2): 221–227. https://doi.org/10.2514/3.5078

14.

Gao

Zhong

Wang

, et al. (2023) Spectral-tchebychev technique for the free and stochastic vibration analysis of functionally graded plates with piezoelectric patches. Engineering Analysis with Boundary Elements 152: 688–703. https://doi.org/10.1016/j.enganabound.2023.04.042

15.

Givi

Ghorbanpour Arani

Khoddami Maraghi

, et al. (2025) Free vibration and supersonic flutter analyses of a sandwich cylindrical shell with cnt-reinforced honeycomb core integrated with piezoelectric layers. Mechanics Based Design of Structures and Machines 53(5): 3225–3253. https://doi.org/10.1080/15397734.2024.2422413

16.

Goldman

(1974) Mode shapes and frequencies of clamped-clamped cylindrical shells. AIAA Journal 12(12): 1755–1756. https://doi.org/10.2514/3.49599

17.

Holmes

(1989) The mechanics of vibrations of cylindrical shells (stefan markus). SIAM Review 31(2): 346–347. https://doi.org/10.1137/1031080

18.

Huang

Zhang

(2024) Investigation of the free vibrations of bi-directional functionally graded material beams using differential quadrature method. Journal of Vibration and Control 30(21–22): 4828–4851. https://doi.org/10.1177/10775463231214314

19.

Kamenskikh

Lekomtsev

Matveenko

(2024) Free vibration of prestressed plates and shallow shells with piezoelectric elements. International Journal of Applied Mechanics 16(07): 2450072. https://doi.org/10.1142/s1758825124500728

20.

Khalifa

(2011) Exact solutions for the vibration of circumferentially stepped orthotropic circular cylindrical shells. Comptes Rendus Mecanique 339(11): 708–718. https://doi.org/10.1016/j.crme.2011.08.001

21.

Khorshidi

Saber-Samandari

Salehi

(2020) Free vibration response of functionally graded carbon nanotube double curved shells and panels with piezoelectric layers in a thermal environment. Scientia Iranica 27(5): 2391–2408.

22.

Luo

Sun

, et al. (2018) The influence of elastic boundary on modal parameters of thin cylindrical shell. Journal of Acoustics and Vibration 23: 93–105. https://doi.org/10.20855/ijav.2018.23.11253

23.

Zhang

Zhou

(2022) Free vibration analysis of a spinning composite laminated truncated conical shell under hygrothermal environment. Symmetry 14(7): 1369. https://doi.org/10.3390/sym14071369

24.

Hao

, et al. (2025) Analysis of free vibration and acoustic responses of fg-cntrc sandwich plates with different piezoelectric material lay-up layers. Journal of Vibration Engineering & Technologies 13(2): 139. https://doi.org/10.1007/s42417-024-01607-6

25.

Love AEH (1888) Xvi. the small free vibrations and deformation of a thin elastic shell. Philosophical Transactions of the Royal Society of London - A 179: 491–546.

26.

Motazedian

Jafari

(2023) Free vibration of joined spherical-cylindrical panel and shells using the 2d gdq method. Journal of Vibration and Sound 12(23): 156–176.

27.

Motazedian

Jafari

Montaseri

(2024) Free vibration of joined spherical–conical panel and shells using the 2d-gdq method. International Journal of Structural Stability and Dynamics 25: 2550182. https://doi.org/10.1142/s0219455425501822

28.

Pang

Liu

Cheng

, et al. (2023) Free vibration analysis of rotating functionally graded joined cylindrical-conical-cylindrical shells with arbitrary boundary conditions. In: ICMD: International Conference on Mechanical Design. Springer, pp. 2179–2191.

29.

Qing

Qiu

Liu

(2006) A semi-analytical solution for static and dynamic analysis of plates with piezoelectric patches. International Journal of Solids and Structures 43(6): 1388–1403. https://doi.org/10.1016/j.ijsolstr.2005.03.048

30.

Reddy

(2006) Theory and Analysis of Elastic Plates and Shells. CRC Press.

31.

Sanders

(1960) An Improved First-approximation Theory for Thin Shells. US Government Printing Office, Vol. 24.

32.

Sayyaadi

Farsangi

MAA

(2014) An analytical solution for dynamic behavior of thick doubly curved functionally graded smart panels. Composite Structures 107: 88–102. https://doi.org/10.1016/j.compstruct.2013.07.039

33.

Sheng

Wang

(2009) Studies on dynamic behavior of functionally graded cylindrical shells with pzt layers under moving loads. Journal of Sound and Vibration 323(3–5): 772–789. https://doi.org/10.1016/j.jsv.2009.01.017

34.

Shu

(1991) Generalized differential-integral Quadrature and Application to the Simulation of Incompressible Viscous Flows Including Parallel Computation. University of Glasgow.

35.

Shu

(2012) Differential Quadrature and its Application in Engineering. Springer Science & Business Media.

36.

Smits

Choi

(1991) The constituent equations of piezoelectric heterogeneous bimorphs. IEEE Transactions on Ultrasonics, Ferroelectrics, and Frequency Control 38(3): 256–270. https://doi.org/10.1109/58.79611

37.

Smits

Dalke

Cooney

(1991) The constituent equations of piezoelectric bimorphs. Sensors and Actuators A: Physical 28(1): 41–61. https://doi.org/10.1016/0924-4247(91)80007-c

38.

Szekrenyes

(2021) Application of differential quadrature method to delaminated first-order shear deformable composite plates. Thin-Walled Structures 166: 108028. https://doi.org/10.1016/j.tws.2021.108028

39.

Szekrenyes

(2022) Differential quadrature solution for composite flat plates with delamination using higher-order layerwise models. International Journal of Solids and Structures 248: 111621. https://doi.org/10.1016/j.ijsolstr.2022.111621

40.

Tornabene

Viola

Inman

(2009) 2-d differential quadrature solution for vibration analysis of functionally graded conical, cylindrical shell and annular plate structures. Journal of Sound and Vibration 328(3): 259–290. https://doi.org/10.1016/j.jsv.2009.07.031

41.

Tran

Nguyen

(2023) Free vibration and dynamic responses of bidirectional functionally graded plates with piezoelectric layers under moving mass. In: The International Conference on Sustainable Civil Engineering and Architecture. Springer, pp. 1394–1402.

42.

Umeda

Nakamura

Ueha

(1997) Energy storage characteristics of a piezo-generator using impact induced vibration. Japanese Journal of Applied Physics 36(5S): 3146. https://doi.org/10.1143/jjap.36.3146

43.

Wang

Abbaspour

Hosseini

(2025) A chebyshev–ritz formulation for vibration analysis of piezoelectric composite doubly curved micro panels with arbitrary boundary conditions. Acta Mechanica 236(10): 6283–6298. https://doi.org/10.1007/s00707-025-04459-7

44.

Weiwei

Siqi

Xiaojun

, et al. (2024) Vibration analysis of porous functionally graded material truncated conical shells in axial motion. Mechanics 30(2): 123–134.

45.

, et al. (2025) The vibration characteristics analysis of fiber-reinforced thin-walled conical–cylindrical composite shells with artificial spring technique. International Journal of Structural Stability and Dynamics 25(19): 2550200. https://doi.org/10.1142/s0219455425502001

46.

Zhou

Zheng

Wang

, et al. (2024) Dynamic characteristics of composite damping sandwich open conical shell. Mechanics of Advanced Materials and Structures 31(5): 1152–1166. https://doi.org/10.1080/15376494.2022.2132435