# Research on Spillover Effects for Vibration Control of Piezoelectric Smart Structures by ANSYS

Authors: Xingjian Dong; Zhike Peng; Wenming Zhang; HongXing Hua; Guang Meng
Journal: Mathematical Problems in Engineering
 (2014)

---

## Abstract

To control vibration of a piezoelectric smart structure, a controller is usually designed based on a reduced order model (ROM) of the system. When such a ROM based controller operates in closed loop with the actual structure, spillover phenomenon occurs because the unmodeled dynamics, which are not included in ROM, will be excited. In this paper, a new approach aiming at investigating spillover effects in ANSYS software is presented. By using the ANSYS parametric design language (APDL), the ROM based controller is integrated into finite element model to provide an accurate representation of what will happen when the controller is connected to the real plant. Therefore, the issues of spillover effects can be addressed in the closed loop simulation. Numerical examples are presented for investigating spillover effects of a cantilever piezoelectric plate subjected to various types of loading. The importance of considering spillover effects in closed loop simulation of piezoelectric smart structures is demonstrated. Moreover, the present study may provide an efficient method especially beneficial for preliminary design of piezoelectric smart structure to evaluate the performance of candidate control laws in finite element environment considering spillover effects.


---

## Body

## 1. Introduction

In the recent years, the active vibration control of piezoelectric smart structures has been intensively investigated due to their potential benefits in shape control and vibration suppression of lightweight flexible structures. The piezoelectric materials, such as lead zirconate titanate (PZT), have the property to generate electrical charge under mechanical load or deformation and the reverse; applying an electrical field to the material results in mechanical strains or stresses. Bonding or embedding piezoelectric patches in a structure can act as sensors to monitor or as actuators to control the response of the structure. It has been observed from the open literature that current analytical and numerical investigations focus on vibration control of thin components such as beams, plates, and shells with embedded or surface bonded piezoelectric patches which are used as sensors or actuators.To control vibration of piezoelectric smart structures, both structural dynamics and control theory need to be considered. It is widely accepted that finite element (FE) method is an efficient tool for analyzing piezoelectric smart structures even with complicated boundary conditions. Comprehensive surveys on FE modeling of piezoelectric smart structures can be found in [1–3]. Among all these works, the simplest and often used model is the equivalent single layer model, which includes the classical laminated plate theory [4, 5], first order shear deformation theory [6, 7], and third order theory [8, 9]. Another model for thick laminated composite structures is layer-wise theory [10, 11]. Besides, some 3D models [12] are also available. These models have their own advantages and disadvantages in terms of accuracy, speed of convergence, and computational cost. Moreover, elements allowing the simulation of the electromechanical coupling properties of piezoelectric materials are routinely implemented in some commercially available FE analysis codes, such as ANSYS. With respect to the control of piezoelectric smart structures, conventional control methods have worked well in the past. A survey on various control algorithms employed for vibration control of smart structures has been presented by Alkhatib and Golnaraghi and Fei et al. [13, 14]. Among various control algorithms, classical constant gain velocity feedback (CGVF) control [15–18] is widely employed for actively vibration control of piezoelectric smart beams and plates. Vasques and Dias Rodrigues [19] present a comparison between the CGVF control algorithm and optimal control strategies (linear quadratic regulator, i.e., LQR, and linear quadratic Gaussian controller, i.e., LQG) to investigate their effectiveness to suppress vibrations in piezoelectric smart structures.However, conventional controller is usually designed from a reduced order model (ROM) of the original system, whereas FE models inevitably have a large number of degrees of freedom (DOFs). FE models are appropriate for dynamic analysis but cannot be used directly as a basis for controller design. There is a gap between dynamic analysis and control design of a piezoelectric smart structure. To this end, some methods or combination of methods for creating a ROM must be employed to compensate the gap. The most widely used model reduction method to try on a large scale FE model is the mode superposition technique [20–22]. Once the ROM is created, a finite dimensional controller can then be designed about it using conventional techniques. Usually a ROM focuses on the low frequency modes, which are extremely important in the overall dynamic characteristics of the structure. Given a ROM representation of the system, the important modes are selected while all others which are not included in the ROM are defined as unmodeled dynamics or residual modes. Obviously, the ROM based finite dimensional controller is designed considering only a few lower modes of a piezoelectric smart structure. When such a ROM based controller is applied to the full order system, actuator forces for reducing the vibration of the lower modes will also influence the residual modes of the structure and produce undesirable vibration due to the unmodeled dynamics. This phenomenon is known as control spillover. Similarly, the sensor will sense the deflection not only from the lower modes but also from the other modes, giving rise to the so-called observation spillover. Spillover effects are undesirable and may cause performance degradation and even system instability. To address the issues of spillover, one popular approach is to use the residual mode filter (RMF) technique [23–28] to block the residual modes from entering the feedback loop and eliminate the spillover effects. In most of these researches, if not all, the design of RMF and the ROM based controller are conducted based on modal space model. This process needs a high order modal state space model to effectively represent the original system. As a result, dynamic analysis of the piezoelectric smart structure is based on the FE model while controller evaluation process is based on modal model. Considering that FE model can accurately predict the dynamic behaviour of the smart structure, it is desired to investigate spillover effects in FE environment. Recently, some articles [29–31] have shown the feasibility and efficiency of integration of controller into ANSYS to perform closed loop simulations.The present study is a continuation of the above mentioned works [30, 31], with focus on demonstrating the spillover effects in FE environment. First, a so-called subspace identification method is employed to develop a state space ROM for control law design from the output of ANSYS. An optimal control scheme coupled with a Luenberger state estimator is then designed based on the developed ROM. Next, ANSYS parametric design language (APDL) is then used to implement the ROM based controller inside the ANSYS software package by updating the actuator voltages according to the controller at each time instant during the transient analysis of FE simulations. In this way, the spillover effects can be taken into account because the unmodeled dynamics are already included in the FE model of piezoelectric structure. Finally, the proposed approach is used for investigating spillover effects of a cantilever piezoelectric plate subjected to various types of loading. The spillover phenomena, which are a crucial problem for active vibration control of piezoelectric smart structures, are illustrated via power spectral density (PSD) of the displacement response. The importance of including spillover effects in closed loop simulation of piezoelectric smart structures is discussed herein. Simulation results have shown the feasibility and efficiency of the proposed approach in investigating spillover effects in FE environment and its potential benefits in vibration control of large scale flexible structures.


## 2. Formulation of a ROM for a Piezoelectric Plate

### 2.1. Finite Element Model

The constitutive equations of a deformable piezoelectric material in Voigt notation, coupling the elastic, and the electric fields can be expressed as(1)
σ
=
c
ε
-
e
T
E
,
D
=
e
ε
+
ϵ
ε
E
,
where σ and ε represent stress vector and strain vector, respectively. E is the electric field vector, D is the electric displacement vector, c is the elastic coefficient matrix, ϵ
ε is the dielectric coefficient matrix at constant ε, and e is the piezoelectric stress coefficient matrix.By applying element discretization and using the Hamilton’s principle, the global dynamic system equation for a piezoelectric smart structure is derived as(2)
[
M
u
u
0
0
0
]
[
q
¨
𝚽
¨
]
+
[
C
u
u
0
0
0
]
[
q
˙
𝚽
˙
]
+
[
K
u
u
K
u
ϕ
K
u
ϕ
T
K
ϕ
ϕ
]
[
q
𝚽
]
=
{
F
u
F
ϕ
}
,
where q denotes the global nodal displacement vector, 𝚽 denotes electric potential vector, and a dot above a variable denotes a time derivative. F
u and F
ϕ denote the nodal external force vector and the nodal externally applied charge vector, respectively. M
u
u and K
u
u are the structural mass and stiffness matrix, respectively, and K
u
ϕ is the piezoelectric coupling matrices; K
ϕ
ϕ is the dielectric permittivity matrix and C
u
u is the damping matrix, which is usually assumed to be constructed as the weighted sum of mass and stiffness matrices
(3)
C
u
u
=
α
M
u
u
+
β
K
u
u
,
where α and β are Rayleigh’s coefficients which are nonnegative and small in the sense that the effective damping coefficient is ≪1.A rectangular aluminum cantilever plate which consists of a host aluminum plate (700 mm × 150 mm × 1.2 mm) and 3 pairs of PZT patches (70 mm × 15 mm × 0.5 mm), as shown in Figure1, is considered. The 3 PZT patches bonded to the top surface are to be used as actuators and 3 other PZT patches bonded symmetrically to the bottom surface are to be used as sensors. The material properties of the host plate and PZT patches are given in Table 1. ANSYS software package is used to model the cantilever plate with PZT patches. In this study, the 8-node SOLID5 coupled field element and 8-node SOLID45 element are used to model the PZT patches and the host structure, respectively. The cantilever boundary conditions are simulated by specifying zero displacement DOF for the x, y, and z directions at the fixed end of the plate. The electrical DOFs of the nodes on the outside surface of a PZT patch are coupled by the CP command so that the nodes on the outside surface represent an electrode. The electrode on the inner surface of a PZT patch is defined in the same manner and a fixed zero voltage potential is defined on the inner surface of the PZT patch to model the ground electrode. A total of 17248 three-dimensional solid elements are used to mesh the structure. Among those elements, 1008 SOLID5 elements are utilized to model piezoelectric patches, and 16240 SOLID45 elements are employed to model host structure. Rayleigh type damping is considered with coefficients α
=
0.0 rad s−1 and β
=
0.003 s.Table 1
Material properties.
PZT
Host plate
ρ (kg/m3)
7600.0
2700
E
11
=
E
22  (GPa)
63.0
69.0
E
33 (GPa)
54.0
69.0
G
12 (GPa)
24.6
27.0
G
23
=
G
13  (GPa)
30.6
27.0
ν
12
0.3
0.32
ν
13
=
ν
23
0.4
0.32
e
31, e
32, e
33, e
15 (C/m2)
−5.4, −5.4, 15.8, 12.3
/
ϵ
11
ε, ϵ
22
ε, ϵ
33
ε (×10−8 F/m)
1.53, 1.53, 1.50
/Figure 1
Configuration of a cantilever plate with PZT S/A pairs.


### 2.2. The Reduced Representation of the System

FE model is appropriate for dynamic analysis but cannot be used directly as a basis for controller design. Considering that FE simulations are analogous to performing experimental investigations where the only direct outputs are the time histories, a so-called subspace identification method is employed to develop a state space ROM for the piezoelectric smart structure. In this way, the multiple-input multiple-output state space model is identified, which is convenient for the following controller design procedure.Suppose that the FE model of a piezoelectric smart structure can be represented by the following state space model:(4)
x
k
+
1
=
A
x
k
+
B
u
k
,
y
k
=
C
x
k
+
D
u
k
.
Here, x
k
∈
ℝ
n is the state vector at discrete time instant k and contains the numerical values of n states. u
k
∈
ℝ
m and y
k
∈
ℝ
l are the measurements at time instant k of, respectively, the m inputs and l outputs. As to the piezoelectric smart structure mentioned above, the inputs are the externally exerted mechanical force and the applied voltages over the actuators and the outputs are the sensor voltages. A
∈
ℝ
n
×
n, B
∈
ℝ
n
×
m, C
∈
ℝ
l
×
n, and D
∈
ℝ
l
×
m are state matrix, input matrix, output matrix, and direct influence matrix, respectively. It is assumed that the matrix pair (
A
,
C
) is observable and (
A
,
B
) is controllable. The objective of the subspace identification method is to determine the order of the unknown system and the system matrices A, B, C, and D from given input-output data of the FE model.The input-output relations of the linear control system described by (4) can be expressed compactly in a matrix form as [32, 33]
(5)
Y
=
Γ
X
+
H
U
,
where the extended observability matrix Γ is defined as
(6)
Γ
=
[
C
C
A
⋯
C
A
i
-
1
]
∈
ℝ
l
i
×
n
.
H is the lower block triangular Toeplitz matrix of impulse response from u to y given by
(7)
H
=
[
D
0
0
⋯
0
C
B
D
0
⋯
0
C
A
B
C
B
D
⋯
0
⋯
⋯
⋯
⋯
⋯
C
A
i
-
2
B
C
A
i
-
3
B
C
A
i
-
4
B
⋯
D
]
∈
ℝ
l
i
×
m
i
.
U and Y are, respectively, the input and output block Hankel matrices, which are defined as
(8)
U
=
[
u
0
u
1
u
2
⋯
u
j
-
1
u
1
u
2
u
3
⋯
u
j
⋯
⋯
⋯
⋯
⋯
u
i
-
1
u
i
u
i
+
1
⋯
u
i
+
j
-
2
]
∈
ℝ
m
i
×
j
,
Y
=
[
y
0
y
1
y
2
⋯
y
j
-
1
y
1
y
2
y
3
⋯
y
j
⋯
⋯
v
⋯
⋯
y
i
-
1
y
i
y
i
+
1
⋯
y
i
+
j
-
2
]
∈
ℝ
l
i
×
j
.
The state sequence X is defined as
(9)
X
=
[
x
0
x
1
⋯
x
j
-
1
]
∈
ℝ
n
×
j
.
The number of block rows (i) is a user specified parameter which is required to be larger than the maximum order of the system one wants to identify. The number of columns (j) is assumed to approach ∞.The key idea of subspace identification method is to estimate the extended observability matrixΓ by computing the singular value decomposition (SVD) of U in the first step:
(10)
U
=
[
P
u
1
P
u
2
]
[
Σ
u
2
0
0
0
]
[
V
u
1
T
V
u
2
T
]
.
It is a property of the SVD that the product V
u
1
T
V
u
2
=
0. Thus, if we multiply (5) on the right by V
u
2 it yields
(11)
Y
V
u
2
=
Γ
X
V
u
2
.
Taking a SVD of Y
V
u
2 results in
(12)
Y
V
u
2
=
P
Σ
V
T
=
[
P
1
P
2
]
[
Σ
1
0
0
0
]
[
V
1
T
V
2
T
]
.
The order of the unknown system is determined by inspecting the singular values in Σ and the SVD is partitioned accordingly to obtain P
1, P
2, and Σ
1. Comparing (11) and (12), the matrix Γ is then identified as
(13)
Γ
=
P
1
with Γ given by (6); the matrix C is the first block row of the observability matrix. The state matrix A is obtained by
(14)
A
=
Γ
-
†
Γ
_
,
where † represents the pseudo inverse and computes a least squares solution. Γ
- and Γ
_ are obtained, respectively, by dropping the last and first block row of Γ. For the determination of B and D, (5) is multiplied on the right by the pseudo inverse of U and on the left by P
u
2
T. Taking into account that U
U
†
=
I and P
u
2
T
Γ
i
=
0, it can be deduced that
(15)
P
u
2
Y
p
U
†
=
P
u
2
H
.
Equation (15) can be solved for B and D using the least squares method. In this way, the system matrices of the unknown system are identified using the subspace identification method. Clearly, special attention should be paid to determining the order of the system from (12). Fortunately, there always exists a significant drop in singular values for an increasing position index that signals the order of the system.Band limited white noise signals in a frequency range varying from 0.5 Hz to 150 Hz, which cover 8 low-frequency modes of the plate, are applied to the 3 PZT actuators as input data and the corresponding voltage output data of the 3 PZT sensors are obtained by performing open loop FE analysis using the software package ANSYS. Transient simulations are performed for 6 s with a time step of 0.0025 s. The subspace identification method is then used to create a state space ROM. The order of the state space model is determined upon inspecting the singular values of the matrixY
V
u
2 mentioned in (12), as shown in Figure 2. For the purpose of comparison, two models are identified. The model of order 6 is called model 1 and the model of order 8 is called model 2. Table 2 lists the natural frequencies of the FE model and the two identified ROMs. It can be seen from Table 2 that model 1 represents the first 4 modes of the FE model very well while only 3 modes are included in model 2. Note that the 3rd mode of the FE model is not identified in model 2. Moreover, by checking the controllability matrices and observability matrices, it is found that both model 1 and model 2 are controllable and observable.Table 2
Natural frequencies of the piezoelectric plate (Hz).
FEM
Model 1
Model 2
Mode 1
2.178
2.180
2.182
Mode 2
13.041
13.089
13.089
Mode 3
20.392
20.427
35.795
Mode 4
35.421
35.787
Mode 5
61.081Figure 2
The distribution of singular values.


## 2.1. Finite Element Model

The constitutive equations of a deformable piezoelectric material in Voigt notation, coupling the elastic, and the electric fields can be expressed as(1)
σ
=
c
ε
-
e
T
E
,
D
=
e
ε
+
ϵ
ε
E
,
where σ and ε represent stress vector and strain vector, respectively. E is the electric field vector, D is the electric displacement vector, c is the elastic coefficient matrix, ϵ
ε is the dielectric coefficient matrix at constant ε, and e is the piezoelectric stress coefficient matrix.By applying element discretization and using the Hamilton’s principle, the global dynamic system equation for a piezoelectric smart structure is derived as(2)
[
M
u
u
0
0
0
]
[
q
¨
𝚽
¨
]
+
[
C
u
u
0
0
0
]
[
q
˙
𝚽
˙
]
+
[
K
u
u
K
u
ϕ
K
u
ϕ
T
K
ϕ
ϕ
]
[
q
𝚽
]
=
{
F
u
F
ϕ
}
,
where q denotes the global nodal displacement vector, 𝚽 denotes electric potential vector, and a dot above a variable denotes a time derivative. F
u and F
ϕ denote the nodal external force vector and the nodal externally applied charge vector, respectively. M
u
u and K
u
u are the structural mass and stiffness matrix, respectively, and K
u
ϕ is the piezoelectric coupling matrices; K
ϕ
ϕ is the dielectric permittivity matrix and C
u
u is the damping matrix, which is usually assumed to be constructed as the weighted sum of mass and stiffness matrices
(3)
C
u
u
=
α
M
u
u
+
β
K
u
u
,
where α and β are Rayleigh’s coefficients which are nonnegative and small in the sense that the effective damping coefficient is ≪1.A rectangular aluminum cantilever plate which consists of a host aluminum plate (700 mm × 150 mm × 1.2 mm) and 3 pairs of PZT patches (70 mm × 15 mm × 0.5 mm), as shown in Figure1, is considered. The 3 PZT patches bonded to the top surface are to be used as actuators and 3 other PZT patches bonded symmetrically to the bottom surface are to be used as sensors. The material properties of the host plate and PZT patches are given in Table 1. ANSYS software package is used to model the cantilever plate with PZT patches. In this study, the 8-node SOLID5 coupled field element and 8-node SOLID45 element are used to model the PZT patches and the host structure, respectively. The cantilever boundary conditions are simulated by specifying zero displacement DOF for the x, y, and z directions at the fixed end of the plate. The electrical DOFs of the nodes on the outside surface of a PZT patch are coupled by the CP command so that the nodes on the outside surface represent an electrode. The electrode on the inner surface of a PZT patch is defined in the same manner and a fixed zero voltage potential is defined on the inner surface of the PZT patch to model the ground electrode. A total of 17248 three-dimensional solid elements are used to mesh the structure. Among those elements, 1008 SOLID5 elements are utilized to model piezoelectric patches, and 16240 SOLID45 elements are employed to model host structure. Rayleigh type damping is considered with coefficients α
=
0.0 rad s−1 and β
=
0.003 s.Table 1
Material properties.
PZT
Host plate
ρ (kg/m3)
7600.0
2700
E
11
=
E
22  (GPa)
63.0
69.0
E
33 (GPa)
54.0
69.0
G
12 (GPa)
24.6
27.0
G
23
=
G
13  (GPa)
30.6
27.0
ν
12
0.3
0.32
ν
13
=
ν
23
0.4
0.32
e
31, e
32, e
33, e
15 (C/m2)
−5.4, −5.4, 15.8, 12.3
/
ϵ
11
ε, ϵ
22
ε, ϵ
33
ε (×10−8 F/m)
1.53, 1.53, 1.50
/Figure 1
Configuration of a cantilever plate with PZT S/A pairs.


## 2.2. The Reduced Representation of the System

FE model is appropriate for dynamic analysis but cannot be used directly as a basis for controller design. Considering that FE simulations are analogous to performing experimental investigations where the only direct outputs are the time histories, a so-called subspace identification method is employed to develop a state space ROM for the piezoelectric smart structure. In this way, the multiple-input multiple-output state space model is identified, which is convenient for the following controller design procedure.Suppose that the FE model of a piezoelectric smart structure can be represented by the following state space model:(4)
x
k
+
1
=
A
x
k
+
B
u
k
,
y
k
=
C
x
k
+
D
u
k
.
Here, x
k
∈
ℝ
n is the state vector at discrete time instant k and contains the numerical values of n states. u
k
∈
ℝ
m and y
k
∈
ℝ
l are the measurements at time instant k of, respectively, the m inputs and l outputs. As to the piezoelectric smart structure mentioned above, the inputs are the externally exerted mechanical force and the applied voltages over the actuators and the outputs are the sensor voltages. A
∈
ℝ
n
×
n, B
∈
ℝ
n
×
m, C
∈
ℝ
l
×
n, and D
∈
ℝ
l
×
m are state matrix, input matrix, output matrix, and direct influence matrix, respectively. It is assumed that the matrix pair (
A
,
C
) is observable and (
A
,
B
) is controllable. The objective of the subspace identification method is to determine the order of the unknown system and the system matrices A, B, C, and D from given input-output data of the FE model.The input-output relations of the linear control system described by (4) can be expressed compactly in a matrix form as [32, 33]
(5)
Y
=
Γ
X
+
H
U
,
where the extended observability matrix Γ is defined as
(6)
Γ
=
[
C
C
A
⋯
C
A
i
-
1
]
∈
ℝ
l
i
×
n
.
H is the lower block triangular Toeplitz matrix of impulse response from u to y given by
(7)
H
=
[
D
0
0
⋯
0
C
B
D
0
⋯
0
C
A
B
C
B
D
⋯
0
⋯
⋯
⋯
⋯
⋯
C
A
i
-
2
B
C
A
i
-
3
B
C
A
i
-
4
B
⋯
D
]
∈
ℝ
l
i
×
m
i
.
U and Y are, respectively, the input and output block Hankel matrices, which are defined as
(8)
U
=
[
u
0
u
1
u
2
⋯
u
j
-
1
u
1
u
2
u
3
⋯
u
j
⋯
⋯
⋯
⋯
⋯
u
i
-
1
u
i
u
i
+
1
⋯
u
i
+
j
-
2
]
∈
ℝ
m
i
×
j
,
Y
=
[
y
0
y
1
y
2
⋯
y
j
-
1
y
1
y
2
y
3
⋯
y
j
⋯
⋯
v
⋯
⋯
y
i
-
1
y
i
y
i
+
1
⋯
y
i
+
j
-
2
]
∈
ℝ
l
i
×
j
.
The state sequence X is defined as
(9)
X
=
[
x
0
x
1
⋯
x
j
-
1
]
∈
ℝ
n
×
j
.
The number of block rows (i) is a user specified parameter which is required to be larger than the maximum order of the system one wants to identify. The number of columns (j) is assumed to approach ∞.The key idea of subspace identification method is to estimate the extended observability matrixΓ by computing the singular value decomposition (SVD) of U in the first step:
(10)
U
=
[
P
u
1
P
u
2
]
[
Σ
u
2
0
0
0
]
[
V
u
1
T
V
u
2
T
]
.
It is a property of the SVD that the product V
u
1
T
V
u
2
=
0. Thus, if we multiply (5) on the right by V
u
2 it yields
(11)
Y
V
u
2
=
Γ
X
V
u
2
.
Taking a SVD of Y
V
u
2 results in
(12)
Y
V
u
2
=
P
Σ
V
T
=
[
P
1
P
2
]
[
Σ
1
0
0
0
]
[
V
1
T
V
2
T
]
.
The order of the unknown system is determined by inspecting the singular values in Σ and the SVD is partitioned accordingly to obtain P
1, P
2, and Σ
1. Comparing (11) and (12), the matrix Γ is then identified as
(13)
Γ
=
P
1
with Γ given by (6); the matrix C is the first block row of the observability matrix. The state matrix A is obtained by
(14)
A
=
Γ
-
†
Γ
_
,
where † represents the pseudo inverse and computes a least squares solution. Γ
- and Γ
_ are obtained, respectively, by dropping the last and first block row of Γ. For the determination of B and D, (5) is multiplied on the right by the pseudo inverse of U and on the left by P
u
2
T. Taking into account that U
U
†
=
I and P
u
2
T
Γ
i
=
0, it can be deduced that
(15)
P
u
2
Y
p
U
†
=
P
u
2
H
.
Equation (15) can be solved for B and D using the least squares method. In this way, the system matrices of the unknown system are identified using the subspace identification method. Clearly, special attention should be paid to determining the order of the system from (12). Fortunately, there always exists a significant drop in singular values for an increasing position index that signals the order of the system.Band limited white noise signals in a frequency range varying from 0.5 Hz to 150 Hz, which cover 8 low-frequency modes of the plate, are applied to the 3 PZT actuators as input data and the corresponding voltage output data of the 3 PZT sensors are obtained by performing open loop FE analysis using the software package ANSYS. Transient simulations are performed for 6 s with a time step of 0.0025 s. The subspace identification method is then used to create a state space ROM. The order of the state space model is determined upon inspecting the singular values of the matrixY
V
u
2 mentioned in (12), as shown in Figure 2. For the purpose of comparison, two models are identified. The model of order 6 is called model 1 and the model of order 8 is called model 2. Table 2 lists the natural frequencies of the FE model and the two identified ROMs. It can be seen from Table 2 that model 1 represents the first 4 modes of the FE model very well while only 3 modes are included in model 2. Note that the 3rd mode of the FE model is not identified in model 2. Moreover, by checking the controllability matrices and observability matrices, it is found that both model 1 and model 2 are controllable and observable.Table 2
Natural frequencies of the piezoelectric plate (Hz).
FEM
Model 1
Model 2
Mode 1
2.178
2.180
2.182
Mode 2
13.041
13.089
13.089
Mode 3
20.392
20.427
35.795
Mode 4
35.421
35.787
Mode 5
61.081Figure 2
The distribution of singular values.


## 3. Active Vibration Control

### 3.1. ROM Based Controller

The LQR based optimal controller is designed for reducing the vibration of the piezoelectric plate. A state feedback is adopted to minimize the cost function such that the requisite design criteria are achieved. The cost function given in a quadratic form is defined as(16)
J
=
E
{
∑
k
=
1
∞
[
x
k
T
Q
x
k
+
u
k
T
R
u
k
]
}
,
where Q is a positive semidefinite state weighting matrix and R is a symmetric positive-definite input weighting matrix. The input voltages at the PZT actuators are defined as the multiple of the feedback gain matrix and the system state values as follows:
(17)
u
k
=
-
G
1
x
k
in which G
1 denotes the state feedback gain matrix and is given by
(18)
G
1
=
(
R
+
B
T
P
B
)
-
1
B
T
P
A
,
where the positive definite matrix P is the solution to the following algebraic Riccati equation:
(19)
A
T
P
A
-
P
-
A
T
P
B
(
R
+
B
T
P
B
)
-
1
B
T
P
A
+
Q
=
0
.The LQR control law expressed by (17) demands knowledge of the full state vector x
k. However, the state vector is seldom available in practice. To incorporate the states information in the control system, it is necessary to estimate the states of the system from the sensor outputs. An estimated state x
^
k of the actual state vector x
k can be obtained based on Luenberger observer, which introduces in (4) an additional term of the form gain multiplied by the output error
(20)
x
^
k
+
1
=
A
x
^
k
+
B
u
k
+
G
2
[
y
k
-
C
x
^
k
-
D
u
k
]
,
where G
2 denotes the observer gain matrix which can be chosen based on the observer pole assignment method such that the matrix A
-
G
2
C is asymptotically stable. The actual state x
k in the feedback control law given by (17) is now replaced by the estimated state
(21)
u
k
=
-
G
1
x
^
k
.In general, higher values forQ correspond to higher priority for vibration suppression, whereas large values of R mean higher interest in reducing the control effort. An efficient way of choosing Q is to compute it from C with
(22)
Q
=
C
T
C
.
Weighting matrix R can be set as λ
I with λ as a scalar design parameter. Therefore, the challenging task of choosing Q and R reduces to choosing one parameter λ. In this study, λ is selected as λ
=
0.01 to most effectively control the structure.


### 3.2. Closed Loop Simulation in ANSYS

In this section, we begin our discussion on implementing the LQR controller coupled with a Luenberger state estimator into ANSYS. The block diagram of the analysis is shown in Figure3. Here the actuators pump energy into the FE model and unmodeled dynamics, which are not included in ROM, can be excited. The sensor outputs, which include unmodeled dynamics, are then processed by the controller to provide control feedback to piezoelectric actuators for active vibration control. Consequently, the high order system (FE model) is driven by the ROM based controller. In this way, the closed loop system shown in Figure 3 can adequately represent the behaviour of the actual closed loop system and the issues of spillover can be addressed through numerical simulations. However, in some traditional vibration control methods, the ROM based controller outputs are fed back to the ROM, but not the FE model, so that the spillover effects cannot be captured in this way.Figure 3
Block diagram of the closed loop simulation in ANSYS.To implement the controller into ANSYS, a macro that involves utilizing the commands for elementary matrix operations and basic program flow controls has to be developed by using APDL. Supposing that the ROM expressed by (4) is available and the state feedback gain matrix G
1 together with the observer gain matrix G
2 have been calculated, the *DIM command and *SET command are used to store all these variables in a matrix form. The Luenberger observer and the control law are implemented according to (20) and (21) by using the elementary matrix manipulation commands in APDL, for example, *MOPER for multiplication, *VFUN for copying, and *VOPER for adding. The main control algorithm is realized by the *DO and *ENDDO commands in which the transient analysis is performed with a time step of 0.0025 s and the actuator voltages are updated at each time instant. To explain in more detail, given that the instantaneous values of the sensor signals are calculated at a time instant under the excitation of actuator voltages and the vibration generating forces, the states of the system at current time instant are estimated according to (20). Then, the feedback voltages, which will be used as the inputs to the PZT actuators during the next time instant, are determined by the state feedback control law according to (21). After that, the sensor signals at the next time instant are calculated under the updated actuator voltages. This process will continue for the predefined time duration of the computed response.


## 3.1. ROM Based Controller

The LQR based optimal controller is designed for reducing the vibration of the piezoelectric plate. A state feedback is adopted to minimize the cost function such that the requisite design criteria are achieved. The cost function given in a quadratic form is defined as(16)
J
=
E
{
∑
k
=
1
∞
[
x
k
T
Q
x
k
+
u
k
T
R
u
k
]
}
,
where Q is a positive semidefinite state weighting matrix and R is a symmetric positive-definite input weighting matrix. The input voltages at the PZT actuators are defined as the multiple of the feedback gain matrix and the system state values as follows:
(17)
u
k
=
-
G
1
x
k
in which G
1 denotes the state feedback gain matrix and is given by
(18)
G
1
=
(
R
+
B
T
P
B
)
-
1
B
T
P
A
,
where the positive definite matrix P is the solution to the following algebraic Riccati equation:
(19)
A
T
P
A
-
P
-
A
T
P
B
(
R
+
B
T
P
B
)
-
1
B
T
P
A
+
Q
=
0
.The LQR control law expressed by (17) demands knowledge of the full state vector x
k. However, the state vector is seldom available in practice. To incorporate the states information in the control system, it is necessary to estimate the states of the system from the sensor outputs. An estimated state x
^
k of the actual state vector x
k can be obtained based on Luenberger observer, which introduces in (4) an additional term of the form gain multiplied by the output error
(20)
x
^
k
+
1
=
A
x
^
k
+
B
u
k
+
G
2
[
y
k
-
C
x
^
k
-
D
u
k
]
,
where G
2 denotes the observer gain matrix which can be chosen based on the observer pole assignment method such that the matrix A
-
G
2
C is asymptotically stable. The actual state x
k in the feedback control law given by (17) is now replaced by the estimated state
(21)
u
k
=
-
G
1
x
^
k
.In general, higher values forQ correspond to higher priority for vibration suppression, whereas large values of R mean higher interest in reducing the control effort. An efficient way of choosing Q is to compute it from C with
(22)
Q
=
C
T
C
.
Weighting matrix R can be set as λ
I with λ as a scalar design parameter. Therefore, the challenging task of choosing Q and R reduces to choosing one parameter λ. In this study, λ is selected as λ
=
0.01 to most effectively control the structure.


## 3.2. Closed Loop Simulation in ANSYS

In this section, we begin our discussion on implementing the LQR controller coupled with a Luenberger state estimator into ANSYS. The block diagram of the analysis is shown in Figure3. Here the actuators pump energy into the FE model and unmodeled dynamics, which are not included in ROM, can be excited. The sensor outputs, which include unmodeled dynamics, are then processed by the controller to provide control feedback to piezoelectric actuators for active vibration control. Consequently, the high order system (FE model) is driven by the ROM based controller. In this way, the closed loop system shown in Figure 3 can adequately represent the behaviour of the actual closed loop system and the issues of spillover can be addressed through numerical simulations. However, in some traditional vibration control methods, the ROM based controller outputs are fed back to the ROM, but not the FE model, so that the spillover effects cannot be captured in this way.Figure 3
Block diagram of the closed loop simulation in ANSYS.To implement the controller into ANSYS, a macro that involves utilizing the commands for elementary matrix operations and basic program flow controls has to be developed by using APDL. Supposing that the ROM expressed by (4) is available and the state feedback gain matrix G
1 together with the observer gain matrix G
2 have been calculated, the *DIM command and *SET command are used to store all these variables in a matrix form. The Luenberger observer and the control law are implemented according to (20) and (21) by using the elementary matrix manipulation commands in APDL, for example, *MOPER for multiplication, *VFUN for copying, and *VOPER for adding. The main control algorithm is realized by the *DO and *ENDDO commands in which the transient analysis is performed with a time step of 0.0025 s and the actuator voltages are updated at each time instant. To explain in more detail, given that the instantaneous values of the sensor signals are calculated at a time instant under the excitation of actuator voltages and the vibration generating forces, the states of the system at current time instant are estimated according to (20). Then, the feedback voltages, which will be used as the inputs to the PZT actuators during the next time instant, are determined by the state feedback control law according to (21). After that, the sensor signals at the next time instant are calculated under the updated actuator voltages. This process will continue for the predefined time duration of the computed response.


## 4. Numerical Example and Discussion

The spillover effects on the controller performance will be investigated through numerical simulations. For all the cases considered following, the plate is subjected to time dependent mechanical loads applied at the free end and the tip displacement responses are used to illustrate the spillover effects. For the purpose of comparison, ROM based controllers of order 8 and 6, respectively, are integrated into the same FE model to perform the closed loop simulations.


### 4.1. Spillover Effects under Sinusoidal Loading

The plate is subjected to a harmonic load0.015
sin
⁡
(
2
π
f
t
+
π
/
6
) N with frequency f
=
2.2 Hz, which is close to the first natural frequency. Figure 4 shows the controlled and uncontrolled tip displacements in the frequency domain. In Figure 4 the three axes are frequency, load cases, and power spectral density (PSD), respectively, and each frequency response is obtained from the time domain data during 15 s. As can be seen from Figure 4, both the controller of order 6 and order 8 reduce the peak value of the first mode. However, compared with the results using the controller of order 8, the effect of the controller of order 6 on the active vibration control is not so significant indicating that the performance of the closed loop system is very much influenced by the order of the controller.Figure 4
Tip displacement PSD under a sinusoidal load of 2.2 Hz.Next, a harmonic load0.3
sin
⁡
(
2
π
f
t
+
π
/
3
) N with frequency f
=
20.1 Hz, which is close to the third natural frequency, is applied to the plate at the free end. Note that the ROM of order 6 contains no information of the third mode. Figure 5 shows the controlled and uncontrolled PSDs. As shown in this figure, the first and the third modes are excited and then are damped out effectively by the controller of order 8. However, the controller of order 6 cannot suppress the third mode because the controller of order 6 has no information of the third order. Moreover, the controller of order 6 pumps energy into the third mode so the peak value of the third mode is much larger than the open loop simulation. This phenomenon is a typical example of spillover and the spillover effects significantly degrade the controller performance.Figure 5
Tip displacement PSD under a sinusoidal load of 20.1 Hz.


### 4.2. Spillover Effects under Multisinusoidal Loading

A final investigation of spillover effects on the controller performance is conducted by applying a mechanical load of the form(23)
0.01
sin
(
2
π
f
1
t
+
π
4
)
+
0.01
sin
(
2
π
f
2
t
+
2
π
3
)
+
0.2
sin
(
2
π
f
3
t
-
π
3
)
with f
1
=
2.0 Hz, f
2
=
13.1 Hz, and f
3
=
20.34 Hz. Figure 6 shows the controlled and uncontrolled results. As expected, the controller of order 6 reduces the peak values of the first two modes while the third mode cannot be damped out by this controller. Once again, spillover effects for the third mode are observed when using the controller of order 6. However, the first three modes are effectively suppressed when the controller of order 8 is adopted, which indicates that a potential method to prevent the spillover problem is to use a higher order controller.Figure 6
Tip displacement PSD under a multisinusoidal load.


## 4.1. Spillover Effects under Sinusoidal Loading

The plate is subjected to a harmonic load0.015
sin
⁡
(
2
π
f
t
+
π
/
6
) N with frequency f
=
2.2 Hz, which is close to the first natural frequency. Figure 4 shows the controlled and uncontrolled tip displacements in the frequency domain. In Figure 4 the three axes are frequency, load cases, and power spectral density (PSD), respectively, and each frequency response is obtained from the time domain data during 15 s. As can be seen from Figure 4, both the controller of order 6 and order 8 reduce the peak value of the first mode. However, compared with the results using the controller of order 8, the effect of the controller of order 6 on the active vibration control is not so significant indicating that the performance of the closed loop system is very much influenced by the order of the controller.Figure 4
Tip displacement PSD under a sinusoidal load of 2.2 Hz.Next, a harmonic load0.3
sin
⁡
(
2
π
f
t
+
π
/
3
) N with frequency f
=
20.1 Hz, which is close to the third natural frequency, is applied to the plate at the free end. Note that the ROM of order 6 contains no information of the third mode. Figure 5 shows the controlled and uncontrolled PSDs. As shown in this figure, the first and the third modes are excited and then are damped out effectively by the controller of order 8. However, the controller of order 6 cannot suppress the third mode because the controller of order 6 has no information of the third order. Moreover, the controller of order 6 pumps energy into the third mode so the peak value of the third mode is much larger than the open loop simulation. This phenomenon is a typical example of spillover and the spillover effects significantly degrade the controller performance.Figure 5
Tip displacement PSD under a sinusoidal load of 20.1 Hz.


## 4.2. Spillover Effects under Multisinusoidal Loading

A final investigation of spillover effects on the controller performance is conducted by applying a mechanical load of the form(23)
0.01
sin
(
2
π
f
1
t
+
π
4
)
+
0.01
sin
(
2
π
f
2
t
+
2
π
3
)
+
0.2
sin
(
2
π
f
3
t
-
π
3
)
with f
1
=
2.0 Hz, f
2
=
13.1 Hz, and f
3
=
20.34 Hz. Figure 6 shows the controlled and uncontrolled results. As expected, the controller of order 6 reduces the peak values of the first two modes while the third mode cannot be damped out by this controller. Once again, spillover effects for the third mode are observed when using the controller of order 6. However, the first three modes are effectively suppressed when the controller of order 8 is adopted, which indicates that a potential method to prevent the spillover problem is to use a higher order controller.Figure 6
Tip displacement PSD under a multisinusoidal load.


## 5. Conclusions

In this paper, a preliminary research is presented for investigation of the spillover effects on the vibration control of piezoelectric smart structure in FE environment, and therefore an efficient method especially beneficial for preliminary design of piezoelectric smart structure is provided to evaluate the performance of candidate control laws in FE environment considering spillover effects. The basic idea of this paper is that the FE model is more accurate than ROM for modeling of physical systems such that the unmodeled dynamics can be included in FE model. By using APDL, a ROM based controller is integrated into the ANSYS software package to provide an accurate representation of what will happen when the controller is connected to the real plant. In this way, the actuators pump energy into the FE model and unmodeled dynamics, which are not included in ROM, and can be excited so that the issues of spillover effects can be addressed in the closed loop FE environment. The ROM based controller is designed to tailor the output of the ROM modes to meet performance specifications; however, there is no prior guarantee that a ROM based controller satisfying the control requirements still works well in closed loop with the FE model because of spillover effects. The spillover effects are illustrated via the PSD of the tip displacement response. It is demonstrated that by using a higher order controller the spillover problem can be prevented. Future research calls for investigation of spillover induced instability in vibration control of smart structures using the proposed method.



---

*Source: 870940-2014-04-02.xml*