# Numerical and Experimental Analysis of the Dynamic Behavior of Piezoelectric Stiffened Composite Plates Subjected to Airflow

Authors: Nguyen Thai Chung; Nguyen Ngoc Thuy; Duong Thi Ngoc Thu; Le Hai Chau
Journal: Mathematical Problems in Engineering
 (2019)

---

## Abstract

In this paper, the authors present results on dynamic behavior analysis of the stiffened composite plate with piezoelectric patches under airflow by finite element method and experimental study. The first-order shear deformation plate theory and nine-noded isoparametric piezoelectric laminated plate finite element with five elastic degrees of freedom at each node and one electric degree of freedom per element per piezoelectric layer were used in the dynamic analysis of plates by finite element method. The modern equipment was used in the dynamic behaviors analysis of plates subjected to airflow load by experimental method. In this study, the results of the theoretical method have been compared with experimental studies.


---

## Body

## 1. Introduction

The research and development of smart stiffened composite structures represent one of the most significant recent trends in the mechanics of structures, especially in the aerospace industry. These structures can provide significant advantages over traditional structures, notably active vibration control. In 1957 and 1959, the first investigations of piezoelectric structures were published by Haskins and Walsh [1] and Toupin [2] who considered cylindrical and spherical shells, respectively. Tiersten’s [3] analysis was in linear vibrations of piezoelectric plates. Adelman and Stavsky [4] studied vibrations of composite cylinders by piezoceramic and metallic layers. V. Balamurugan, S. Narayanon and Rajan L.Wankhade, Kamal M.B [5, 6] studied consistency technique for vibration control of smart stiffened plates using distributed piezoelectric sensors and actuators subjected to cyclic loads. P.veera Sanjeeva Kumar, B.Chandra Mohana Reddy [7] have presented a paper on the vibration of smart composite laminate plates using higher order theory. Sangamesh B. Herakal, Sai Kumar Dathrika, P. Giriraj Goud, and P. Sravan [8] analyzed bending of the smart composite plate under thermal environment, and Kanjuro Makihara, Junjiro Onoda, and Kenji Minesugi [9] studied flutter of cantilevered plate wing with piezoelectric material layers by finite element method. In 2016, Zafer K and Zahit M [10] presented a paper on flutter analysis of a laminated composite plate with temperature dependent material properties. N. T. Chung, H. X. Luong, and N. T. T. Xuan [11] studied the dynamic stability of laminated composite plate with piezoelectric layers. N. N. Thuy and N. T. Chung [12] have presented a paper on the dynamic analysis of smart stiffened composite plates using higher order theory.In this paper, in order to have more studies about the dynamic responses of piezoelectric composite plates, the authors examine the problems with piezoelectric stiffened composite plates subjected to airflow by finite element and experimental method.


## 2. Finite Element Formulation and the Governing Equations

Consider isoparametric piezoelectric laminated stiffened plate with the general coordinate system (x, y, z), in which the x, y plane coincides with the neutral plane of the plate. The top surface and lower surface of the plate are bonded to the piezoelectric patches (actuator and sensor). The plate is subjected to the airflow load acting (Figure1).Figure 1
Smart stiffened plate and coordinate system (a), and lamina details (b).
(a) (b)


### 2.1. Piezoelectric Equations

Consider a laminated composite plate with integrated sensors and actuators as shown in Figure1. It is assumed that each layer of the plate possesses a plane of elastic symmetry parallel to the x-y plane; the kth layer’s lamina constitutive equations coupling the direct and converse piezoelectric equations can be expressed as [13, 14](1)σ1σ2τ12k=Q-11Q-120Q-12Q-22000Q-66kε1ε2γ12k-00e-3100e-32000kE1E2E3k,(2)D1D2D3k=000000e-31e-320kε1ε2γ12k+p-11000p-22000p-33kE1E2E3k,where Q-ij,e-ij,p-ij are, respectively, the plane-stress reduced elastic constants, the piezoelectric constants, and the permittivity coefficients of the kth lamina in its material coordinate system; σk=σ1σ2τ12kT,  εk=ε1ε2γ12kT,  Ek=E1E2E3kT,  Dk=D1D2D3kT are the stress, strain, electric field, and electric displacement matric components, respectively, to the material coordinate system. If the voltage is applied to the actuator in the thickness only, then Ek=00-1/hkTϕk=Bϕϕk, Vk is the applied voltage across the kth ply, and hk is the thickness of the kth piezoelectric layer. The plane-stress reduced elastic constants Q-ij are given as(3)Q-11=E11-ν12ν21,Q-12=ν12E21-ν12ν21,Q-22=E21-ν12ν21,Q-66=G12,Upon transformation, the lamina piezoelectric equations can be expressed in terms of the stress, strains, and electric displacements in the plate coordinates as [13, 14](4)σxσyτxyk=Q11Q12Q16Q12Q22Q26Q16Q26Q66kεxεyγxyk-00e3100e3200e36kExEyEzk,(5)DxDyDzk=000000e31e32e36kεxεyγxyk+p11p120p12p22000p33kExEyEzk.Equations (4) and (5) can also be written as(6)σk=Qkεk-ekTEk,(7)Dk=ekεk+pkEk,TheQij are the components of the transformed lamina stiffness matrix which are defined as follows:(8)Q11=Q-11cos4θ+Q-22sin4θ+2Q-12+2Q-66sin2θcos2θ,Q12=Q-11+Q-22-4Q-66sin2θcos2θ+Q-12sin4θ+cos4θ,Q22=Q-11sin4θ+Q-22cos4θ+2Q-12+2Q-66sin2θcos2θ,Q16=Q-11-Q-12-2Q-66cos3θsin⁡θ-Q-22-Q-12-2Q-66sin3θcos⁡θ,Q26=Q-11-Q-12-2Q-66sin3θcos⁡θ-Q-22-Q-12-2Q-66sin⁡θcos3θ,Q66=Q-11+Q-12-2Q-12-2Q-66sin2θcos2θ+Q-66sin4θ+cos4θ,where θ is the lamina orientation angle.The piezoelectric constants matrix [e] is unavailable and it can be expressed by the more commonly available piezoelectric strain constant matrix which is [d] as [13](9)e=dQ,where d=000d15000d2400d31d32000.


### 2.2. Smart Piezoelectric Composite Plate Element Formulation

A smart piezoelectric composite plate element is considered with coordinates x, y along the in-plane direction and z along the thickness direction as shown in Figure2.Figure 2
Smart piezoelectric composite plate element.Using the Mindlin formulation, the displacements u, v, and w at a point (x,y,z) form the median surface and are expressed as functions of mid-plane displacements u0, v0, and w0 and independent rotations θx and θy of the normal in xz and yz planes, respectively, as [11, 15, 16](10)ux,y,z,t=u0x,y,t+zθyx,y,t,vx,y,z,t=v0x,y,t-zθxx,y,t,wx,y,z,t=w0x,y,t,The components of the strain vector corresponding to the displacement field (10) are [7, 11](11)εxεyγxy=εbln+εnl=εb,where {εbln},{εbnl} are the linear and nonlinear strain vector, respectively, and defined as follows:(12)εbln=εxεyγxy=εxoεyoγxyo+zκxκyκxy=∂∂x00∂∂y∂∂y∂∂xu0v0+z-∂∂y00-∂∂x-∂∂y∂∂xθxθy=ε0+zκ,εbnl=12∂w0∂x00∂w0∂y∂w0∂y∂w0∂x∂∂x∂∂yw0.In the finite formulation, the displacement fieldu and the electric potential ϕ over an element are related to the corresponding node values ue of the element and the electric charges of the piezo layer ϕe by the mean of the shape functions [Nu], [Nϕ], as follows [15]:(13)u=Nuqe,(14)ϕ=Nϕϕe.Equation (11) can be expressed as(15)ε=Bulnue+Bunlue.where Buln and Bunl are linear and nonlinear strain matric, respectively, Nu are the shape functions, and Nϕ are the shape functions for the electric potential.The dynamic equations of a finite smart laminated composite plate can be derived by using Hamilton’s principle [15, 17]:(16)δ∫t1t2Te-Ue+Wedt=0,where Te is the kinetic energy, Ue is the strain energy, and We is the work done by the applied forces.The kinetic energy at the element level is defined as(17)Te=12∫Veρu˙Tu˙dV,where Ve is the volume of the plate element.The strain energy can be written as(18)Ue=12∫VeεTσdV,The work done by the external forces(19)We=∫VeuTfbdV+∫SeuTfsdS+uTfc,where {fb} is the body force, Se is the surface area of the plate element, {fs} is the surface force, and {fc} is the concentrated load.Substituting (17), (18), and (19) into (16), noting that the electric field vector E=-Bϕϕe, and using (6), (7), (13), (14), and (15), the dynamic matrix equations can be written as(20)Meu¨e+Kuulne+Kuunleue+Kuϕeϕe=fem,(21)Kϕueue+Kϕϕeϕe=Φepz,Substitute (21) into (20) to obtain(22)Meu¨e+Kuulne+Kuunle-KuϕeKϕϕe-1Kϕueue=fem-KuϕeKϕϕe-1Φepz,where Me=∫VeρNuTNudV, Kuulne=∫VeBulnTQBulndV, Kuunle=∫VeBuunlTQBuunldV, Kuϕe=∫VeBuTeTBϕdV, Kϕϕe=∫VeBϕTpBϕdV, Kϕue=KuϕeT, which are, respectively, the element mass, linear mechanical stiffness, nonlinear mechanical stiffness, mechanical-electrical coupling, piezoelectric permittivity, and electrical-mechanical coupling matrix. The vector fem and Φepz are element external mechanical force vector and external voltage applied to the piezo layer, respectively, which are defined as follows:(23)fem=∫VeNuTfbdV+∫SeNuTfsdS+NuTfc,For the sensor layer, charge sensing is considered. With zero voltage, from (21), the sensed voltage ϕse is given by [13, 16, 18](24)ϕse=-Kϕϕse-1Kϕuseue.The operation of the amplified control loop implies the actuating voltage is(25)ϕae=Gdϕse+Gvϕ˙se,where Gd and Gv are the feedback control gains for displacement and velocity, respectively.From (21), the charge in the actuator due to actuator strain in response to plate vibration modified by control system feedback is(26)Φae=Kϕueue+KϕϕeGdϕse+Gvϕ˙se=Kϕueue-GdKϕϕeKϕϕse-1Kϕuseue-GvKϕϕeKϕϕse-1Kϕuseu˙e.Substitute (26) into (22) to yield(27)Meu¨e-GvKuϕeKϕϕe-1KϕϕeKϕϕse-1Kϕuseu˙e+Kuulne+Kuunle-GdKuϕeKϕϕe-1KϕϕeKϕϕse-1Kϕuseue=fem,where Kϕuse=Kϕuae=KϕueandKϕϕse=Kϕϕae=Kϕϕe are the mechanical-electrical coupling and pizoelectric permittivity stiffness matrices, respectively.Equation (27) can be rewritten as follows:(28)Meu¨e+CAeu˙e+Kuulne+Kuunle+KAeue=fem,where CAe=-GvKuϕeKϕϕe-1KϕϕeKϕϕse-1Kϕuse is the element active damping matrix and KAe=-GdKuϕeKϕϕe-1KϕϕeKϕϕse-1Kϕuse is the element active stiffness matrix.


### 2.3. Formulation of Stiffener

#### 2.3.1. Formulation of x-Stiffener

(29) U x s x , z = u 0 x + z θ x s x , W x s x , z = w x s x ,where x-axis is taken along the stiffener centerline and the z-axis is its upward normal. The plate and stiffener element are shown in Figure 3.Figure 3
Modeling of plate and stiffener element.If we consider that the x-stiffener is attached to the lower side of the plate, conditions of displacement compatibility along their line of connection can be written as(30)upz=-tp/2=uxsz=txs/2,θxpz=-tp/2=θxsz=txs/2,wpz=-tp/2=wxsz=txs/2,where tp is the plate thickness and txs is the x-stiffener depth.The element stiffness and mass matrices are defined as follows [15]:(31)Kxse=∫leBxsTDxsBxsdx,(32)Mxse=∫AePNu0TNu0+NwTNw+IyNθxTNθxdA,where [Bxs] is the strain-displacement relations matrix, [Dxs] is the stress-strain relations matrix, and le is the element length. Nu0,NwandNθx are the shape function matrices relating the primary variables u0, w, x, in terms of nodal unknowns, Iy is the area moment of inertia related to the y-axis, and P=∑k=1n∫hk-1hkρkdz, with ρk being density of kth layer.


#### 2.3.2. Formulation of y-Stiffener

The same as for x-stiffener, the element stiffness and mass matrices of the y-stiffener are defined as follows:(33)Kyse=∫leBysTDysBysdy,(34)Myse=∫AePNu0TNu0+NwTNw+IxNθyTNθydA,


### 2.4. Modeling the Effect of Aerodynamic Pressure and Motion Equations of the Smart Composite Plate-Stiffeners Element

Based on the first-order theory, the aerodynamic pressurelh and moment mθ can be described as [19, 20](35)lw=12ρaUcos⁡α2BkH1∗w˙Ucos⁡α+kH2∗Bθ˙Ucos⁡α+k2H3∗θ+12CpρaUsin⁡α2,mθ=12ρaUcos⁡α2B2kA1∗w˙Ucos⁡α+kA2∗Bθ˙Ucos⁡α+k2A3∗θ,where k=bω/U is defined as the reduced frequency, ω is the circular frequency of osscillation of the airfoil, U is the mean wind velocity, B is the half-chord length of the airfoil or half-width of the plate, ρa is the air density, and α is the angle of attack.The functionsAi∗(K),Hi∗(K) are defined as follows:(36)H1∗K=-πkFk,H2∗K=-π4k1+Fk+2Gkk,H3∗K=-π2k2Fk-kGk2,A1∗K=π4kFk,A2∗K=-π16k1-Fk-2Gkk,A3∗K=π8k2k28+Fk-kGk2,where F(k) and G(k) are defined as(37)Fk=0,500502k3+0,512607k2+0,2104k+0,021573k3+1,035378k2+0,251293k+0,021508,Gk=-0,000146k3+0,122397k2+0,327214k+0,001995k3+2,481481k2+0,93453k+0,089318.Using finite element method, aerodynamic force vector can be described as(38)feair=-Kaireue-Caireu˙e+fen,where [Kair]e,[Cair]eand {f}en are the aerodynamic stiffness, damping matrices, and lift force vector, respectively:(39)Keair=ρaUcos⁡α2Bk2∫AeH3∗kNwTNθx+BA3∗k∂Nθy∂xTNθxdA,(40)Ceair=ρaUcos⁡αBk∫AeH1∗kNwTNw+BH2∗kNwTNθxdA+∫AeBA1∗k∂Nθy∂xTNw+B2A2∗k∂Nθy∂xTNθxdA,(41)fen=CpρaUsin⁡α2∫AeNwTdA,where Ae is the element area and Nw,Nθ are the shape functions.From (28), (31), (32), (33), (34), and (38), the governing equations of motion of the smart composite plate-stiffeners element subjecteced to aerodynamic without damping can be derived as(42)M∗eu¨e+CAeu˙e+K∗e+KAe+Kaireue=f∗em,with  M∗e=Me+Mxse+Myse,K∗e=Kuulne+Kuunle+Kxse+Kyse,f∗em=fem+fen.


### 2.5. Governing Differential Equations for Total System

Finally, the elemental equations of motion are assembled to obtain the open-loop global equation of motion of the overall stiffened composite plate with the PZT patches as follows:(43)M∗u¨+CR+CAu˙+K∗+KA+Kairue=f∗m,where  CR=αMuu+βKuuln+Kuunl,  CA=-∑NePZTGvKuϕeKϕϕe-1KϕϕeKϕϕse-1Kϕuse,  KA=∑nPZT-GdKuϕeKϕϕe-1KϕϕeKϕϕse-1Kϕuse, α and β are Rayleigh’s coefficients, to account for inherent structural damping.The solution of nonlinear equation (43) is carried out using Newmark direct and Newton-Raphson iteration method [11, 21].


### 2.6. Numerical Applications

A rectangle cantilever laminated composite plate is assumed to be[0°/90°s with total thickness of 4mm, length of 600mm, and width of 400mm with three stiffeners along each of the directions x and y. The geometrical dimension of the stiffener is height of 5mm and width of 10mm. The material properties for plate and stiffeners made of graphite/epoxy are E11 = 181GPa, E22 = E33 = 10.3GPa, E12 = 7.17GPa, ν12 = 0.35, ν23 = ν32 = 0.38, are ρ = 1,600kg·m−3. Material properties for piezoelectric layer made of PZT-5A are d31 = d32 = -171×10−12m/V, d33 = 374×10−12m/V, d15 = d24 = -584×10−12m/V, G12 = 7.17GPa, G23 = 2.87GPa, G32 = 7.17GPa, νPZT = 0.3, ρPZT = 7,600kg·m−3 and thickness tPZT = 0.15876mm, ξ = 0.05, Gv = 0.5, Gd = 15. The effects of the excitation frequency and location of the actuators are presented through a parametric study to investigate the vibration shape of the composite plate activated by the surface bonded piezoelectric actuators. The iterative error of the load εD = 0.02% is chosen.The stiffened plate is subjected to the airflow in the positive x direction as shown in Figure1(a).Figure4 presents the time history response of the plate at critical airflow velocity and α = 0°.Figure 4
History of the plate at a critical airflow velocityUcr = 30.5m/s.Figure7 shows the history of plate dynamics at α = 0° and critical airflow velocity Ucr = 30.5m/s. In this case, at t = 15s, tip displacement (point A in Figure 1(a)) reached 1.75cm. The plate begins to be instable and, at this time, piezoelectric voltage was 125V due to the piezoelectric effect.


## 2.1. Piezoelectric Equations

Consider a laminated composite plate with integrated sensors and actuators as shown in Figure1. It is assumed that each layer of the plate possesses a plane of elastic symmetry parallel to the x-y plane; the kth layer’s lamina constitutive equations coupling the direct and converse piezoelectric equations can be expressed as [13, 14](1)σ1σ2τ12k=Q-11Q-120Q-12Q-22000Q-66kε1ε2γ12k-00e-3100e-32000kE1E2E3k,(2)D1D2D3k=000000e-31e-320kε1ε2γ12k+p-11000p-22000p-33kE1E2E3k,where Q-ij,e-ij,p-ij are, respectively, the plane-stress reduced elastic constants, the piezoelectric constants, and the permittivity coefficients of the kth lamina in its material coordinate system; σk=σ1σ2τ12kT,  εk=ε1ε2γ12kT,  Ek=E1E2E3kT,  Dk=D1D2D3kT are the stress, strain, electric field, and electric displacement matric components, respectively, to the material coordinate system. If the voltage is applied to the actuator in the thickness only, then Ek=00-1/hkTϕk=Bϕϕk, Vk is the applied voltage across the kth ply, and hk is the thickness of the kth piezoelectric layer. The plane-stress reduced elastic constants Q-ij are given as(3)Q-11=E11-ν12ν21,Q-12=ν12E21-ν12ν21,Q-22=E21-ν12ν21,Q-66=G12,Upon transformation, the lamina piezoelectric equations can be expressed in terms of the stress, strains, and electric displacements in the plate coordinates as [13, 14](4)σxσyτxyk=Q11Q12Q16Q12Q22Q26Q16Q26Q66kεxεyγxyk-00e3100e3200e36kExEyEzk,(5)DxDyDzk=000000e31e32e36kεxεyγxyk+p11p120p12p22000p33kExEyEzk.Equations (4) and (5) can also be written as(6)σk=Qkεk-ekTEk,(7)Dk=ekεk+pkEk,TheQij are the components of the transformed lamina stiffness matrix which are defined as follows:(8)Q11=Q-11cos4θ+Q-22sin4θ+2Q-12+2Q-66sin2θcos2θ,Q12=Q-11+Q-22-4Q-66sin2θcos2θ+Q-12sin4θ+cos4θ,Q22=Q-11sin4θ+Q-22cos4θ+2Q-12+2Q-66sin2θcos2θ,Q16=Q-11-Q-12-2Q-66cos3θsin⁡θ-Q-22-Q-12-2Q-66sin3θcos⁡θ,Q26=Q-11-Q-12-2Q-66sin3θcos⁡θ-Q-22-Q-12-2Q-66sin⁡θcos3θ,Q66=Q-11+Q-12-2Q-12-2Q-66sin2θcos2θ+Q-66sin4θ+cos4θ,where θ is the lamina orientation angle.The piezoelectric constants matrix [e] is unavailable and it can be expressed by the more commonly available piezoelectric strain constant matrix which is [d] as [13](9)e=dQ,where d=000d15000d2400d31d32000.


## 2.2. Smart Piezoelectric Composite Plate Element Formulation

A smart piezoelectric composite plate element is considered with coordinates x, y along the in-plane direction and z along the thickness direction as shown in Figure2.Figure 2
Smart piezoelectric composite plate element.Using the Mindlin formulation, the displacements u, v, and w at a point (x,y,z) form the median surface and are expressed as functions of mid-plane displacements u0, v0, and w0 and independent rotations θx and θy of the normal in xz and yz planes, respectively, as [11, 15, 16](10)ux,y,z,t=u0x,y,t+zθyx,y,t,vx,y,z,t=v0x,y,t-zθxx,y,t,wx,y,z,t=w0x,y,t,The components of the strain vector corresponding to the displacement field (10) are [7, 11](11)εxεyγxy=εbln+εnl=εb,where {εbln},{εbnl} are the linear and nonlinear strain vector, respectively, and defined as follows:(12)εbln=εxεyγxy=εxoεyoγxyo+zκxκyκxy=∂∂x00∂∂y∂∂y∂∂xu0v0+z-∂∂y00-∂∂x-∂∂y∂∂xθxθy=ε0+zκ,εbnl=12∂w0∂x00∂w0∂y∂w0∂y∂w0∂x∂∂x∂∂yw0.In the finite formulation, the displacement fieldu and the electric potential ϕ over an element are related to the corresponding node values ue of the element and the electric charges of the piezo layer ϕe by the mean of the shape functions [Nu], [Nϕ], as follows [15]:(13)u=Nuqe,(14)ϕ=Nϕϕe.Equation (11) can be expressed as(15)ε=Bulnue+Bunlue.where Buln and Bunl are linear and nonlinear strain matric, respectively, Nu are the shape functions, and Nϕ are the shape functions for the electric potential.The dynamic equations of a finite smart laminated composite plate can be derived by using Hamilton’s principle [15, 17]:(16)δ∫t1t2Te-Ue+Wedt=0,where Te is the kinetic energy, Ue is the strain energy, and We is the work done by the applied forces.The kinetic energy at the element level is defined as(17)Te=12∫Veρu˙Tu˙dV,where Ve is the volume of the plate element.The strain energy can be written as(18)Ue=12∫VeεTσdV,The work done by the external forces(19)We=∫VeuTfbdV+∫SeuTfsdS+uTfc,where {fb} is the body force, Se is the surface area of the plate element, {fs} is the surface force, and {fc} is the concentrated load.Substituting (17), (18), and (19) into (16), noting that the electric field vector E=-Bϕϕe, and using (6), (7), (13), (14), and (15), the dynamic matrix equations can be written as(20)Meu¨e+Kuulne+Kuunleue+Kuϕeϕe=fem,(21)Kϕueue+Kϕϕeϕe=Φepz,Substitute (21) into (20) to obtain(22)Meu¨e+Kuulne+Kuunle-KuϕeKϕϕe-1Kϕueue=fem-KuϕeKϕϕe-1Φepz,where Me=∫VeρNuTNudV, Kuulne=∫VeBulnTQBulndV, Kuunle=∫VeBuunlTQBuunldV, Kuϕe=∫VeBuTeTBϕdV, Kϕϕe=∫VeBϕTpBϕdV, Kϕue=KuϕeT, which are, respectively, the element mass, linear mechanical stiffness, nonlinear mechanical stiffness, mechanical-electrical coupling, piezoelectric permittivity, and electrical-mechanical coupling matrix. The vector fem and Φepz are element external mechanical force vector and external voltage applied to the piezo layer, respectively, which are defined as follows:(23)fem=∫VeNuTfbdV+∫SeNuTfsdS+NuTfc,For the sensor layer, charge sensing is considered. With zero voltage, from (21), the sensed voltage ϕse is given by [13, 16, 18](24)ϕse=-Kϕϕse-1Kϕuseue.The operation of the amplified control loop implies the actuating voltage is(25)ϕae=Gdϕse+Gvϕ˙se,where Gd and Gv are the feedback control gains for displacement and velocity, respectively.From (21), the charge in the actuator due to actuator strain in response to plate vibration modified by control system feedback is(26)Φae=Kϕueue+KϕϕeGdϕse+Gvϕ˙se=Kϕueue-GdKϕϕeKϕϕse-1Kϕuseue-GvKϕϕeKϕϕse-1Kϕuseu˙e.Substitute (26) into (22) to yield(27)Meu¨e-GvKuϕeKϕϕe-1KϕϕeKϕϕse-1Kϕuseu˙e+Kuulne+Kuunle-GdKuϕeKϕϕe-1KϕϕeKϕϕse-1Kϕuseue=fem,where Kϕuse=Kϕuae=KϕueandKϕϕse=Kϕϕae=Kϕϕe are the mechanical-electrical coupling and pizoelectric permittivity stiffness matrices, respectively.Equation (27) can be rewritten as follows:(28)Meu¨e+CAeu˙e+Kuulne+Kuunle+KAeue=fem,where CAe=-GvKuϕeKϕϕe-1KϕϕeKϕϕse-1Kϕuse is the element active damping matrix and KAe=-GdKuϕeKϕϕe-1KϕϕeKϕϕse-1Kϕuse is the element active stiffness matrix.


## 2.3. Formulation of Stiffener

### 2.3.1. Formulation of x-Stiffener

(29) U x s x , z = u 0 x + z θ x s x , W x s x , z = w x s x ,where x-axis is taken along the stiffener centerline and the z-axis is its upward normal. The plate and stiffener element are shown in Figure 3.Figure 3
Modeling of plate and stiffener element.If we consider that the x-stiffener is attached to the lower side of the plate, conditions of displacement compatibility along their line of connection can be written as(30)upz=-tp/2=uxsz=txs/2,θxpz=-tp/2=θxsz=txs/2,wpz=-tp/2=wxsz=txs/2,where tp is the plate thickness and txs is the x-stiffener depth.The element stiffness and mass matrices are defined as follows [15]:(31)Kxse=∫leBxsTDxsBxsdx,(32)Mxse=∫AePNu0TNu0+NwTNw+IyNθxTNθxdA,where [Bxs] is the strain-displacement relations matrix, [Dxs] is the stress-strain relations matrix, and le is the element length. Nu0,NwandNθx are the shape function matrices relating the primary variables u0, w, x, in terms of nodal unknowns, Iy is the area moment of inertia related to the y-axis, and P=∑k=1n∫hk-1hkρkdz, with ρk being density of kth layer.


### 2.3.2. Formulation of y-Stiffener

The same as for x-stiffener, the element stiffness and mass matrices of the y-stiffener are defined as follows:(33)Kyse=∫leBysTDysBysdy,(34)Myse=∫AePNu0TNu0+NwTNw+IxNθyTNθydA,


## 2.3.1. Formulation of x-Stiffener

(29) U x s x , z = u 0 x + z θ x s x , W x s x , z = w x s x ,where x-axis is taken along the stiffener centerline and the z-axis is its upward normal. The plate and stiffener element are shown in Figure 3.Figure 3
Modeling of plate and stiffener element.If we consider that the x-stiffener is attached to the lower side of the plate, conditions of displacement compatibility along their line of connection can be written as(30)upz=-tp/2=uxsz=txs/2,θxpz=-tp/2=θxsz=txs/2,wpz=-tp/2=wxsz=txs/2,where tp is the plate thickness and txs is the x-stiffener depth.The element stiffness and mass matrices are defined as follows [15]:(31)Kxse=∫leBxsTDxsBxsdx,(32)Mxse=∫AePNu0TNu0+NwTNw+IyNθxTNθxdA,where [Bxs] is the strain-displacement relations matrix, [Dxs] is the stress-strain relations matrix, and le is the element length. Nu0,NwandNθx are the shape function matrices relating the primary variables u0, w, x, in terms of nodal unknowns, Iy is the area moment of inertia related to the y-axis, and P=∑k=1n∫hk-1hkρkdz, with ρk being density of kth layer.


## 2.3.2. Formulation of y-Stiffener

The same as for x-stiffener, the element stiffness and mass matrices of the y-stiffener are defined as follows:(33)Kyse=∫leBysTDysBysdy,(34)Myse=∫AePNu0TNu0+NwTNw+IxNθyTNθydA,


## 2.4. Modeling the Effect of Aerodynamic Pressure and Motion Equations of the Smart Composite Plate-Stiffeners Element

Based on the first-order theory, the aerodynamic pressurelh and moment mθ can be described as [19, 20](35)lw=12ρaUcos⁡α2BkH1∗w˙Ucos⁡α+kH2∗Bθ˙Ucos⁡α+k2H3∗θ+12CpρaUsin⁡α2,mθ=12ρaUcos⁡α2B2kA1∗w˙Ucos⁡α+kA2∗Bθ˙Ucos⁡α+k2A3∗θ,where k=bω/U is defined as the reduced frequency, ω is the circular frequency of osscillation of the airfoil, U is the mean wind velocity, B is the half-chord length of the airfoil or half-width of the plate, ρa is the air density, and α is the angle of attack.The functionsAi∗(K),Hi∗(K) are defined as follows:(36)H1∗K=-πkFk,H2∗K=-π4k1+Fk+2Gkk,H3∗K=-π2k2Fk-kGk2,A1∗K=π4kFk,A2∗K=-π16k1-Fk-2Gkk,A3∗K=π8k2k28+Fk-kGk2,where F(k) and G(k) are defined as(37)Fk=0,500502k3+0,512607k2+0,2104k+0,021573k3+1,035378k2+0,251293k+0,021508,Gk=-0,000146k3+0,122397k2+0,327214k+0,001995k3+2,481481k2+0,93453k+0,089318.Using finite element method, aerodynamic force vector can be described as(38)feair=-Kaireue-Caireu˙e+fen,where [Kair]e,[Cair]eand {f}en are the aerodynamic stiffness, damping matrices, and lift force vector, respectively:(39)Keair=ρaUcos⁡α2Bk2∫AeH3∗kNwTNθx+BA3∗k∂Nθy∂xTNθxdA,(40)Ceair=ρaUcos⁡αBk∫AeH1∗kNwTNw+BH2∗kNwTNθxdA+∫AeBA1∗k∂Nθy∂xTNw+B2A2∗k∂Nθy∂xTNθxdA,(41)fen=CpρaUsin⁡α2∫AeNwTdA,where Ae is the element area and Nw,Nθ are the shape functions.From (28), (31), (32), (33), (34), and (38), the governing equations of motion of the smart composite plate-stiffeners element subjecteced to aerodynamic without damping can be derived as(42)M∗eu¨e+CAeu˙e+K∗e+KAe+Kaireue=f∗em,with  M∗e=Me+Mxse+Myse,K∗e=Kuulne+Kuunle+Kxse+Kyse,f∗em=fem+fen.


## 2.5. Governing Differential Equations for Total System

Finally, the elemental equations of motion are assembled to obtain the open-loop global equation of motion of the overall stiffened composite plate with the PZT patches as follows:(43)M∗u¨+CR+CAu˙+K∗+KA+Kairue=f∗m,where  CR=αMuu+βKuuln+Kuunl,  CA=-∑NePZTGvKuϕeKϕϕe-1KϕϕeKϕϕse-1Kϕuse,  KA=∑nPZT-GdKuϕeKϕϕe-1KϕϕeKϕϕse-1Kϕuse, α and β are Rayleigh’s coefficients, to account for inherent structural damping.The solution of nonlinear equation (43) is carried out using Newmark direct and Newton-Raphson iteration method [11, 21].


## 2.6. Numerical Applications

A rectangle cantilever laminated composite plate is assumed to be[0°/90°s with total thickness of 4mm, length of 600mm, and width of 400mm with three stiffeners along each of the directions x and y. The geometrical dimension of the stiffener is height of 5mm and width of 10mm. The material properties for plate and stiffeners made of graphite/epoxy are E11 = 181GPa, E22 = E33 = 10.3GPa, E12 = 7.17GPa, ν12 = 0.35, ν23 = ν32 = 0.38, are ρ = 1,600kg·m−3. Material properties for piezoelectric layer made of PZT-5A are d31 = d32 = -171×10−12m/V, d33 = 374×10−12m/V, d15 = d24 = -584×10−12m/V, G12 = 7.17GPa, G23 = 2.87GPa, G32 = 7.17GPa, νPZT = 0.3, ρPZT = 7,600kg·m−3 and thickness tPZT = 0.15876mm, ξ = 0.05, Gv = 0.5, Gd = 15. The effects of the excitation frequency and location of the actuators are presented through a parametric study to investigate the vibration shape of the composite plate activated by the surface bonded piezoelectric actuators. The iterative error of the load εD = 0.02% is chosen.The stiffened plate is subjected to the airflow in the positive x direction as shown in Figure1(a).Figure4 presents the time history response of the plate at critical airflow velocity and α = 0°.Figure 4
History of the plate at a critical airflow velocityUcr = 30.5m/s.Figure7 shows the history of plate dynamics at α = 0° and critical airflow velocity Ucr = 30.5m/s. In this case, at t = 15s, tip displacement (point A in Figure 1(a)) reached 1.75cm. The plate begins to be instable and, at this time, piezoelectric voltage was 125V due to the piezoelectric effect.


## 3. Experimental Validation

### 3.1. Experimental Model

A rectangle cantilever laminated composite plate is assumed to be45°/0°/-45°s with total thickness of 4mm, length of 500mm, and width of 300mm with three stiffeners along each of the directions x and y (45°/0°/-45°s). The geometrical dimension of the stiffener is height of 8mm and width of 10mm (Figure 5). The piezoelectric slice is 200 mm by 90mm dimension and placed in the location as in Figure 6. The plate and stiffeners made of T300/2500 graphite/epoxy are produced by Torayca Co., Japan. The mean values and coefficients of variation of the experimentally determined material constants are given as E1 = 124.68GPa, E2 = 9.6GPa, G12 = 8.64GPa, G23 = 2.32GPa, and ν12 = 0.33.Figure 5
The preparation of experimental model.Figure 6
The real PZT model.Figure 7
The diagram of experiment.


### 3.2. The Diagram of Experiment

The diagram of an experiment is shown in Figure7. Alternating current from the electric grid running through transmitter is linearly amplified to generate voltage V with frequency f imposing on PZT slice. The PZT slices attaching on the shell are playing a role to arouse the vibration of the whole system. Through two acceleration sensors attached to the shell, the dynamic response of the shell at the location of the measurement points will be captured, displaying on the oscilloscope display and then storing in the computer.The main experiments are shown herein:(i)
Measurement of the acceleration response of the structure at the measuring points arranged on the top of the plate corresponding to different voltage and frequency stimulation levels placed on two piezoelectric plates.(ii)
Measurement of the first free vibration of the structure.


### 3.3. Measuring Devices

#### 3.3.1. The Signal Generator

The signal generator has a function to generate randomly the sinusoidal alternating current with voltage V and frequency f placing on 2 piezoelectric plates attached to the shell to establish vibration of the structure. In this experiment, the signal generator is used by combining the speaker test oscillator and a linear amplifier as in Figure8. The speaker test oscillator is used in this experiment is Onsoku branch (Japan) model OG-422A. The specifications of the machine are as follows: the maximum output voltage: ± 200 V, the response frequency range (output): 10÷20 kHz, capacity compliant amplifier: 20 (W), compatible power supply: 100V, 115V, 220V and 240 V, frequency 50/60Hz.Figure 8
The signal generator.
(a)
The speaker test oscillator (b)
The linear amplifier


#### 3.3.2. Generating Load Equipment

The airflow load is created through an open wind tunnel with capacity 11kW, the area of test cross section is 1000mm×1000mm, and airflow speed can be changed from 0 to 40 m/s (Figure9). The experimental model is placed on a bracket with a diameter ϕ500 mm; this bracket can rotate with different angles compared to the wind direction.Figure 9
The open wind tunnel and the bracket.
(a) (b)


#### 3.3.3. The Oscilloscopes

The oscilloscope has a function which is capturing the acceleration response of the plate at the point of measurement, then displaying on screen and storing the data in the computer. The devices to measure the vibrations used in this experiment include accelerometer sensor (2 sensors), two piezoelectric data receivers HnB75B, oscilloscope display Tektronix TDS-1012, and a computer. The accelerometer sensors used in this experiment are type ACH01-02 (Figure9(a)). The specifications of these accelerometer sensors include parameters: the sensitivity 10mV/g, frequency band 1,0Hz ÷ 20kHz, dynamic range: ± 250g, resolution 40μgHz, resonant frequency: > 35kHz, linearity: 0.1%, and the maximum pick: 1000g (g: gravity acceleration).Two piezoelectric channels’ data acquisition HnB75B has a function to amplify the electrical signals from the sensor to the oscilloscope and the computer (Figure9(b)). The specifications of this device include: number of channels - 2 simultaneous channels, voltage range measurement: 0 ÷ 5V, accuracy: +/- 0.2%, sampling frequency: maximum 10 kHz, connection port: 1 com port, 2 AC ports.Through the simulation program set up already in devices and computer, the vibrations of the plate at measurement points are displayed on the screen of Tektronix TDS-1012 machine with the vibration parameters at each oscillation cycle (Figure10(c)).Figure 10
The oscilloscopes.
(a)
The acceleration sensor (b)
HnB75B data receiver (c)
The oscilloscope display


### 3.4. The Results of Measuring Acceleration Response Experiment

In this experiment, the authors conduct three cases:α = 0°, α = 22.5°, and α = 45°.For α = 0°: The speed of the airflow U = 5m/s and 10 m/s; the voltage Vin = 9.30V, and Kamp = 20. Hence the voltage level will be Vin × Kamp = 9.30 V × 20 = 186.0V; sampling frequency fsamp = 1000 Hz; frequency of stimulationf = 6.944 Hz.For α = 22.5°: The speed of the airflow U = 5m/s and 10m/s; the voltage Vin = 9.50V, and Kamp = 20. Hence the voltage level will be Vin × Kamp = 9.50 V × 20 = 190.0V; sampling frequency fsamp = 1000 Hz; frequency of stimulationf = 6.981 Hz.For α = 45°: The speed of the airflow U = 5m/s and 10m/s; the voltage Vin = 9.20V, and Kamp = 20. Hence the voltage level will be Vin × Kamp = 9.20 V × 20 = 184.0 V; sampling frequency fsamp = 1000 Hz; frequency of stimulation f = 6.993 Hz.The acceleration-time and amplitude-frequency response at the measuring point of the above stiffened composite plate with PZT patches are shown in Figures11 and 12.Figure 11
The acceleration response in one measurement withα = 0°, U = 5m/s.Figure 12
The acceleration response and frequency spectrum with f = 13.672Hz, U = 10m/s,α = 45°.We calculate the acceleration response for the above stiffened composite plate with PZT patches by our computer program. The acceleration response at the measuring point of the plate will be compared with those of experimental ones and is given in Figure13 and Table 1.Table 1
The maximum acceleration at U = 5m/s.
Method a m a x m / s 2 α [ d e g ] Error [%] Experiment 0.050 45 11.8 FEM 0.0559  Experiment 0.061 22.5 10.7 FEM 0.0675  Experiment 0.072 0 10.1 FEM 0.0793 Comment: The theoretical and experimental results for the acceleration response of the plate at the point of measurement over time are fairly uniform, with maximum error of 11.80% for acceleration being acceptable. Thus, the results of the experiment showed the relevance of algorithms and programs in published theoretical research.Figure 13
The acceleration response of plate at the measuring position (U = 5m/s).


## 3.1. Experimental Model

A rectangle cantilever laminated composite plate is assumed to be45°/0°/-45°s with total thickness of 4mm, length of 500mm, and width of 300mm with three stiffeners along each of the directions x and y (45°/0°/-45°s). The geometrical dimension of the stiffener is height of 8mm and width of 10mm (Figure 5). The piezoelectric slice is 200 mm by 90mm dimension and placed in the location as in Figure 6. The plate and stiffeners made of T300/2500 graphite/epoxy are produced by Torayca Co., Japan. The mean values and coefficients of variation of the experimentally determined material constants are given as E1 = 124.68GPa, E2 = 9.6GPa, G12 = 8.64GPa, G23 = 2.32GPa, and ν12 = 0.33.Figure 5
The preparation of experimental model.Figure 6
The real PZT model.Figure 7
The diagram of experiment.


## 3.2. The Diagram of Experiment

The diagram of an experiment is shown in Figure7. Alternating current from the electric grid running through transmitter is linearly amplified to generate voltage V with frequency f imposing on PZT slice. The PZT slices attaching on the shell are playing a role to arouse the vibration of the whole system. Through two acceleration sensors attached to the shell, the dynamic response of the shell at the location of the measurement points will be captured, displaying on the oscilloscope display and then storing in the computer.The main experiments are shown herein:(i)
Measurement of the acceleration response of the structure at the measuring points arranged on the top of the plate corresponding to different voltage and frequency stimulation levels placed on two piezoelectric plates.(ii)
Measurement of the first free vibration of the structure.


## 3.3. Measuring Devices

### 3.3.1. The Signal Generator

The signal generator has a function to generate randomly the sinusoidal alternating current with voltage V and frequency f placing on 2 piezoelectric plates attached to the shell to establish vibration of the structure. In this experiment, the signal generator is used by combining the speaker test oscillator and a linear amplifier as in Figure8. The speaker test oscillator is used in this experiment is Onsoku branch (Japan) model OG-422A. The specifications of the machine are as follows: the maximum output voltage: ± 200 V, the response frequency range (output): 10÷20 kHz, capacity compliant amplifier: 20 (W), compatible power supply: 100V, 115V, 220V and 240 V, frequency 50/60Hz.Figure 8
The signal generator.
(a)
The speaker test oscillator (b)
The linear amplifier


### 3.3.2. Generating Load Equipment

The airflow load is created through an open wind tunnel with capacity 11kW, the area of test cross section is 1000mm×1000mm, and airflow speed can be changed from 0 to 40 m/s (Figure9). The experimental model is placed on a bracket with a diameter ϕ500 mm; this bracket can rotate with different angles compared to the wind direction.Figure 9
The open wind tunnel and the bracket.
(a) (b)


### 3.3.3. The Oscilloscopes

The oscilloscope has a function which is capturing the acceleration response of the plate at the point of measurement, then displaying on screen and storing the data in the computer. The devices to measure the vibrations used in this experiment include accelerometer sensor (2 sensors), two piezoelectric data receivers HnB75B, oscilloscope display Tektronix TDS-1012, and a computer. The accelerometer sensors used in this experiment are type ACH01-02 (Figure9(a)). The specifications of these accelerometer sensors include parameters: the sensitivity 10mV/g, frequency band 1,0Hz ÷ 20kHz, dynamic range: ± 250g, resolution 40μgHz, resonant frequency: > 35kHz, linearity: 0.1%, and the maximum pick: 1000g (g: gravity acceleration).Two piezoelectric channels’ data acquisition HnB75B has a function to amplify the electrical signals from the sensor to the oscilloscope and the computer (Figure9(b)). The specifications of this device include: number of channels - 2 simultaneous channels, voltage range measurement: 0 ÷ 5V, accuracy: +/- 0.2%, sampling frequency: maximum 10 kHz, connection port: 1 com port, 2 AC ports.Through the simulation program set up already in devices and computer, the vibrations of the plate at measurement points are displayed on the screen of Tektronix TDS-1012 machine with the vibration parameters at each oscillation cycle (Figure10(c)).Figure 10
The oscilloscopes.
(a)
The acceleration sensor (b)
HnB75B data receiver (c)
The oscilloscope display


## 3.3.1. The Signal Generator

The signal generator has a function to generate randomly the sinusoidal alternating current with voltage V and frequency f placing on 2 piezoelectric plates attached to the shell to establish vibration of the structure. In this experiment, the signal generator is used by combining the speaker test oscillator and a linear amplifier as in Figure8. The speaker test oscillator is used in this experiment is Onsoku branch (Japan) model OG-422A. The specifications of the machine are as follows: the maximum output voltage: ± 200 V, the response frequency range (output): 10÷20 kHz, capacity compliant amplifier: 20 (W), compatible power supply: 100V, 115V, 220V and 240 V, frequency 50/60Hz.Figure 8
The signal generator.
(a)
The speaker test oscillator (b)
The linear amplifier


## 3.3.2. Generating Load Equipment

The airflow load is created through an open wind tunnel with capacity 11kW, the area of test cross section is 1000mm×1000mm, and airflow speed can be changed from 0 to 40 m/s (Figure9). The experimental model is placed on a bracket with a diameter ϕ500 mm; this bracket can rotate with different angles compared to the wind direction.Figure 9
The open wind tunnel and the bracket.
(a) (b)


## 3.3.3. The Oscilloscopes

The oscilloscope has a function which is capturing the acceleration response of the plate at the point of measurement, then displaying on screen and storing the data in the computer. The devices to measure the vibrations used in this experiment include accelerometer sensor (2 sensors), two piezoelectric data receivers HnB75B, oscilloscope display Tektronix TDS-1012, and a computer. The accelerometer sensors used in this experiment are type ACH01-02 (Figure9(a)). The specifications of these accelerometer sensors include parameters: the sensitivity 10mV/g, frequency band 1,0Hz ÷ 20kHz, dynamic range: ± 250g, resolution 40μgHz, resonant frequency: > 35kHz, linearity: 0.1%, and the maximum pick: 1000g (g: gravity acceleration).Two piezoelectric channels’ data acquisition HnB75B has a function to amplify the electrical signals from the sensor to the oscilloscope and the computer (Figure9(b)). The specifications of this device include: number of channels - 2 simultaneous channels, voltage range measurement: 0 ÷ 5V, accuracy: +/- 0.2%, sampling frequency: maximum 10 kHz, connection port: 1 com port, 2 AC ports.Through the simulation program set up already in devices and computer, the vibrations of the plate at measurement points are displayed on the screen of Tektronix TDS-1012 machine with the vibration parameters at each oscillation cycle (Figure10(c)).Figure 10
The oscilloscopes.
(a)
The acceleration sensor (b)
HnB75B data receiver (c)
The oscilloscope display


## 3.4. The Results of Measuring Acceleration Response Experiment

In this experiment, the authors conduct three cases:α = 0°, α = 22.5°, and α = 45°.For α = 0°: The speed of the airflow U = 5m/s and 10 m/s; the voltage Vin = 9.30V, and Kamp = 20. Hence the voltage level will be Vin × Kamp = 9.30 V × 20 = 186.0V; sampling frequency fsamp = 1000 Hz; frequency of stimulationf = 6.944 Hz.For α = 22.5°: The speed of the airflow U = 5m/s and 10m/s; the voltage Vin = 9.50V, and Kamp = 20. Hence the voltage level will be Vin × Kamp = 9.50 V × 20 = 190.0V; sampling frequency fsamp = 1000 Hz; frequency of stimulationf = 6.981 Hz.For α = 45°: The speed of the airflow U = 5m/s and 10m/s; the voltage Vin = 9.20V, and Kamp = 20. Hence the voltage level will be Vin × Kamp = 9.20 V × 20 = 184.0 V; sampling frequency fsamp = 1000 Hz; frequency of stimulation f = 6.993 Hz.The acceleration-time and amplitude-frequency response at the measuring point of the above stiffened composite plate with PZT patches are shown in Figures11 and 12.Figure 11
The acceleration response in one measurement withα = 0°, U = 5m/s.Figure 12
The acceleration response and frequency spectrum with f = 13.672Hz, U = 10m/s,α = 45°.We calculate the acceleration response for the above stiffened composite plate with PZT patches by our computer program. The acceleration response at the measuring point of the plate will be compared with those of experimental ones and is given in Figure13 and Table 1.Table 1
The maximum acceleration at U = 5m/s.
Method a m a x m / s 2 α [ d e g ] Error [%] Experiment 0.050 45 11.8 FEM 0.0559  Experiment 0.061 22.5 10.7 FEM 0.0675  Experiment 0.072 0 10.1 FEM 0.0793 Comment: The theoretical and experimental results for the acceleration response of the plate at the point of measurement over time are fairly uniform, with maximum error of 11.80% for acceleration being acceptable. Thus, the results of the experiment showed the relevance of algorithms and programs in published theoretical research.Figure 13
The acceleration response of plate at the measuring position (U = 5m/s).


## 4. Conclusion

The nonlinear dynamics analysis of the piezoelectric stiffened composite plate subjected to airflow using the finite element and experimental method has been presented. In this paper, we have presented a nine-noded stiffened rectangular composite plate element with PZT patches for the nonlinear vibration analysis of the piezoelectric stiffened composite plates subjected to airflow. The critical velocity of the airflow is determined by numerical calculations. The finite element results compare well with experimental ones. It is recommended that the present formulation can be used to determine the characteristics of the vibration and stability in the analysis and design of the piezoelectric stiffened composite plate structures subjected to airflow applied to the flying instruments.



---

*Source: 2697242-2019-01-29.xml*