# Application of Energy Finite Element Method in Active Vibration Control of Piezoelectric Intelligent Beam

Authors: Jinhua Xie; Rui Huo; Yanfeng Guan; Zhen Zhou
Journal: Advances in Acoustics and Vibration
 (2012)

---

## Abstract

Based on the transmission and equilibrium relationship of vibration energy in beam-like structures, the Galerkin weighted residual method was applied to equation discretization. An equivalent transformation of feedback element was suggested to develop the Energy Finite Element model of a composite piezoelectric cantilever beam driven by harmonic excitation on lateral direction, with both systems with and without time delay being studied and the power input estimation of harmonic excitation was discussed for the resolution of Energy Finite Element function. Then the energy density solutions of the piezoelectric coupling beam through Energy Finite Element Method (EFEM) and classical wave theory were compared to verify the EFEM model, which presented a good accordance. Further investigation was undertaken about the influence of control parameters including the feedback gain and  arrangement of piezoelectric patches on characteristics of system energy density distribution.


---

## Body

## 1. Introduction

The birth of intelligent structure system is due to the mutual penetration of information science and engineering and materials science. Piezoelectric intelligent structure, as one of intelligent structures, is interpreted as light damping system attached with piezoelectric materials on the surface. Based on various advantages of piezoelectric materials, such as easy measurement and control of input and output of the signals, wide range of frequency response, both use of sensor and actuator and so on, it has been largely applied in aerospace, robot, information technology, IT, materials science, and other high technologies. Many scholars have made deep researches and gained rich achievements, and as a result, active vibration control systems using piezoelectric materials as distributed sensors and actuators (S/As) are becoming more and more popular in recent years.As the object of study, piezoelectric intelligent structure’s dynamic modeling and optimal design of structure parameters are the problems we concerned about most. External perturbation is always treated as pulse excitation to discuss the suppression problems of system vibration. Obtaining the system’s state-space equation through modal theory, the optimal locations of piezoelectric sensors and actuators were searched taking the smallest vibration energy as the objective function [1]. In  [2], governing equation was established based on steady stimulation and then made an experimental research on frequency response characteristic. In   [3], the boring bar vibration was controlled by using piezoelectric materials; simulation results suggest that the control effect was well.The application of Finite Element Method (FEM) will be involved when it needs to solve the analysis and design problems of complex engineering structures. In recent years, the FEM analysis for vibration control of piezoelectric intelligent structures has become some research hotspot. A piezoelectric smart structure model had been discussed in   [4], using the finite element theory.  The finite element dynamics equation of piezoelectric intelligent structure was deduced in   [5]. However, the application of FEM is generally confined in relatively lower frequency domain, which is resulted from a high sensitivity of the numerical calculation to tiny alteration of structure parameters when the number of elements in the FEM model becomes very great. Although Statistical Energy Analysis (SEA) has become a mature method in energy flow analysis, and related commercial software such as AutoSEA has been in use, it would face the problem of insufficiency in modal density in relatively lower frequency band. In addition, there is not mature theory of structure in loss factor and coupling loss factor [6]. There exists some middle frequency band in which either FEM or SEA would not be credible. Compared with SEA, as an energy analysis method based on wave theory, EFEM does not have any request for structural modal density, and it is more rigorous in describing the energy of system balance. Thus, the EFEM would overcome both the disadvantages of FEM and SEA in solving middle frequency problems and will have a great development in resolving the structure vibration and acoustic radiation problems [7].In this paper, based on the transmission and equilibrium relationship of vibration energy in beam-like structures, the Galerkin weighted residual method was applied to equation discretization. An equivalent transformation of feedback element was suggested to develop the Energy Finite Element model of a composite piezoelectric cantilever beam driven by harmonic excitation on lateral direction, both systems with and without time delay being studied. And the power input estimation of harmonic excitation was discussed for the resolution of Energy Finite Element function. Then the energy density solutions of the piezoelectric coupling beam through Energy Finite Element Method (EFEM) and classical wave theory were compared to verify the EFEM model, which presented a good accordance. Further investigation was undertaken about the influence of control parameters including the feedback gain and arrangement of piezoelectric patches on characteristics of system energy density distribution.


## 2. Energy Finite Element Model of Piezoelectric Cantilever System

### 2.1. System Model

Composite piezoelectric cantilever beam is one of the most typical piezoelectric intelligent structures, which is presented in Figure1(a). In the demonstrated system, series of piezoelectric material patches are grouped in couples and attached symmetrically to a host beam, one patch of each couple being used as an actuator and the other a sensor. The host beam is excited by a harmonically varying force with magnitude F0 which would cause a transverse bending vibration. hp and hb represent the thickness of piezoelectric patch and cantilever beam, respectively.Piezoelectric cantilever system model and energy flow balance.
(a)
(b)Figure1(b) shows the vibration energy inflow and outflow of cubic element, in which q and πdiss represent power flow through cubic element and the energy dissipated by damping, respectively, and Qin is the energy input from external. At each point of the structure, the instantaneous energy balance equation is given by  [8]
(1)∂e∂t=-∂q∂x-πdiss+Qin,
where ∂e/∂t is the variation ratio of energy density to time.The time-averaged value ofq, πdiss, and e, that is, 〈q〉, 〈πdiss〉, and 〈e〉, could be related as
(2)〈q〉=(cg2ηω)·d〈e〉dx,〈πdiss〉=ηω〈e〉.By substituting expression (2) into  (1), the energy balance equation is written as
(3)(cg2ηω)·(d2〈e〉dx2)+ηω〈e〉=〈Qin〉,
where cg=2(ω2EI/ρA) is the group velocity of bending wave and η is the damping loss factor.


### 2.2. Energy Density Control Equations without Feedback

When the feedback does not exist, the harmonic excitation is the only energy source in the system. Equation (3) could be discretized by Galerkin weighted residual method:
(4)K(n)e(n)=P(n)+Q(n),(5)K(n)=∫xL(n)xR(n)(-cg2ωηdNTdxdNdx+ωηNTN)dx,(6)P(n)=N(xR(n))〈q〉|xR(n)-N(xL(n))〈q〉|xL(n),Q(n)=∫xixjN〈Qin〉dx,
where  the superscript (n) denotes the nth discretizednodes, K(n) is a coefficient matrix related to stiffness and mass, P(n) is the power flow at corresponding node, Q(n) is the external energy input, and N(x) is the lagrange linear interpolation basic function.For solving expression (4), coupling relationship among elements needs to be solved. Suppose the numbers of two adjacent elements are i and i+1 and  τ(i,i+1) and γ(i,i+1) represent transmission and reflection coefficient between coupling nodes:
(7)〈q-〉R(i)=γ(i,i+1)〈q+〉R(i)+τ(i,i+1)〈q-〉L(i+1),〈q+〉L(i+1)=τ(i,,i+1)〈q+〉R(i)+γ(i,i+1)〈q-〉L(i+1).The node input power flow for each element is(8)〈q〉|xR(i)=τ(i,i+1)〈q-〉L(i+1)-[1-γ(i,i+1)]〈q+〉R(i),〈q〉|xL(i+1)=τ(i,i+1)〈q+〉R(i)-[1-γ(i,i+1)]〈q-〉L(i+1).Combine expressions (7) and (8):
(9)[〈q〉|xR(i)〈q〉|xL(i+1)]=α(i,i+1)[cg(i)-cg(i+1)-cg(i)cg(i+1)][〈e〉R(i)〈e〉L(i+1)],
where α(i,i+1)=0.5τ(i,i+1)/γ(i,i+1).As an example, the piezoelectric beam is divided into four units as shown in Figure2 according to the classification rule of energy finite element, in which unit 2 is a piezoelectric beam element, and others are base beam elements. Substitute expression (7) into (9):(10)[k11(1)k12(1)00000k21(1)k22(1)-α(1,2)cg(1)α(1,2)cg(2)00000α(1,2)cg(1)k11(2)-α(1,2)cg(2)k12(2)00000k21(2)k22(2)-α(2,3)cg(2)α(2,3)cg(3)00000α(2,3)cg(2)k11(3)-α(2,3)cg(3)k12(3)00000k21(3)k22(3)+k11(4)k12(4)00000k21(4)k22(4)][e1e2e3e4e5e6e7]=[000000Qin].Figure 2
Element Finite Element model of cantilever beam system.From expression (10), if the material parameters of the base beam and piezoelectric beam elements are known, the energy density of each node could be obtained after measuring or calculating the input power of system.


### 2.3. Transformation of Feedback Element

#### 2.3.1. Ideal System without Time Delay

Equation (10) is deduced without consideration of the feedback element performed by the piezoelectric actuator and sensor. In the negative feedback control model, the piezoelectric actuator would exert a resistance moment to suppress the vibration of host beam according to the vibration response signal picked up by the piezoelectric sensor. In this situation, the active control moment would become a second energy source besides the harmonic exciting force. However, assuming an ideal system in which the transmission of feedback signal and the rapid reaction of the actuator were all instant, that is, there were not any time delay in the feedback control element; an equivalent transformation could be adopted to match the effect of feedback control with some simple modification of system parameters.In most conditions, the attached piezoelectric patches should not be very long in size, and the electric chargeQa(x,t) generated by the piezoelectric sensor under a bending deformation of the host beam could be expressed as
(11)Qa(t)=∫xx+lbd31(Epε1)dx=12d31Ephbb[w2,x(x+l,t)-w2,x(x,t)]≐12d31Ephbblw2,xx(t),
where b,d31,ε1, and Ep represent width, piezoelectric strain constant, and elasticity modulus of piezoelectric element, x and x+l  are coordinates of the piezoelectric patch’s two ends, w2,xx is the curvature of the composite beam.For a proportional feedback control, the bending momentma(t) supplied by the piezoelectric actuator is decided by
(12)ma(t)=KmUa(t)=KmCpQa(t)=Kw2,xx(t),
where Km is a transfer coefficient, Cp is the capacitance of the piezoelectric sensor, and K=0.5d31EphbblKm/Cp is total feedback gain.Taking into account the dynamic balance of moments, when a feedback control went into action, the internal bending momentM in the composite beam would be
(13)M(jω)=Ma(jω)+EI(1+jη)W2,xx(jω)=[K+EI(1+jη)]W2,xx(jω).The above equation is expressed in frequency domain, whereMa and W2,xx are the Fourier transforms of ma and w2,xx, respectively.The revelation of (13) is that an attendance of negative feedback control would be equivalent to an increase of the flexural rigidity of the composite beam by K.


#### 2.3.2. Discussion of Systems with Time Delay

Time delay, which is caused by the transmission and processing of feedback signals and reaction inertia of actuators, is inevitable in all kinds of active control systems. Time delay would result in complicated dynamic characteristics. The following is a qualitative analysis of the influence of time delay to the modeling method investigated here.When a time delay is taken into consideration, the expression (12) should be written in the following description:
(14)ma(t)=Kw2,xx(t-τ).Assuming the system is in a harmonic vibration style, making a Fourier transform to the above equation would result in(15)Ma(jω)=Ke-jωτW2,xx(ω)=KW2,xx(ω)[cos(ωτ)-jsin(ωτ)].Substitute the above equation into (13),
(16)M(jω)={[EI+Kcos(ωτ)]+EI·j[η-Ksin(ωτ)EI]}×W2,xx(jω).Comparing (16) with (13), one could conclude that(1)
ifωτ≪1, that means cos(ωτ)≈1 and sin(ωτ)≈0, (16) and (13) are actually the same and time delay could be neglected;(2)
ifωτ is great enough to be counted, a feedback control accompanying time delay would be equal to a change of flexural rigidity by Kcos(ωτ), and a change of damping loss factor by -Ksin(ωτ)/EI;(3)
from the viewpoint of rigidity improvement, time delay would not be expected. The feedback control would more or less produce a positive effect when0<ωτ<π/2. But if ωτ is too great to exceed π/2, in the range of π/2<ωτ<3π/2, the feedback control would produce a worsen effect than it comes into action before. And this might lead to a negative rigidity if the feedback gain K was great enough, which would be the worst situation that the system would become unstable;(4)
if there was any possibility that one might benefit from a time delay, the benefit might lie in an equivalent increase of damping loss factor whenπ<ωτ<2π, which has been presented in the second part of (16), and this might be helpful to vibration damping. But there would also be a risk that negative damping might appear if improperly controlled, and the system stability would be destroyed.


### 2.4. Input Power of the Piezoelectric Cantilever Beam Using Two Different Calculation Methods

In order to solve expression (10), the input energy Qin should be determined at first. Qin here denotes the energy intensity input by the external excitation to the system and is specifically defined as the time-averaged power input in a vibration period.For steady harmonic vibration,Qin can be expressed as
(17)Qin(ω)=ω2π∫02π/ωRe(F)·Re(V)dt=12|F||V|cosφ=12|F|2Re(1Z)=12|V|2Re(Z),
where F=F(t)ejωt,V=V(t)ejω(t+φ), and Z is the impedance of the structure at the excitation point.The estimation ofQin could be performed by modal method or infinite beam method.


#### 2.4.1. Modal Method

For beam structures as shown in Figure1, the displacement response function is
(18)w(x,t)=F0ejωt∑r=1∞Yr(x0)Yr(x)Mr[ωr2(1+jη)-ω2],
where x0 are the coordinates of the exciting force, and the beam end and Yr(x) and Mr are the modal shape function and modal mass of the beam, respectively.Thus, the input power could be calculated by(19)Qin=12|F|2Re(1Z)=12|F|2Re(jω∑r=1∞Yr(x0)Yr(x)M0[ωr2(1+jη)-ω2]).


#### 2.4.2. Infinite Beam Method

Considering a bilateral infinite beam excited by a harmonic force at pointx0=0, the displacement response in frequency domain could be obtained according to the wave theory:
(20)Wi(x,ω)=Ci1e-jλix+Ci2e-λix+Ci3ejλix+Ci4eλix,
where C1, C2, C3, and C4 are determined by the shear equilibrium conditions as below:
(21)F02-EI∂3w∂x3|x=0=F02+EIKB3(-jC3+C4)=0(x=0-),(22)F02+EI∂3w∂x3|x=0=F02+EIKB3(-jC1+C2)=0(x=0+),(23)KB(jC3+C4)=0,-KB(jC1+C2)=0(x=0).Solving the simultaneous equations (21)~(23),
(24)C1=jC2=C3=jC4=-jF04EIKB3.Thus, the input impedance of bilateral infinite beam at the excitation point is(25)Zs=F0ejωt∂w(0,t)/∂t=2ρACB(1+j)≈2.67ρAClhf(1+j).Similarly, the input impedance of unilateral infinite beam could also be deduced:(26)Zd=F0ejωt∂w(0,t)/∂t=2ρACB(1+j)≈0.67ρAClhf(1+j).The input power of infinite beam could then be acquired by (17).Figure3 shows the results of modal method and infinite beam method. It could be found that the latter is the average value of modal solution over modal density. Figure 4 presents the influence of the damping loss factor on the input power. Increasing the damping will reduce the size of each peak in the input power curve, but the sizes of other regions will increase. Figure 5 compares input power curves under different elastic modulus. It can be observed that reducing elastic modulus will lead to the movement of peaks to lower frequency region. This is because the beam’s vibration modals have been excited earlier after the drop of the stiffness. Reducing the beam density would bring about a similar result as that of increasing the elastic modulus, which is shown in Figure 6.The input power of modal method and infinite beam method.
(a)
(b)Figure 4
Influence of damping on the input power.Figure 5
Influence of stiffness on the input power.Figure 6
Influence of density on the input power.


## 2.1. System Model

Composite piezoelectric cantilever beam is one of the most typical piezoelectric intelligent structures, which is presented in Figure1(a). In the demonstrated system, series of piezoelectric material patches are grouped in couples and attached symmetrically to a host beam, one patch of each couple being used as an actuator and the other a sensor. The host beam is excited by a harmonically varying force with magnitude F0 which would cause a transverse bending vibration. hp and hb represent the thickness of piezoelectric patch and cantilever beam, respectively.Piezoelectric cantilever system model and energy flow balance.
(a)
(b)Figure1(b) shows the vibration energy inflow and outflow of cubic element, in which q and πdiss represent power flow through cubic element and the energy dissipated by damping, respectively, and Qin is the energy input from external. At each point of the structure, the instantaneous energy balance equation is given by  [8]
(1)∂e∂t=-∂q∂x-πdiss+Qin,
where ∂e/∂t is the variation ratio of energy density to time.The time-averaged value ofq, πdiss, and e, that is, 〈q〉, 〈πdiss〉, and 〈e〉, could be related as
(2)〈q〉=(cg2ηω)·d〈e〉dx,〈πdiss〉=ηω〈e〉.By substituting expression (2) into  (1), the energy balance equation is written as
(3)(cg2ηω)·(d2〈e〉dx2)+ηω〈e〉=〈Qin〉,
where cg=2(ω2EI/ρA) is the group velocity of bending wave and η is the damping loss factor.


## 2.2. Energy Density Control Equations without Feedback

When the feedback does not exist, the harmonic excitation is the only energy source in the system. Equation (3) could be discretized by Galerkin weighted residual method:
(4)K(n)e(n)=P(n)+Q(n),(5)K(n)=∫xL(n)xR(n)(-cg2ωηdNTdxdNdx+ωηNTN)dx,(6)P(n)=N(xR(n))〈q〉|xR(n)-N(xL(n))〈q〉|xL(n),Q(n)=∫xixjN〈Qin〉dx,
where  the superscript (n) denotes the nth discretizednodes, K(n) is a coefficient matrix related to stiffness and mass, P(n) is the power flow at corresponding node, Q(n) is the external energy input, and N(x) is the lagrange linear interpolation basic function.For solving expression (4), coupling relationship among elements needs to be solved. Suppose the numbers of two adjacent elements are i and i+1 and  τ(i,i+1) and γ(i,i+1) represent transmission and reflection coefficient between coupling nodes:
(7)〈q-〉R(i)=γ(i,i+1)〈q+〉R(i)+τ(i,i+1)〈q-〉L(i+1),〈q+〉L(i+1)=τ(i,,i+1)〈q+〉R(i)+γ(i,i+1)〈q-〉L(i+1).The node input power flow for each element is(8)〈q〉|xR(i)=τ(i,i+1)〈q-〉L(i+1)-[1-γ(i,i+1)]〈q+〉R(i),〈q〉|xL(i+1)=τ(i,i+1)〈q+〉R(i)-[1-γ(i,i+1)]〈q-〉L(i+1).Combine expressions (7) and (8):
(9)[〈q〉|xR(i)〈q〉|xL(i+1)]=α(i,i+1)[cg(i)-cg(i+1)-cg(i)cg(i+1)][〈e〉R(i)〈e〉L(i+1)],
where α(i,i+1)=0.5τ(i,i+1)/γ(i,i+1).As an example, the piezoelectric beam is divided into four units as shown in Figure2 according to the classification rule of energy finite element, in which unit 2 is a piezoelectric beam element, and others are base beam elements. Substitute expression (7) into (9):(10)[k11(1)k12(1)00000k21(1)k22(1)-α(1,2)cg(1)α(1,2)cg(2)00000α(1,2)cg(1)k11(2)-α(1,2)cg(2)k12(2)00000k21(2)k22(2)-α(2,3)cg(2)α(2,3)cg(3)00000α(2,3)cg(2)k11(3)-α(2,3)cg(3)k12(3)00000k21(3)k22(3)+k11(4)k12(4)00000k21(4)k22(4)][e1e2e3e4e5e6e7]=[000000Qin].Figure 2
Element Finite Element model of cantilever beam system.From expression (10), if the material parameters of the base beam and piezoelectric beam elements are known, the energy density of each node could be obtained after measuring or calculating the input power of system.


## 2.3. Transformation of Feedback Element

### 2.3.1. Ideal System without Time Delay

Equation (10) is deduced without consideration of the feedback element performed by the piezoelectric actuator and sensor. In the negative feedback control model, the piezoelectric actuator would exert a resistance moment to suppress the vibration of host beam according to the vibration response signal picked up by the piezoelectric sensor. In this situation, the active control moment would become a second energy source besides the harmonic exciting force. However, assuming an ideal system in which the transmission of feedback signal and the rapid reaction of the actuator were all instant, that is, there were not any time delay in the feedback control element; an equivalent transformation could be adopted to match the effect of feedback control with some simple modification of system parameters.In most conditions, the attached piezoelectric patches should not be very long in size, and the electric chargeQa(x,t) generated by the piezoelectric sensor under a bending deformation of the host beam could be expressed as
(11)Qa(t)=∫xx+lbd31(Epε1)dx=12d31Ephbb[w2,x(x+l,t)-w2,x(x,t)]≐12d31Ephbblw2,xx(t),
where b,d31,ε1, and Ep represent width, piezoelectric strain constant, and elasticity modulus of piezoelectric element, x and x+l  are coordinates of the piezoelectric patch’s two ends, w2,xx is the curvature of the composite beam.For a proportional feedback control, the bending momentma(t) supplied by the piezoelectric actuator is decided by
(12)ma(t)=KmUa(t)=KmCpQa(t)=Kw2,xx(t),
where Km is a transfer coefficient, Cp is the capacitance of the piezoelectric sensor, and K=0.5d31EphbblKm/Cp is total feedback gain.Taking into account the dynamic balance of moments, when a feedback control went into action, the internal bending momentM in the composite beam would be
(13)M(jω)=Ma(jω)+EI(1+jη)W2,xx(jω)=[K+EI(1+jη)]W2,xx(jω).The above equation is expressed in frequency domain, whereMa and W2,xx are the Fourier transforms of ma and w2,xx, respectively.The revelation of (13) is that an attendance of negative feedback control would be equivalent to an increase of the flexural rigidity of the composite beam by K.


### 2.3.2. Discussion of Systems with Time Delay

Time delay, which is caused by the transmission and processing of feedback signals and reaction inertia of actuators, is inevitable in all kinds of active control systems. Time delay would result in complicated dynamic characteristics. The following is a qualitative analysis of the influence of time delay to the modeling method investigated here.When a time delay is taken into consideration, the expression (12) should be written in the following description:
(14)ma(t)=Kw2,xx(t-τ).Assuming the system is in a harmonic vibration style, making a Fourier transform to the above equation would result in(15)Ma(jω)=Ke-jωτW2,xx(ω)=KW2,xx(ω)[cos(ωτ)-jsin(ωτ)].Substitute the above equation into (13),
(16)M(jω)={[EI+Kcos(ωτ)]+EI·j[η-Ksin(ωτ)EI]}×W2,xx(jω).Comparing (16) with (13), one could conclude that(1)
ifωτ≪1, that means cos(ωτ)≈1 and sin(ωτ)≈0, (16) and (13) are actually the same and time delay could be neglected;(2)
ifωτ is great enough to be counted, a feedback control accompanying time delay would be equal to a change of flexural rigidity by Kcos(ωτ), and a change of damping loss factor by -Ksin(ωτ)/EI;(3)
from the viewpoint of rigidity improvement, time delay would not be expected. The feedback control would more or less produce a positive effect when0<ωτ<π/2. But if ωτ is too great to exceed π/2, in the range of π/2<ωτ<3π/2, the feedback control would produce a worsen effect than it comes into action before. And this might lead to a negative rigidity if the feedback gain K was great enough, which would be the worst situation that the system would become unstable;(4)
if there was any possibility that one might benefit from a time delay, the benefit might lie in an equivalent increase of damping loss factor whenπ<ωτ<2π, which has been presented in the second part of (16), and this might be helpful to vibration damping. But there would also be a risk that negative damping might appear if improperly controlled, and the system stability would be destroyed.


## 2.3.1. Ideal System without Time Delay

Equation (10) is deduced without consideration of the feedback element performed by the piezoelectric actuator and sensor. In the negative feedback control model, the piezoelectric actuator would exert a resistance moment to suppress the vibration of host beam according to the vibration response signal picked up by the piezoelectric sensor. In this situation, the active control moment would become a second energy source besides the harmonic exciting force. However, assuming an ideal system in which the transmission of feedback signal and the rapid reaction of the actuator were all instant, that is, there were not any time delay in the feedback control element; an equivalent transformation could be adopted to match the effect of feedback control with some simple modification of system parameters.In most conditions, the attached piezoelectric patches should not be very long in size, and the electric chargeQa(x,t) generated by the piezoelectric sensor under a bending deformation of the host beam could be expressed as
(11)Qa(t)=∫xx+lbd31(Epε1)dx=12d31Ephbb[w2,x(x+l,t)-w2,x(x,t)]≐12d31Ephbblw2,xx(t),
where b,d31,ε1, and Ep represent width, piezoelectric strain constant, and elasticity modulus of piezoelectric element, x and x+l  are coordinates of the piezoelectric patch’s two ends, w2,xx is the curvature of the composite beam.For a proportional feedback control, the bending momentma(t) supplied by the piezoelectric actuator is decided by
(12)ma(t)=KmUa(t)=KmCpQa(t)=Kw2,xx(t),
where Km is a transfer coefficient, Cp is the capacitance of the piezoelectric sensor, and K=0.5d31EphbblKm/Cp is total feedback gain.Taking into account the dynamic balance of moments, when a feedback control went into action, the internal bending momentM in the composite beam would be
(13)M(jω)=Ma(jω)+EI(1+jη)W2,xx(jω)=[K+EI(1+jη)]W2,xx(jω).The above equation is expressed in frequency domain, whereMa and W2,xx are the Fourier transforms of ma and w2,xx, respectively.The revelation of (13) is that an attendance of negative feedback control would be equivalent to an increase of the flexural rigidity of the composite beam by K.


## 2.3.2. Discussion of Systems with Time Delay

Time delay, which is caused by the transmission and processing of feedback signals and reaction inertia of actuators, is inevitable in all kinds of active control systems. Time delay would result in complicated dynamic characteristics. The following is a qualitative analysis of the influence of time delay to the modeling method investigated here.When a time delay is taken into consideration, the expression (12) should be written in the following description:
(14)ma(t)=Kw2,xx(t-τ).Assuming the system is in a harmonic vibration style, making a Fourier transform to the above equation would result in(15)Ma(jω)=Ke-jωτW2,xx(ω)=KW2,xx(ω)[cos(ωτ)-jsin(ωτ)].Substitute the above equation into (13),
(16)M(jω)={[EI+Kcos(ωτ)]+EI·j[η-Ksin(ωτ)EI]}×W2,xx(jω).Comparing (16) with (13), one could conclude that(1)
ifωτ≪1, that means cos(ωτ)≈1 and sin(ωτ)≈0, (16) and (13) are actually the same and time delay could be neglected;(2)
ifωτ is great enough to be counted, a feedback control accompanying time delay would be equal to a change of flexural rigidity by Kcos(ωτ), and a change of damping loss factor by -Ksin(ωτ)/EI;(3)
from the viewpoint of rigidity improvement, time delay would not be expected. The feedback control would more or less produce a positive effect when0<ωτ<π/2. But if ωτ is too great to exceed π/2, in the range of π/2<ωτ<3π/2, the feedback control would produce a worsen effect than it comes into action before. And this might lead to a negative rigidity if the feedback gain K was great enough, which would be the worst situation that the system would become unstable;(4)
if there was any possibility that one might benefit from a time delay, the benefit might lie in an equivalent increase of damping loss factor whenπ<ωτ<2π, which has been presented in the second part of (16), and this might be helpful to vibration damping. But there would also be a risk that negative damping might appear if improperly controlled, and the system stability would be destroyed.


## 2.4. Input Power of the Piezoelectric Cantilever Beam Using Two Different Calculation Methods

In order to solve expression (10), the input energy Qin should be determined at first. Qin here denotes the energy intensity input by the external excitation to the system and is specifically defined as the time-averaged power input in a vibration period.For steady harmonic vibration,Qin can be expressed as
(17)Qin(ω)=ω2π∫02π/ωRe(F)·Re(V)dt=12|F||V|cosφ=12|F|2Re(1Z)=12|V|2Re(Z),
where F=F(t)ejωt,V=V(t)ejω(t+φ), and Z is the impedance of the structure at the excitation point.The estimation ofQin could be performed by modal method or infinite beam method.


### 2.4.1. Modal Method

For beam structures as shown in Figure1, the displacement response function is
(18)w(x,t)=F0ejωt∑r=1∞Yr(x0)Yr(x)Mr[ωr2(1+jη)-ω2],
where x0 are the coordinates of the exciting force, and the beam end and Yr(x) and Mr are the modal shape function and modal mass of the beam, respectively.Thus, the input power could be calculated by(19)Qin=12|F|2Re(1Z)=12|F|2Re(jω∑r=1∞Yr(x0)Yr(x)M0[ωr2(1+jη)-ω2]).


### 2.4.2. Infinite Beam Method

Considering a bilateral infinite beam excited by a harmonic force at pointx0=0, the displacement response in frequency domain could be obtained according to the wave theory:
(20)Wi(x,ω)=Ci1e-jλix+Ci2e-λix+Ci3ejλix+Ci4eλix,
where C1, C2, C3, and C4 are determined by the shear equilibrium conditions as below:
(21)F02-EI∂3w∂x3|x=0=F02+EIKB3(-jC3+C4)=0(x=0-),(22)F02+EI∂3w∂x3|x=0=F02+EIKB3(-jC1+C2)=0(x=0+),(23)KB(jC3+C4)=0,-KB(jC1+C2)=0(x=0).Solving the simultaneous equations (21)~(23),
(24)C1=jC2=C3=jC4=-jF04EIKB3.Thus, the input impedance of bilateral infinite beam at the excitation point is(25)Zs=F0ejωt∂w(0,t)/∂t=2ρACB(1+j)≈2.67ρAClhf(1+j).Similarly, the input impedance of unilateral infinite beam could also be deduced:(26)Zd=F0ejωt∂w(0,t)/∂t=2ρACB(1+j)≈0.67ρAClhf(1+j).The input power of infinite beam could then be acquired by (17).Figure3 shows the results of modal method and infinite beam method. It could be found that the latter is the average value of modal solution over modal density. Figure 4 presents the influence of the damping loss factor on the input power. Increasing the damping will reduce the size of each peak in the input power curve, but the sizes of other regions will increase. Figure 5 compares input power curves under different elastic modulus. It can be observed that reducing elastic modulus will lead to the movement of peaks to lower frequency region. This is because the beam’s vibration modals have been excited earlier after the drop of the stiffness. Reducing the beam density would bring about a similar result as that of increasing the elastic modulus, which is shown in Figure 6.The input power of modal method and infinite beam method.
(a)
(b)Figure 4
Influence of damping on the input power.Figure 5
Influence of stiffness on the input power.Figure 6
Influence of density on the input power.


## 2.4.1. Modal Method

For beam structures as shown in Figure1, the displacement response function is
(18)w(x,t)=F0ejωt∑r=1∞Yr(x0)Yr(x)Mr[ωr2(1+jη)-ω2],
where x0 are the coordinates of the exciting force, and the beam end and Yr(x) and Mr are the modal shape function and modal mass of the beam, respectively.Thus, the input power could be calculated by(19)Qin=12|F|2Re(1Z)=12|F|2Re(jω∑r=1∞Yr(x0)Yr(x)M0[ωr2(1+jη)-ω2]).


## 2.4.2. Infinite Beam Method

Considering a bilateral infinite beam excited by a harmonic force at pointx0=0, the displacement response in frequency domain could be obtained according to the wave theory:
(20)Wi(x,ω)=Ci1e-jλix+Ci2e-λix+Ci3ejλix+Ci4eλix,
where C1, C2, C3, and C4 are determined by the shear equilibrium conditions as below:
(21)F02-EI∂3w∂x3|x=0=F02+EIKB3(-jC3+C4)=0(x=0-),(22)F02+EI∂3w∂x3|x=0=F02+EIKB3(-jC1+C2)=0(x=0+),(23)KB(jC3+C4)=0,-KB(jC1+C2)=0(x=0).Solving the simultaneous equations (21)~(23),
(24)C1=jC2=C3=jC4=-jF04EIKB3.Thus, the input impedance of bilateral infinite beam at the excitation point is(25)Zs=F0ejωt∂w(0,t)/∂t=2ρACB(1+j)≈2.67ρAClhf(1+j).Similarly, the input impedance of unilateral infinite beam could also be deduced:(26)Zd=F0ejωt∂w(0,t)/∂t=2ρACB(1+j)≈0.67ρAClhf(1+j).The input power of infinite beam could then be acquired by (17).Figure3 shows the results of modal method and infinite beam method. It could be found that the latter is the average value of modal solution over modal density. Figure 4 presents the influence of the damping loss factor on the input power. Increasing the damping will reduce the size of each peak in the input power curve, but the sizes of other regions will increase. Figure 5 compares input power curves under different elastic modulus. It can be observed that reducing elastic modulus will lead to the movement of peaks to lower frequency region. This is because the beam’s vibration modals have been excited earlier after the drop of the stiffness. Reducing the beam density would bring about a similar result as that of increasing the elastic modulus, which is shown in Figure 6.The input power of modal method and infinite beam method.
(a)
(b)Figure 4
Influence of damping on the input power.Figure 5
Influence of stiffness on the input power.Figure 6
Influence of density on the input power.


## 3. Numerical Analysis

A concrete example of the system in Figure2 is investigated. The physical parameters are listed in Table 1.Table 1
Parameters of piezoelectric cantilever beam.
Length
Section areaA
Moment of inertiaI
Elasticity modulusE
Densityρ
Beam
L
b
=
0.5
A
b
=
4
e
-
5
I
b
=
1.33
e
-
11
E
b
=
2.06
e
11
ρ
b
=
7.8
e
3
PZ (Case1)
L
a
=
0.1
A
a
=
2
e
-
5
I
a
=
1.67
e
-
12
E
a
=
0.7
e
11
ρ
a
=
7.75
e
3
PZ (Case2)
L
a
=
0.05
A
a
=
1
e
-
5
I
a
=
2.08
e
-
13
E
a
=
0.7
e
11
ρ
a
=
7.75
e
3The EFEM model would be verified through a comparison with the result from a wave theory which has been deduced in the appendix. As a first stage, Table2 gives a comparison between natural frequencies of the composite piezoelectric cantilever beam calculated by wave theory and ANSYS, respectively. It could also be found that the natural frequencies are   little higher than those of the host beam after the path of piezoelectric elements.Table 2
Natural frequency of the system before and after the patch of piezoelectric elements.
Mode number
1
2
3
4
5
Natural frequency of host beam (wave theory)/Hz
13.55
84.94
237.8
466.1
770.4
Natural frequency of host beam (FEM)/Hz
13.55
84.94
237.8
465.9
770.2
Natural frequency of piezoelectric beam (FEM)/Hz
15.25
92.19
250.7
480.4
784.2


### 3.1. Comparison of EFEM and Exact Solutions by Wave Theory

Figure7 gives a comparison between the energy density calculated by EFEM and the wave theory in the appendix. The average value of wave solution over the beam length, or in other words, EFEM solution is a time-averaged and space-averaged energy density. Figure 7 also indicates an excellent identity in the results of EFEM and wave method. In Figure 8 different feedback gain is applied in the wave method and the effect of feedback control is evaluated through the variation of vibration amplitude. A careful comparison of Figures 8 and 9(a) could further verify the coincidence of EFEM and wave theory.Figure 7
Comparison between solutions of EFEM and wave theory.Figure 8
Vibration amplitude under different feedback gain solved by wave theory.Influence of variation of system parameters on energy density distribution by EFEM.
(a)
Comparison of different feedback gains
(b)
Comparison of different piezoelectric patch positions
(c)
Comparison of different piezoelectric patch lengths
(d)
Comparison of different piezoelectric patch pairs


### 3.2. The Influence of Variation of System Parameter on the Distribution of Energy Density

Figure9 is an investigation about the influence of variation of system parameters on distribution of system energy density along the piezoelectric cantilever beam by FFEM solution. As can be seen in Figure 9(a), increasing the feedback gain can reduce the energy density value of the beam. By comprehensive comparison of Figure 8 with Figure 9(a), one can still conclude that the energy density and vibration amplitude are equivalent in describing vibration level. However, the energy density, on account of having been space averaged, is more recapitulative in vibration evaluation. And the EFEM could be utilized to deal with more complex engineering structures than classical wave method.Figure9(b) shows the energy density distribution when the piezoelectric patch is located at different positions. It seems that the active control would be more effective if the piezoelectric patches were fixed near the clamped end of the cantilever beam. Figure 9(c) shows the alteration of energy density distribution when the length of piezoelectric element is altered, and it presents that a bit longer piezoelectric patch might be more beneficial to vibration suppression. Figure 9(d) compares the situation that two pairs of piezoelectric patches are attached discretely to the host beam with only one pair of piezoelectric patches being used. Obviously, lower energy density could be achieved in the former situation. And it could be realized that the beam segments having been patched with piezoelectric elements have lower energy density than other unpatched beam segments, and a more satisfied control effect would be expected if more pairs of piezoelectric patches were applied.


## 3.1. Comparison of EFEM and Exact Solutions by Wave Theory

Figure7 gives a comparison between the energy density calculated by EFEM and the wave theory in the appendix. The average value of wave solution over the beam length, or in other words, EFEM solution is a time-averaged and space-averaged energy density. Figure 7 also indicates an excellent identity in the results of EFEM and wave method. In Figure 8 different feedback gain is applied in the wave method and the effect of feedback control is evaluated through the variation of vibration amplitude. A careful comparison of Figures 8 and 9(a) could further verify the coincidence of EFEM and wave theory.Figure 7
Comparison between solutions of EFEM and wave theory.Figure 8
Vibration amplitude under different feedback gain solved by wave theory.Influence of variation of system parameters on energy density distribution by EFEM.
(a)
Comparison of different feedback gains
(b)
Comparison of different piezoelectric patch positions
(c)
Comparison of different piezoelectric patch lengths
(d)
Comparison of different piezoelectric patch pairs


## 3.2. The Influence of Variation of System Parameter on the Distribution of Energy Density

Figure9 is an investigation about the influence of variation of system parameters on distribution of system energy density along the piezoelectric cantilever beam by FFEM solution. As can be seen in Figure 9(a), increasing the feedback gain can reduce the energy density value of the beam. By comprehensive comparison of Figure 8 with Figure 9(a), one can still conclude that the energy density and vibration amplitude are equivalent in describing vibration level. However, the energy density, on account of having been space averaged, is more recapitulative in vibration evaluation. And the EFEM could be utilized to deal with more complex engineering structures than classical wave method.Figure9(b) shows the energy density distribution when the piezoelectric patch is located at different positions. It seems that the active control would be more effective if the piezoelectric patches were fixed near the clamped end of the cantilever beam. Figure 9(c) shows the alteration of energy density distribution when the length of piezoelectric element is altered, and it presents that a bit longer piezoelectric patch might be more beneficial to vibration suppression. Figure 9(d) compares the situation that two pairs of piezoelectric patches are attached discretely to the host beam with only one pair of piezoelectric patches being used. Obviously, lower energy density could be achieved in the former situation. And it could be realized that the beam segments having been patched with piezoelectric elements have lower energy density than other unpatched beam segments, and a more satisfied control effect would be expected if more pairs of piezoelectric patches were applied.


## 4. Summary and Conclusions

The research in this paper is focused on the modeling and analysis of active vibration control system of the composite piezoelectric beam based on Energy Finite Element Method. It should be pointed that environment noise, signal transmission error, and other uncontrollable factors might cause more or less departure of theoretical model from the practice, which is difficult to forecast   now. However, through the theoretical and numerical analysis in this work, it could be proved that the EFEM is of application value in dealing with similar systems and it could be concluded that(1)
the EFEM gives a time-averaged and space-averaged value of energy density at the conjunction node between discretized beam elements. It could be equally efficient as classical wave method or other mature method such as FEM and SEA in vibration evaluation, but the method would be more suitable for dealing with complex engineering structures associating broad frequency band;(2)
in modeling the piezoelectric intelligent beam, the feedback gain could be equivalent to an extra increase of the flexural rigidity of the composite beam segment; that is, greater feedback gain would provide greater reduction of system vibration level, whereas the introduction of time delay would bring about complicated situations. In general, both rigidity modification and damping modification should be taken into consideration, and attention should be paid to potential problems on negative rigidity and damping;(3)
numerical analysis indicates that proper design of system configuration parameters would be essential in order to achieve the best control efficiency, and how to find the optimal system parameters would be of significant value to be explored further. For the piezoelectric cantilever beam as studied in this paper, the control efficiency would be raised by arranging the piezoelectric patches near-clamped-end in position, a bit longer in size and multipled in number.



---

*Source: 819725-2012-12-20.xml*