Structural-Electrical-Coupled Formulation for the Free Vibration of a Piezoelectric-Laminated Plate Using the Analytical Arbitrary Quadrilateral p Element

Abstract

An analytical quadrilateral p element is developed for solving the free vibrations of piezoelectric-laminated plates. The formulations of the displacement and strain fields are based on first-order shear deformation plate theory. The coupling effect between the electrical and stress fields is also considered. The Legendre orthogonal polynomials are used as the element interpolation functions, and the analytical integration technique is adopted. It is found that the present p element method gives high numerical precision results, fast and monotonic convergence rate. In the numerical cases, the effects of the number of hierarchical terms and mesh size on the convergence rate are investigated. Examples of square plates with different displacement and potential boundary conditions are studied. In the comparisons, the solutions of the present element are in good agreement with those obtained from other classical and finite element methods.

1. Introduction

Various beam/plate problems have been hot research topics for decades (e.g., [1–15]), and many of them are solved using the finite element method (e.g., [16–21]). Beams/plates embedded with piezoelectric materials are widely used in various intelligent material systems for active control because of the electromechanical-coupling properties. For cases of simple geometric shapes, the two-dimensional (2D) and three-dimensional (3D) analytical methods that were developed by Benjeddou and Deü [22], Correia et al. [23], Heyliger and Saravanos [24], and Ding and Chen [25] are effective in analyzing the free vibrations of piezoelectric-laminated plates and provide clear physical meanings in the analytical and symbolical formulations. For cases of complex geometries, boundary conditions, and loadings, the finite element method (FEM) is more powerful and versatile than classical solutions. Most of the classical continuum solutions of plates have been limited to single- or two-mode approximations. This is due to the difficulties in obtaining the general multiple mode governing equations using Galerkin’s approach, especially for plates with complex boundary conditions. Hence, Saravanos et al., Lam et al., He et al., and Loja et al. [26–30] developed the finite-element and finite-strip methods for complicated piezoelectric-laminated plates. There are two main kinds of theories used for finite-element formulations. One is the classical plate theory (CPT), and the other one is the shear deformation theory, which branches out into first-order shear deformation theory (FSDT) and higher order shear deformation theory (HSDT). In general, theories that consider the shear deformation effect are more accurate for the analysis of piezoelectric-laminated plates. Thus, one such theory is used in this paper. The performance of an FEM model can be improved by a finer mesh (h-version) or refined integrated piezoelectric sensor and actuator interpolation (p-version). Zienkiewicz and Taylor [31] mentioned that the convergence rates of p-version elements would be more rapid than those of h-version elements for the same degrees of freedom (DOFs). The p-version elements that were developed by Houmat [32] and Woo et al. [33] were demonstrated to have fast monotonic convergence for cases of triangular plates and thick skew plates. According to the free-vibration analysis of skew Mindlin plates that Leung and Zhu [34] conducted, if the p-version approach is employed, then analytical integration is preferred in the procedures of setting up the element formulation; otherwise, the numerical integration errors may adversely affect the accuracy, and monotonic convergence of the natural frequencies cannot be guaranteed. In the axial free vibration analyses of beams, Leung and Chan [35] adopted trigonometric interpolation functions and the analytical integration technique to improve the accuracy. The analytical p element approach has been further applied to various vibrating plate problems [36–39] such as membrane vibration and thick plates and laminated plates. Although this kind of element can be applied to the analysis of plates with triangular and other complicated shapes, it is tedious to set up the mesh and assemble the elements when compared with quadrilateral elements.

Based on FSDT, the quadrilateral hierarchical element is developed for the free vibrations of piezoelectric-laminated plates in this paper. Legendre orthogonal polynomials are adopted in the interpolation functions to improve the accuracy. Comparisons with published results and case studies show the good accuracy and efficiency of the present quadrilateral p element.

2. Theoretical Formulation

2.1. Analytical Model

Figure 1 shows a laminated plate and the coordinate system defined at the midplane. The fiber direction is indicated by an angle 𝜃, which is the positive-rotation angle of the principal material axes from the arbitrary 𝑥𝑦 axes. The elasticity modulus for a layer parallel to the fibers is 𝐸₁₁, and perpendicular to the fibers is 𝐸₂₂. The layers can be piezoelectric or of composite materials. It is assumed that the potential of the piezoelectric layer varies linearly through its thickness.

Figure 1 

The configuration of a laminated composite plate.

2.2. Displacements and Strains

According to FSDT [40], the displacement fields of the plate are expressed as:𝑢𝑥,𝑦,𝑧=𝑢𝑥,𝑦+𝑧𝜓_𝑥𝑥,𝑦,𝑣𝑥,𝑦,𝑧=𝑣𝑥,𝑦+𝑧𝜓_𝑦𝑥,𝑦,𝑤𝑥,𝑦,𝑧=𝑤𝑥,𝑦, (2.1) where 𝑢,𝑣, and 𝑤 are the displacement components along the 𝑥, 𝑦, and 𝑧 axes, respectively, 𝑢, 𝑣, and 𝑤 are the associated midplane displacements, and 𝜓_𝑥, 𝜓_𝑦 are the normal rotations about the 𝑦 and 𝑥 axes, respectively.

The linear strain-displacement relations are given by𝜀_𝑥=𝜕𝑢𝜕𝑥+𝑧𝜕𝜓_𝑥,𝜀𝜕𝑥_𝑦=𝜕𝑣𝜕𝑦+𝑧𝜕𝜓_𝑦,𝛾𝜕𝑦_𝑥𝑦=𝜕𝑢+𝜕𝑦𝜕𝑣𝜕𝑥+𝑧𝜕𝜓_𝑥𝜕𝑦+𝑧𝜕𝜓_𝑦,𝛾𝜕𝑥_𝑦𝑧=𝜕𝑤𝜕𝑦+𝜓_𝑦,𝛾_𝑥𝑧=𝜕𝑤𝜕𝑥+𝜓_𝑥. (2.2)

The electrical potential is assumed to be linear through the thickness in each piezoelectric layer. Hence, the potential of the kth layer 𝜙 can be expressed as𝜙=ℎ𝑥,𝑦,𝑧^𝑘+1−𝑧ℎ^𝑘+1−ℎ^𝑘𝜙_𝑘𝑥,𝑦+𝑧−ℎ^𝑘ℎ^𝑘+1−ℎ^𝑘𝜙_𝑘+1𝑥,𝑦, (2.3) where 𝜙_𝑘,𝜙_𝑘+1 are the electrical potentials of the bottom and top surfaces of the kth layer, respectively, and ℎ^𝑘 and ℎ^𝑘+1 are the 𝑧-coordinates corresponding to the bottom and top surfaces.

Therefore, the electric field in the 𝑘th layer can be written as: 𝐄^𝑘1=−∇𝜙=−ℎ^𝑘+1−ℎ^𝑘ℎ^𝑘+1𝜕−𝑧𝜕𝑥𝑧−ℎ^𝑘𝜕ℎ𝜕𝑥_𝑘+1𝜕−𝑧𝜕𝑦𝑧−ℎ^𝑘𝜕𝜙𝜕𝑦−11^𝑘𝜙_𝑘+1. (2.4)

In two-dimensional problems, the transverse displacement 𝑤 can be interpolated using the following equation (2.5):𝑤𝜉,𝜂=^𝑝+2_𝑘=1^𝑞+2_𝑙=1𝑤_𝑘,𝑙𝑁_𝑘,𝑙𝜉,𝜂, (2.5) where 𝑁_𝑘,𝑙(𝜉,𝜂)=𝑓_𝑘(𝜉)𝑓_𝑙(𝜂); 𝑓_𝑘(𝜉) and 𝑓_𝑙(𝜂) are 𝐶⁰ Legendre orthogonal polynomials [40] (which are also used in the interpolations of displacements 𝑢, 𝑣, 𝜓_𝑥, 𝜓_𝑦).

When 𝑘, 𝑙=1 and 2, (2.5) represents a typical FEM interpolation. The hierarchical shape functions, when 𝑘 or 𝑙>2, lead to zero displacement at the corner nodes. Additional DOFs appear along the four edges and in the interior of the element. The DOFs at the four corner nodes are represented by 𝑘 and 𝑙≤2, and the DOFs along the four edges are represented by 𝑘 or 𝑙>2. Finally, the DOFs in the interior are represented by 𝑘 and 𝑙>2.

2.3. Constitutive Relationships

The constitutive relationships for the lamina-oriented arbitrarily, taking into account the piezoelectric effects are given as follows [25]:𝝈=𝐐𝜺−𝐞^𝐓𝐄,𝐃=𝐞𝜺+𝐩𝐄, (2.6) where 𝝈 is the stress matrix, 𝜺 is the strain vector, 𝐃 is the electrical displacement vector, 𝐄 is the electric field vector, 𝐐 is the transformed elastic stiffness coefficient matrix, the coefficients of which are explicitly given by Reddy [41], 𝐩 is the permittivity matrix, and 𝐞 is the transformed piezoelectric stress coefficient matrix, which is given by𝑒₃₁=𝑒₃₁cos²𝜃+𝑒₃₂sin²𝜃,𝑒₃₂=𝑒₃₁sin²𝜃+𝑒₃₂cos²𝜃,𝑒₁₄=𝑒₁₅−𝑒₂₄sin𝜃cos𝜃,𝑒₁₅=𝑒₁₅cos²𝜃+𝑒₂₄sin²𝜃,𝑒₂₄=𝑒₂₄cos²𝜃+𝑒₁₅sin²𝜃,𝑒₂₅=𝑒₁₅−𝑒₂₄sin𝜃cos𝜃,𝑒₃₆=𝑒₃₁−𝑒₃₂sin𝜃cos𝜃, (2.7) where 𝑒_𝑖𝑗 are the piezoelectric stress coefficients in the directions parallel and perpendicular to the fibers, and 𝜃 is the fiber direction.

2.4. Finite Element Model

The coordinate systems that are used to define an arbitrary quadrilateral plate element are shown in Figures 1 and 2 depict the Cartesian plane coordinate system of the element and the mapped 𝜉−𝜂 square plane region. The Jacobian matrix is expressed in terms of the Cartesian coordinates at the four corner nodes: 𝐉=𝜕𝑥𝜕𝜉𝜕𝑦𝜕𝜉𝜕𝑥𝜕𝜂𝜕𝑦=𝜕𝜂𝑎+𝑏𝜂𝑑+𝑒𝜂𝑐+𝑏𝜉𝑓+𝑒𝜉, (2.8) where 𝑎=0.25(𝑥₂+𝑥₃−𝑥₄), 𝑏=0.25(−𝑥₂+𝑥₃−𝑥₄), 𝑐=0.25(−𝑥₂+𝑥₃+𝑥₄), 𝑑=0.25(𝑦₂+𝑦₃−𝑦₄), 𝑒=0.25(−𝑦₂+𝑦₃−𝑦₄), and 𝑓=0.25(−𝑦₂+𝑦₃+𝑦₄).

(a)

(b)

Figure 2 

The coordinate transformation of a quadrilateral element.

The determinant of the Jacobian matrix is |𝐉|=(𝑎𝑒−𝑏𝑑)𝜉+(𝑏𝑓−𝑐𝑒)𝜂+𝑎𝑓−𝑐𝑑, and𝐉⁻¹=1|𝐉|−()𝑓+𝑒𝜉−𝑑+𝑒𝜂𝑐+𝑏𝜉𝑎+𝑏𝜂. (2.9)

The displacement fields 𝑢,𝑣,𝑤,𝜓_𝑥, and 𝜓_𝑦 and the electric potential 𝜙^𝑘 are interpolated using the 𝐶⁰ Legendre orthogonal polynomials as𝑢𝑣𝑤𝜓𝐮=_𝑥𝜓_𝑦=𝑁_1,1𝐈,𝑁_1,2𝐈,…,𝑁_1,𝑞+2𝐈,𝑁_2,1𝐈,…,𝑁_{𝑝+2,𝑞+2}𝐈𝜹^𝑒=𝐍_𝑢⋅𝜹^𝑒,𝜙_𝑘=𝑁_1,1,𝑁_1,2,…,𝑁_1,𝑞+2,𝑁_2,1,…,𝑁_{𝑝+2,𝑞+2}𝜙^𝑒_𝑘=𝐍_𝜙⋅𝜙^𝑒_𝑘, (2.10) where 𝜹^𝑒 is the vector of generalized DOFs, 𝜙^𝑘_𝑒 is the nodal electric potential in the bottom surface of the kth layer, and 𝐈 is the 5×5 identity matrix. Substituting the expressions of the displacement fields into (2.2) and (2.4) gives𝜀𝜺=^𝑥𝜀_𝑦𝛾_𝑥𝑦𝛾_𝑦𝑧𝛾_𝑥𝑧=𝜕𝜕𝑥00𝑧𝜕00𝜕𝜕𝑥𝜕𝑦00𝑧𝜕𝜕𝜕𝑦𝜕𝜕𝑦0𝜕𝑥𝑧𝜕𝜕𝑦𝑧𝜕𝜕𝜕𝑥00𝜕𝜕𝑥1000𝜕𝑦01𝐮=𝐁_𝑢⋅𝜹^𝑒,𝐄^𝑘=𝐸_𝑥𝐸_𝑦𝐸_𝑧1=−ℎ^𝑘+1−ℎ^𝑘ℎ^𝑘+1𝜕−𝑧𝜕𝑥𝑧−ℎ^𝑘𝜕ℎ𝜕𝑥_𝑘+1𝜕−𝑧𝜕𝑦𝑧−ℎ^𝑘𝜕𝜙𝜕𝑦−11^𝑘𝜙_𝑘+1=−𝐁_𝜙⋅𝚽^𝑒, (2.11) where Φ^𝑒=[𝜙^𝑒₁,…,𝜙^𝑒_𝑘,𝜙^𝑒_𝑘+1,…,𝜙^𝑒_𝑛+1]^𝐓, 𝑛 is the total number of layers.

The potential energies of the element include the elastic strain energy, the piezoelectric energy, and the electrical energy, which are given as 1𝐻=2^𝑛_𝑘=1^ℎ𝑘+1_ℎ𝑘𝜺^𝑇𝐐𝜺−2𝜺^𝑇𝐞^𝑇𝐄−𝐄^𝑇=1𝐩𝐄𝑑𝑥𝑑𝑦𝑑𝑧2𝜹^𝑛𝑒𝑇_𝑘=1^ℎ𝑘+1_ℎ𝑘𝐁^𝑇_𝑢𝐐𝐁_𝑢𝑑𝑥𝑑𝑦𝑑𝑧⋅𝜹^𝑒+𝜹^𝑛𝑒𝑇_𝑘=1^ℎ𝑘+1_ℎ𝑘𝐁^𝑇_𝑢𝐞^𝑇𝐁_𝜙𝑑𝑥𝑑𝑦𝑑𝑧⋅𝚽^𝑒−12𝚽^𝑛𝑒𝑇_𝑘=1^ℎ𝑘+1_ℎ𝑘𝐁^𝑇_𝜙𝐩𝐁_𝜙𝑑𝑥𝑑𝑦𝑑𝑧⋅𝚽^𝑒. (2.12)

The kinetic energy of the element is given by𝑇^𝑒=12^𝑛_𝑘=1^ℎ𝑘+1_ℎ𝑘𝜌𝑢̇²+𝑣²̇+𝑤̇²+𝑧²𝜌𝜓̇²_𝑥+𝜓̇²_𝑦𝑑𝑥𝑑𝑦𝑑𝑧. (2.13)

The Hamilton principle is adopted in the derivation of the equation of the eigenvalue problem, which is as follows:𝐊_𝑢𝑢𝐊_𝑢𝜙𝐊_𝜙𝑢𝐊_𝜙𝜙−𝜔²𝐌_𝑢𝑢𝟎𝜹𝟎𝟎𝜙=𝟎, (2.14) where𝐊_𝑢𝑢=_𝑒𝐊^𝑒_𝑢𝑢,𝐊_𝑢𝜙=_𝑒𝐊^𝑒_𝑢𝜙,𝐌_𝑢𝑢=_𝑒𝐌^𝑒_𝑢𝑢,𝐊_𝜙𝑢=𝐊^𝑇_𝑢𝜙,𝐊^𝑒_𝑢𝑢=^𝑛_𝑘=1^ℎ𝑘+1_ℎ𝑘𝐁^𝑇_𝑢𝐐𝐁_𝑢𝑑𝑥𝑑𝑦𝑑𝑧=^𝑛_𝑘=1^ℎ𝑘+1_ℎ𝑘¹₋₁¹₋₁𝐁^𝑇_𝑢𝐐𝐁_𝑢⋅|𝐉|𝐊𝑑𝜉𝑑𝜂𝑑𝑧,_𝑒𝑢𝜙=^𝑛_𝑘=1^ℎ𝑘+1_ℎ𝑘𝐁^𝑇_𝑢𝐞^𝑇𝐁_𝜙𝑑𝑥𝑑𝑦𝑑𝑧=^𝑛_𝑘=1^ℎ𝑘+1_ℎ𝑘¹₋₁¹₋₁𝐁^𝑇_𝑢𝐞^𝑇𝐁_𝜙⋅|𝐉|𝐊d𝜉𝑑𝜂𝑑𝑧,_𝑒𝜙𝜙=^𝑛_𝑘=1^ℎ𝑘+1_ℎ𝑘𝐁^𝑇_𝜙𝐩𝐁_𝜙𝑑𝑥𝑑𝑦𝑑𝑧=^𝑛_𝑘=1^ℎ𝑘+1_ℎ𝑘¹₋₁¹₋₁𝐁^𝑇_𝜙𝐩𝐁_𝜙⋅|𝐉|𝐌𝑑𝜉𝑑𝜂𝑑𝑧,_𝑒𝑢𝑢=^𝑛_𝑘=1^ℎ𝑘+1_ℎ𝑘𝜌⋅𝐍^𝑇diag1,1,1,𝑧²,𝑧²=𝐍𝑑𝑥𝑑𝑦𝑑𝑧,_{𝑛𝑘=1ℎ𝑘+1ℎ𝑘1−11−1}𝜌⋅𝐍^𝑇diag1,1,1,𝑧²,𝑧²𝐍⋅|𝐉|𝑑𝜉𝑑𝜂𝑑𝑧. (2.15)(2.16) In the procedure of coordinate mapping, the derivative in 𝐁_𝑢 and 𝐁_𝜙 should be replaced by𝜕𝜕𝜕𝑥𝜕𝑦=𝐉⁻¹𝜕𝜕𝜕𝜉𝜕𝜂. (2.17) Note that analytic integration is adopted in (2.16). Hence, the analytic integrations of the stiffness and mass matrices are in the form of 𝜉^𝑖𝜂^𝑗/(𝐴𝜉+𝐵𝜂+𝐶), where 𝐴, 𝐵, and 𝐶 are constants.

3. Numerical Results and Discussions

In this section, case studies are carried out to show the performance of the quadrilateral p element. The elastic, piezoelectric, and dielectric material properties that are adopted in the case studies are given in Table 1. There are two sets of electric boundary conditions for the outer surface of the piezoelectric layers: a closed circuit condition Φ=0, (i.e., the potential is grounded) and an open circuit condition 𝐷_𝑧=0, (i.e., the electric displacement is zero). Unless otherwise stated, the number of the hierarchical term p is 4, and the mesh in Figure 3(c) is used.

Table 1 

Elastic, piezoelectric, and dielectric properties of the materials.

(a)

(b)

(c)

Figure 3 

The 3 meshes adopted for the square plate.

3.1. Convergence Study of a Square Plate

The free vibration of a simply supported hybrid sandwich square plate is investigated here. The accuracy of the present element can be seen in this convergence study. The lamina configuration is PZT-4/GE/PZT-4. The thicknesses of the layers are 0.1 h, 0.8 h, 0.1 h, respectively. The close circuit condition is considered. The plate is meshed into one, two, or four quadrilateral elements, as shown in Figure 3. The fundamental frequency parameter is shown in Table 2 for different numbers of hierarchical terms. It can be observed that the monotonic convergence rate is very fast with respect to the number of hierarchical terms, and the results are in good agreement with those obtained from the 2D analytical method [22] and the 3D state space method [42]. It should be noted that the electrical potentials in [22, 42] are assumed to be quadratic through the plate thickness, whereas they are linear in the present method. Thus, there is a small difference between the results.

Table 2 

Fundamental frequency parameter 𝜆(=𝜔𝑎²𝜌^1/2/2000𝜋ℎ) of a simply supported sandwich square plate (Φ=0,𝑎/ℎ=10,𝜌=𝜌_𝑖ℎ_𝑖/ℎ).

3.2. Square Plates with Different Displacement Boundary Conditions

The first ten frequencies of a three-ply piezoelectric-laminated square plate (G-1195N/Ti-6A1-4V/G-1195N, 0.1 mm/0.48 mm/0.1 mm) are computed for three boundary conditions: simply supported (SSSS), fully clamped (CCCC), and cantilevered (CFFF). The electric boundary condition is a closed circuit. He et al. [28] derived finite element solutions using a traditional laminated plate element and CPT. A comparison of the results of our element and those of the traditional elements is shown in Table 3. Besides, a comparison for a static case is shown in Figure 4. The transverse normal stress, transverse displacement, and in-plane stress results of the sensor plate from the proposed method reasonably agree with those from [8].

Table 3 

Natural frequencies (Hz) of a three-ply square plate (𝑎=0.4𝑚,Φ=0).

(a)

(b)

(c)

Figure 4 

The comparisons between (a) the transverse normal stress, (b) transverse displacement, and (c) in-plane stress results of a sensor plate from the proposed method and [8].

3.3. Square and Trapezoidal Plates with Different Potential Boundary Conditions

To study the validity of the present element for different potential boundary conditions, a five-ply simply supported plate, which is made of three plies of G/E with the laminate sequence (0°, 90°, 0°)and two surface-bonded PZT-4 piezoelectric layers, is considered. The thickness of both piezoelectric layers is 0.1 h, and all G/E layers have the same thickness. The unit density is adopted for fair comparison. The solutions of the present method are presented in Table 4 and are compared with the 3D exact solution. Besides, Figure 5 shows the trapezoidal plates with PZT-5H layers perfectly bonded on both the upper and lower surfaces having various boundary conditions (the number of Fourier terms used are 𝑝=𝑞=5). The thicknesses of the trapezoidal plates and PZT-5H layers are 10 mm and 1 mm, respectively. Table 5 shows the natural frequencies of the aluminum trapezoidal plates with PZT-5H layers.

Table 4 

Nondimensional frequencies (𝜔𝑎²𝜌^1/2/2000𝜋ℎ) of a simple supported square plate with different potential boundary conditions.

Table 5 

Natural frequencies (Hz) of the aluminum trapezoidal plates with a top and a bottom PZT-5H layers.

(a)

(b)

(c)

Figure 5 

The 3 meshes adopted for the 10 mm aluminum trapezoidal plates with 1 mm PZT-5H layers perfectly bonded on both the upper and lower surfaces having various boundary conditions (the number of Fourier terms: 𝑝=𝑞=5, 𝐿𝑎=0.1m).

4. Conclusion

A quadrilateral p element with analytical integration for the free vibration of piezoelectric-laminated composite plates is presented. The Legendre orthogonal polynomials are used as the element interpolation functions, and the analytical integration technique is adopted. The monotonic convergence rate of the present element is very fast with respect to the number of hierarchical terms. Comparisons between the solutions that are obtained from our method and those of other methods show that the element is accurate and efficient for free-vibration analyses of piezoelectric-laminated plates.

Acknowledgment

The research described herein was supported by Grant number CityU SRG7008100.