國 立 交 通 大 學
電機與控制工程學系
碩
士
論 文
使用多項式近似法實現
Preisach
模型
並應用於微壓電致動器
Realization of Preisach Model Using Polynomial Approximation
with Applications for Micro Piezoelectric Actuators
研 究 生
:
溫宏揚
指導教授
:
鄭木火 博士
使用多項式近似法實現
Preisach
模型並應用於
微壓電致動器
研究生
:
溫宏揚
指導教授
:
鄭木火博士
國立交通大學電機與控制學系
中文摘要
近年來由於工程應用上對精準度要求的增加
,磁滯非線性系統的建模已受到相當的注意。 在現今
已發展的數種描述磁滯系統的模型
, Preisach模型是最廣為人知的一種。 磁滯行為的描述可以被
Preisach
模型的無限多組可計算的一階逆曲線
(First Order Reversal Curves, FORC)表示。 實現
方面
,傳統的方法是以資料表格實現
Preisach模型
;其資料相應
FORC的有限數量的取樣點
,而透
過表格資料的線性內插計算模型輸出。 然而
,這種方法有兩個缺點
:為了準確地預測磁滯行為需要
大量的記憶體
,以及發展有效的方法調整表格資料以反映磁滯元件的老化或時間效應是困難的
,或
說甚至是不可能的。 為了克服這些缺點
,本論文提出使用一組多項式取代資料表格來實現
Preisach模型。 這樣的方式明顯地減少所需的記憶體因為只需儲存多項式的少量係數。 再者
,多項式的係數
可以使用最小平方法或是適應性鑑別演算法取得
,因此可以追蹤磁滯系統的參數。 我們提出的方
法已經被驗證在微壓電致動器的位移預測與追跡控制上。 我們應用最小平均平方的方法發展一種
適應性的演算法來鑑別微壓電致動器的多項式係數
;以此結果實現並與表格實現的方式相比較
,在
模型準確度與減少所需記憶體大小方面有顯著的改善。
關鍵詞
:磁滯
, Preisach模型
,壓電致動器
,多項式
Realization of Preisach Model Using Polynomial
Approximation with Applications for Micro
Piezoelectric Actuators
Student: Home-Young Win
Advisor: Dr. Mu-Huo Cheng
Institute of Electrical and Control Engineering
National Chiao-Tung University
Abstract
Modeling of systems with hysteresis nonlinearities has received considerable attention re-cently due to the increasing accuracy requirement in engineering applications. Several models have been developed to characterize systems with hysteresis, and the most popular one is the Preisach model. The Preisach model to characterize the hysteresis behavior can be represented by infinite but countable first order reversal curves (FORC). In practice, the conventional approach to realize the Preisach model uses a data table; the data of which correspond with the samples of a finite number of FORC, then the model output is evaluated via the linear interpolation from the table data. This approach, however, suffers from two drawbacks: it requires a large amount of memory in order to obtain an accurate prediction of hysteresis behavior and it is difficult or even impossible to derive efficient ways to modify the data table in order to reflect the aging or timing effect of elements with hysteresis. To overcome these drawbacks, this thesis proposes to use a set of polynomials instead of a data table for realization of the Preisach model. The proposed approach reduces significantly the required memory because it only requires to store a small num-ber of polynomial coefficients. Furthermore, the polynomial coefficients can be obtained using the least-square approximation or the adaptive identification algorithm such that the tracking of hysteresis model parameters is possible. The proposed approach has been verified by the displace-ment prediction and the tracking control of the micro piezoelectric actuator. We apply the least mean square method to develop an adaptive algorithm for identification of the polynomial coeffi-cients for the micro piezoelectric actuator; the resulting realization compared with the realization via table method yields significant improvement in model accuracy as well as the reduction in the required memory size.
Contents
ABSTRACT IN CHINESE I
ABSTRACT IN ENGLISH II
1 Introduction 1
1.1 Introduction of Hysteresis . . . 1
1.2 Organization of this Thesis . . . 2
2 Preisach Model Representation 3 2.1 Typical Properties of Hysteresis . . . 3
2.1.1 Hysteresis Loop . . . 3
2.1.2 Static Hysteresis Nonlinearity . . . 5
2.2 Preisach Model . . . 6
2.2.1 Hysteresis Operator . . . 6
2.3 Preisach Plane . . . 7
2.3.1 Geometric Interpretation of the Relay Operator . . . 8
2.3.2 Wiping-Out Property . . . 9
2.4 First Order Reversal Curve . . . 11
2.4.1 Experimental Determination of Weighting Function . . . 13
2.4.2 Numerical Realization of the Preisach Model . . . 13
3 Preisach Model Realization and Inverse Preisach Model 16 3.1 Table Method . . . 17
3.2 Polynomial Approximation . . . 18
3.2.1 Polynomial Curve Fitting . . . 18
3.2.2 Polynomial Approximation . . . 19
3.3 Adaptive Polynomial Approximation . . . 21
3.3.1 Adaptive Polynomial Approximation with LMS . . . 22
3.4.1 Cubic Equation . . . 23
3.4.2 Solving Polynomial of the Preisach Model . . . 25
4 Applications for Piezoelectric Actuators 26 4.1 Experimental Setup . . . 26
4.2 Table Method . . . 27
4.3 Polynomial Approximation . . . 28
4.4 Polynomial Approximation with Adaptive Identification . . . 31
4.5 Tracking Control . . . 35
5 Conclusions 39
List of Tables
4.1 Sampling data of FORC. . . 28
4.2 The coefficients of polynomial. . . 30
List of Figures
2.1 Hysteresis behavior: (a) the input of a hysteresis system; (b) the output of a
hys-teresis system. . . 4
2.2 Phase transition of a hysteresis system. . . 4
2.3 Hysteresis behavior: (a) the input of a hysteresis system; (b) the phase transition of a hysteresis system. . . 5
2.4 Hysteresis behavior: (a) the input of a hysteresis system; (b) the input of a hystere-sis system; (c) the phase transition of a hysterehystere-sis system. . . 5
2.5 Block diagram interpretation of the Preisach model. . . 6
2.6 The hysteresis operator: (a) the phase diagram; (b) the operator output when the inputu(t) is increasing; (c) the operator output when the input u(t) is decreasing. . 7
2.7 The Preisach plane. . . 8
2.8 The interpretation of relay operator via Preisach plane: (a) when the input increases from zero tou1; (b) when the input changes direction and decreases fromu1 tou2. 8 2.9 A changing input and its memory formation mechanism of the Preisach plane. . . . 9
2.10 Wiping out property: (a)M2 < M1(b)M2 = M1 (c)M2 > M1. . . 10
2.11 Alternating series of the input extremum. . . 11
2.12 The change on the Preisach plane to yield a FORC: (a) the input increases from 0 to some valueM; (b) then the input decreases from M to some value m. . . 12
2.13 First order reversal on the phase diagram:(a) an ascending branch; (b) an descend-ing branch. . . 12
2.14 Corresponding change on Preisach plane of an input:(a) from 0 to M1; (b) from M1 tom1; (c) fromm1 tou(t); (d) from M2tou(t). . . 15
3.1 Linear interpolation betweena and b. . . 17
3.2 Approximatefαβ(M, m) where (M,m) belongs to a rectangular cell. . . 17
3.3 Approximatefαβ(M, m) where (M,m) belongs to a triangular cell. . . 18
3.5 Block diagram of adaptive system identification. . . 22
4.1 Real system for experiment. . . 27
4.2 The block diagram of the entire experimental system. . . 27
4.3 The input and output of the micro piezoelectric actuator. . . 28
4.4 The phase transition of micro piezoelectric actuator. . . 29
4.5 Prediction of micro piezoelectric actuator using Preisach model realization via ta-ble method. . . 29
4.6 Phase transition of micro piezoelectric actuator using Preisach model realization via table method. . . 30
4.7 Prediction of micro piezoelectric actuator using Preisach model realization via polynomial approximation. . . 31
4.8 Phase transition of micro piezoelectric actuator using Preisach model realization via polynomial approximation. . . 31
4.9 Mean square error of the prediction of micro piezoelectric actuator for varying step-sizeµ using adaptive polynomial approximation. . . 32
4.10 Prediction of micro piezoelectric actuator using Preisach model realization via adaptive polynomial approximation. . . 33
4.11 Phase transition of micro piezoelectric actuator using Preisach model realization via adaptive polynomial approximation. . . 34
4.12 Phase transition of the model error and the input. . . 34
4.13 Prediction of micro piezoelectric actuator using Preisach model realization via modified polynomial approximation. . . 35
4.14 Phase transition of micro piezoelectric actuator using Preisach model realization via modified polynomial approximation. . . 36
4.15 Model error of the Preisach model realization via the table method and the modi-fied polynomial approximation. . . 36
4.16 Block diagram of the tracking control system using inverse Preisach model. . . 37
4.17 Tracking of micro piezoelectric actuator using inverse Preisach model. . . 37
4.18 Block diagram of the tracking control system using inverse Preisach model with PID controller. . . 37
4.19 Tracking of micro piezoelectric actuator using inverse Preisach model with PID controller. . . 38
Chapter 1
Introduction
1.1
Introduction of Hysteresis
Hysteresis is a nonlinear phenomenon occurring in many engineering devices such as the micro piezoelectric actuators (PZA), the shape memory alloys (SMA), and the ferromagnetic ele-ments [1, 3]. A system with hysteresis is usually difficult to describe accurately and may result in unstable behaviors if not controlled appropriately [3]. Therefore, an accurate model of hysteresis is critical in order to develop suitable control algorithms or applications with these systems [4]. The Preisach model has been known as the most popular one to characterize hysteresis behaviors [2, 7]. The model uses double integrals with a relay operator, and weighting parameters to describe the system input/output relation; this relation can be further represented by an infinite but countable first order reversal curves (FORC). The conventional approach to realize the Preisach model for a hysteresis system approximates the representation by the sampling of a finite number of FORC and stores the sample data in a table [10, 13]. The model output is then evaluated from the table data and the given input by simple linear interpolation. This approach, however, suffers from two drawbacks: One is that it requires a large amount of memory in order to obtain an accurate pre-diction of hysteresis behavior because the number of FORC and the sampling resolution should be increased. The other is that it is difficult or even impossible to derive efficient ways to adjust the data table to reflect the parameter variations with aging or timing effect. Hence, the realization accuracy may degrade with time and eventually needs to rebuild the table.
To overcome these drawbacks, in this thesis we propose to use a set of polynomials instead of a data table for realization of the Preisach model. The proposed approach reduces the required mem-ory because it only requires to store a small number of polynomial coefficients. Furthermore, the polynomial coefficients can be obtained using the least-square approximation [11] or the adaptive identification algorithm [14] such that the tracking of hysteresis model parameters is possible. The
proposed approach has been applied to model the displacement of a micro piezoelectric actuator. We apply the least mean square (LMS) adaptive algorithm to develop an adaptive algorithm for identification of the polynomial coefficients of the micro piezoelectric actuator; the resulting re-alization compared with the rere-alization via table method yields significant improvement in model accuracy as well as the reduction in the required memory size.
1.2
Organization of this Thesis
The remainder of this thesis is divided to four chapters including conclusions. Chapter 2 reviews the hysteresis and the Preisach model. Chapter 3 illustrates the difference between the traditional method and our method to realize the Preisach model. Chapter 4 demonstrates the applications of the micro piezoelectric actuators. The final chapter is the conclusions.
Chapter 2
Preisach Model Representation
In this chapter, the Preisach model to characterize the hysteresis is discussed with the focus on the representation of the model using FORC and its numerical implementation via an infinite number of FORC [2, 7].
2.1
Typical Properties of Hysteresis
Before introducing the Preisach model, we first review some typical hysteresis properties through few simple experiments [9, 10]. One can see that a system with hysteresis is a nonlin-ear system whose output depends not only on the current input but the past input [7]. Due to these facts, it is difficult to express the hysteresis via a simple mathematical form.
2.1.1
Hysteresis Loop
We begin to introduce hysteresis sketchily by a simple illustration. Fig. 2.1(a) shows an input of a system with hysteresis which increases to one value and then decreases to the initial value; Fig. 2.1(b) is the system response output and the system input/output relation, defined as the phase transition, plotted on Fig. 2.2. As the input increasing, the ascending branch comes up in the phase transition shown in Fig. 2.2; whereas the decreasing input results in descending branch. A hysteresis loop is formed because of the ascending branch disagrees with the descending branch. Hysteresis loops appear in the phase transition is one of the major characteristics of a system with hysteresis; it is obvious that a system with hysteresis is a nonlinear system according to this characteristic. Another major characteristic of a system with hysteresis is that the output depends on the past input. For instance, Fig. 2.3(b) has different descending branches corresponding to the same current input but the past input is different shown in Fig. 2.3(a).
40 45 50 55 60 0 0.5 1 1.5 2 2.5 3 3.5 4 4.5 time (sec) input signal 40 45 50 55 60 0 5 10 15 20 25 30 35 40 45 time (sec) output (a) (b)
Figure 2.1: Hysteresis behavior: (a) the input of a hysteresis system; (b) the output of a hysteresis system. 0 0.5 1 1.5 2 2.5 3 3.5 4 4.5 0 5 10 15 20 25 30 35 40 45 input signal output ascending branch descending branch
40 45 50 55 60 65 70 75 0 0.5 1 1.5 2 2.5 3 3.5 4 4.5 time (sec) input signal 0 1 2 3 4 0 5 10 15 20 25 30 35 40 45 input signal output (a) (b)
Figure 2.3: Hysteresis behavior: (a) the input of a hysteresis system; (b) the phase transition of a hysteresis system. 30 35 40 45 0 0.5 1 1.5 2 2.5 3 3.5 4 4.5 time (sec) input signal of case 1
30 40 50 60 70 80 90 0 0.5 1 1.5 2 2.5 3 3.5 4 4.5 time (sec) input signal of case 2
0 1 2 3 4 0 5 10 15 20 25 30 35 40 45 input output phase transition case 1 case 2 (a) (b) (c)
Figure 2.4: Hysteresis behavior: (a) the input of a hysteresis system; (b) the input of a hysteresis system; (c) the phase transition of a hysteresis system.
2.1.2
Static Hysteresis Nonlinearity
Static hysteresis nonlinearity is also an important quality of standard hysteresis systems [7]. What the term “static” means is that the output depends only on the past extremum values of input. In another word, how fast does the input change doesn’t influence on the phase transition. Fig. 2.4 is an experiment which is used to illustrate this phenomenon. Fig. 2.4(a) and Fig. 2.4(b) are different inputs but they has the same extremum values. Both phase transitions are plotted on Fig. 2.4(c), and one can find that their phase transition are almost overlapping.
2.2
Preisach Model
The Preisach model represents the input and output relations of a system with hysteresis by a double integral formula as below
f (t) = Z Z
α≥β
µ(α, β)γαβ[u(·)](t)dαdβ (2.1)
whereu(t) denotes the input signal, f (t) is the system response output, µ(α, β) is the weighting
function, andγαβ[u(·)](t) is an operator with two distinct output values, as discussed later below.
The ideal of the Preisach model, therefore, is to regard a hysteresis transducer as the sum of an
infinite set of weighted simplest hysteresis operators µ(α, β)γαβ[u(·)] shown in Fig. 2.5 which
indicates that the input u(t) is applied to a system of parallel connected weighted two-position
relays, and then their individual outputs integrated over all the values ofα and β [7].
u(t) f (t) P γα1β1[u(·)](t)µ(α1, β1) γα2β1[u(·)](t)µ(α2, β1) γα2β2[u(·)](t)µ(α1, β1) γαnβn[u(·)](t)µ(αn, βn)
Figure 2.5: Block diagram interpretation of the Preisach model.
2.2.1
Hysteresis Operator
The hysteresis operatorγαβ is also called the relay operator where the valuesα and β
corre-spond to “up” and “down” switching values of input, respectively [7, 8]. In this thesis, we assume
thatα ≥ β. The operator is characterized by its switching values α and β, and illustrated by a
rect-angular loop on the input-output phase diagram, depicted in Fig. 2.6(a). As shown in the figure,
remains to be 1 until the input is below the valueβ. To demonstrate further the operator function,
Fig. 2.6(b) shows the operator output when the input is increasing and Fig. 2.6(c) illustrates the operator output when the input is decreasing.
t 1 0 t 1 0 1 0 (a) (b) (c) α α α β β β γαβ[u(·)](t) γαβ[u(·)](t) γαβ[u(·)](t) u(t) u(t) u(t) u(t) u(t)
Figure 2.6: The hysteresis operator: (a) the phase diagram; (b) the operator output when the input
u(t) is increasing; (c) the operator output when the input u(t) is decreasing.
2.3
Preisach Plane
The relay operator is the core of the Preisach model; its function can be clarified via the
geometric interpretation, called the Preisach plane. The triangle region P shown in Fig. 2.7 is
referred to as the Preisach plane, where the region is defined as P = {(α, β) | 0 ≤ β ≤ α, 0 ≤
α ≤ umax, 0 ≤ β ≤ umax}, with umax denoting the maximum value of inputu(t). Because the
operator output is zero for(α, β) outside the region P, we can assume that the weighting function
µ(α, β) also equals zero outside the region P. The relay operator output for a given input can be
described as to paint the Preisach plane into two subregionsP+andP−where
P+(t) = {(α, β) ∈
P | γαβ[u(·)](t) is + 1} and P−(t) = {(α, β) ∈ P | γαβ[u(·)](t) is 0}; hence P+(t) ∪ P−(t) = P
at all times. For a given input u(t), there exits a one-to-one mapping between operators γαβ[u(·)]
and each point(α, β) of the Preisach plane. The point (α, β) corresponding to the relay operator
output, with switching valuesα and β, is at state +1 if it lies in P+and at state 0 if it is inP−. This
geometrical interpretation of the relay operator enables a realization of the Preisach model using the first-order reversal curves (FORC), as discussed later in section 2.4.
(α , ) α α β β β γαβ[u(·)](t) u(t) umax P
Figure 2.7: The Preisach plane.
2.3.1
Geometric Interpretation of the Relay Operator
Assume first the system is initially at rest; that is, the system input and output are both zero and
the values of allγαβ[u(·)] are also zero. As the input increases from zero to a value u1, γαβ[u(·)]
will be switched to the +1 state if the switching value α are less than the current input value u1.
In this case, the Preisach plane is divided into two regions, shown in Fig. 2.8(a). When the input
changes direction and decreases from u1 to a valueu2, then the operator outputγαβ[u(·)] will be
switched to the 0 state for the value β greater than u2. This direction change will result in the
change in the Preisach plane, as depicted in Fig. 2.8(b). Therefore, as the input decreases, the
regionP−widens or equivalently the region
P+narrows. ( , ) ( , ) ( , ) ( , ) (a) (b) α α β β α1 α1 α1 α2 α2 α2 β1 β1 β1 β1 β2 β2 γαβ[u(·)](t) γαβ[u(·)](t) γαβ[u(·)](t) γαβ[u(·)](t) u(t) u(t) u(t) u(t) u1 u1 u1 u2 u2 u2 P− P− P+ P+
Figure 2.8: The interpretation of relay operator via Preisach plane: (a) when the input increases
From the division of the Preisach plane, it can record ascending and descending changes of the
inputu(t) geometrically. By generalizing the previous analysis, the region P is divided into P+
andP−by the interface link; the interface link moves up asu(t) increases, and moves from right
to left asu(t) decreases as illustrated in Fig. 2.9. Using the results that the relay operator output is
0 inP−and 1 inP+and the disjoint of the two subregions, the Preisach model output (2.1) can be
written as: f (t) = Z Z α≥β µ(α, β)γαβ[u(·)](t)dαdβ = Z Z P µ(α, β)γαβ[u(·)](t)dαdβ = Z Z P+ µ(α, β)γαβ[u(·)](t)dαdβ + Z Z P− µ(α, β)γαβ[u(·)](t)dαdβ = Z Z P+ µ(α, β)dαdβ (2.2)
The integration region is constrained to be inP+.
( , ) ( , ) ( , ) ( , ) α β m0 m0 M1 M1 M1 M2 M2 M2 m1 m1 m1 m2 m2 u(t) t P− P+
Figure 2.9: A changing input and its memory formation mechanism of the Preisach plane.
2.3.2
Wiping-Out Property
Next we discuss the wiping-out property; it clarifies that to determine the current Preisach model output doesn’t have to memorize all the past input extremum but only parts of them [2,
7]. This property is illustrated in Fig. 2.10. Assume the input u(t) have reached to one local
maximum value M1 and then u(t) reaches to another local maximum value M2 later. If M2 is
smaller thanM1, the Preisach plane as shown in Fig. 2.10(a), theM1should be record to determine
the interface link of the Preisach plant as stated in section (2.3.1). If M2 is equal to or greater
thanM1, however, the previous input maximumM1 in the Preisach plane would be wiped out as
illustrated in Fig. 2.10(b)(c).
The wiping-out property means that the input extremum which be wiped out are not needed to determine the interface link of the Preisach plane; what we need to memorize are called the
(a) (b) (c) α α α β β β M1 M1 M1 M1 M1 M1 M2 M2 M2 M2 u(t) u(t) u(t) t t t P− P− P− P+ P+ P+
Figure 2.10: Wiping out property: (a)M2 < M1 (b)M2 = M1 (c)M2 > M1.
alternating series of the input extremum, explained as below. Consider a particular input u(t) in
the time interval [t0, t′] shown in Fig. 2.11. Assume that u(t0) = m0 = 0. We use the notations
Tkfor the time the global maximumMkis reached andtkfor the time the global minimummk is
s(t0, u(t′)) = {M1, m1, . . . , Mk, mk, . . .} (k = 1, 2, . . .), where Mk = max [tk−1,t′] u(t) = u(Tk) (2.3) and mk = max [Tk,t′] u(t) = u(tk) (2.4) u(t) t t1 t2 tk T1 T2 Tk M1 M2 Mk m1 m2 mk t0
Figure 2.11: Alternating series of the input extremum.
2.4
First Order Reversal Curve
The geometric interpretation of the relay operator via the Preisach plane poses a practical
way to determine the weighting parameters µ(α, β) and to realize the Preisach model. The
ap-proach uses a set of first order reversal curves (FORC) which can be obtained experimentally. A FORC is composed of two branches, an ascending branch and a descending branch. When the
input increases from the initial 0 to some value M, the corresponding Preisach plane is shown
in Fig. 2.12(a), yielding an ascending branch on the phase diagram. One practical demonstration of a practical micro piezoelectric actuator for an increasing input is also shown in Fig. 2.13(a).
The input, then, decreases from M to some value m, the corresponding Preisach plane is shown
in Fig. 2.12(b), a descending branch is formed on the input-output phase diagram. One practical measurement is also shown in Fig. 2.13(b). The term “first-order” is used to emphasize the fact that each of these curves is formed after the first reversal of the input.
(a) (b) α α β β m M M P+ P+
Figure 2.12: The change on the Preisach plane to yield a FORC: (a) the input increases from 0 to
some valueM; (b) then the input decreases from M to some value m.
0 0.5 1 1.5 2 2.5 3 3.5 0 5 10 15 20 25 30 35 input voltage (V) displacement ( µ m) an ascending branch 0 0.5 1 1.5 2 2.5 3 3.5 0 5 10 15 20 25 30 35 input voltage (V) displacement ( µ m)
an ascending branchand a descending branch
an ascending branch a descending branch an ascending branch
(a) (b)
M m M
Figure 2.13: First order reversal on the phase diagram:(a) an ascending branch; (b) an descending branch.
In order to connect FORC and the algebraic calculation of model output, define a function
fαβ(M, m) is the output of the plant as the input increases from 0 to M and then decreases to m.
Hence, fαβ(M, m) denotes the ending value of a descending branch of FORC. This value can be
obtained using (2.2) and Fig. 2.12, yielding
fαβ(M, m) = Z Z P+ µ(α, β)dαdβ = Z M 0 Z α 0 µ(α, β)dβdα − Z M m Z α m µ(α, β)dβdα (2.5)
SubstitutingM for m in (2.5), we obtain
fαβ(M, M) = Z M 0 Z α 0 µ(α, β)dβdα − Z M M Z α M µ(α, β)dβdα = Z M 0 Z α 0 µ(α, β)dβdα (2.6)
The valuefαβ(M, M) represents the beginning value of a descending branch and equivalently the
ending value of a ascending branch of FORC. This result enables us to determine the weighting parameters experimentally.
2.4.1
Experimental Determination of Weighting Function
The weighting µ(α, β) function can be determined from the FORC which can be obtained
experimentally. That is, one can collect experimentally a sets offαβ(M, m) with various M and m
values by measuring the system output with an input increasing from 0 toM, and then decreasing
fromM to m. Since we can obtain the relation below by taking partial derivatives of (2.5),
∂2f αβ(M, m) ∂M∂m = − ∂ ∂M Z M m ∂ ∂m Z α m µ(α, β)dβ dα = ∂ ∂M Z M m µ(α, m)dα = µ(M, m) (2.7)
It is clear from (2.7) that if thefαβ(M, m) is known for all points (M, m), the weighting
func-tionµ(α, β) can be determined. Then, the Preisach model output can be evaluated using µ(α, β)
and the model equation. The approach using weighting parameter first to evaluate the model out-put is unappealing because of two difficulties. First, to evaluate the Preisach model outout-put using
µ(α, β) requires numerical double integration, which is a complex and time-consuming procedure.
Second, the determination ofµ(α, β) requires the differentiations of experimental data, which may
amplify errors in the experimental data. One practical and popular approach to evaluate model out-put using FORC is discussed in the following.
2.4.2
Numerical Realization of the Preisach Model
We have seen that the evaluation of the Preisach model output using µ(α, β) is complex and
may amplify errors. An alternative method, called numerical implementation of the Preisach
model, is discussed. This method evaluates the Preisach model using fα,β(M, m) directly. In
order to develop the numerical implementation of the Preisach model, the ascending input u(t)
is discussed first. The change of this input is from 0 to M1 and then tom1, The corresponding
change on Preisach plane is shown in Fig. 2.14(a)(b)(c), and the output of the Preisach model can be written as: f (t) = Z M1 0 Z α 0 µ(α, β)dβdα − Z M1 m1 Z α m1 µ(α, β)dβdα + Z u(t) m1 Z α m1 µ(α, β)dβdα (2.8)
Setm0 = 0, from (2.5), (2.6) and (2.8), we derive the following expression for f (t) in the case of
increasing input:
f (t) = [fαβ(M1, M1) − fαβ(M1, m0)] − [fαβ(M1, M1) − fαβ(M1, m1)]
+ [fαβ(u(t), u(t)) − fαβ(u(t), m1)]
In general, asu(t) is ascending, f (t) =
n
X
k=1
[fαβ(Mk, Mk) − fαβ(Mk, mk−1)] + [fαβ(u(t), u(t)) − fαβ(u(t), mn)] (2.10)
wherem0 = 0 and {M1, m1, . . . , Mn, mn} is the alternating series of the input.
Next consider that the inputu(t) is descending, the change of this input is from 0 to M1and then
tom1, and then toM2, the corresponding change on Preisach plane is shown in Fig. 2.14(a)(b)(c)(d).
Similarly, the output of the Preisach model can be written as:
f (t) = [fαβ(M1, m1) − fαβ(M1, m0)]
+ [fαβ(M2, M2) − fαβ(M2, m1)] − [fαβ(M2, M2) − fαβ(M2, u(t))]
= [fαβ(M1, m1) − fαβ(M1, m0)] + [fαβ(M2, u(t)) − fαβ(M2, m1)] (2.11)
In general, asu(t) is descending,
f (t) =
n
X
k=1
[fαβ(Mk, mk) − fαβ(Mk, mk−1)] + [fαβ(Mn+1, u(t)) − fαβ(Mn+1, mn)] (2.12)
wherem0 = 0 and {M1, m1, . . . , Mn, mn, Mn+1} is the alternating series of the input.
In another aspect, we can also utilize the alternating series to simplify (2.10) and (2.12).
As-sumeu(t) is ascending, s = {M1, m1, . . . , mn} and u(tn) = mn; (2.10) can equivalently written
as:
f (t) = f (tn) + [fαβ(u(t), u(t)) − fαβ(u(t), mn)] (2.13)
Assume u(t) is descending, s = {M1, m1, . . . , Mn−1, mn−1, Mn} and u(Tn+1) = Mn+1; (2.10)
can equivalently written as:
f (t) = f (Tn) + [fαβ(Mn, u(t)) − fαβ(Mn, Mn)] (2.14)
From the above discussion, the numerical implementation of the Preisach model can be
ob-tained by combining the Preiasch plane and fαβ(M, m), which can be used for evaluation of the
(a) (b) (c) (d) α α α α β β β β u(t) u(t) m1 m1 m1 M1 M1 M1 M1 M2
Figure 2.14: Corresponding change on Preisach plane of an input:(a) from 0 toM1; (b) fromM1
Chapter 3
Preisach Model Realization and Inverse
Preisach Model
This chapter investigates practical approaches to find fαβ(M, m), including the polynomial
approximation, our new approach. In Chapter 2, a numerical implementation of the Preisach
model has been developed. The two formulations (2.13) and (2.14) expressed with fα,β(M, m)
can be used for evaluating the Preisach model output. From the experimental FORC measurement,
we can collect the data offαβ(M, m) which can be used to predict the Preisach model output f (t).
In practice, however, it is impossible to obtain all values of fαβ(M, m) because the completed
information offαβ(M, m) requires infinite and uncountable experiments. One realizable method to
determine values offαβ(M, m) is to find finite data of fαβ(M, m) and store it as a table. When we
need a value offαβ(M, m), we can obtain its value by the interpolation using its neighbors stored
in the table. In this chapter, we introduce three methods to approximate fαβ(M, m). They are
the table method [10, 13], the polynomial approximation and adaptive polynomial approximation, respectively.
The remainder of this chapter introduces inverse Preisach model [9, 13]. The control theory is focus on determine an input feed to the system such that desired output has a good tracking with system output. One way to achieve this goal is using inverse system; that is, to generate a system has inverse relation between input and output corresponding with the original system. Realization the Preisach model using the polynomial approximation can describe the system with hysteresis via polynomial equations; it’s possible to find the inverse system with hysteresis by solving polynomial equations.
3.1
Table Method
Linear interpolation is a very intuitive way to estimate values between two known values [11].
As an example, if you want to estimate x shown in Fig. 3.1, a fundamental formula for finding
internal point of division is:
x = (1 − r)a + rb (3.1)
a x b
r 1 − r
Figure 3.1: Linear interpolation betweena and b.
If a point (M, m) belongs to a rectangular cell formed by (M1, m1), (M1, m2), (M2, m1),
(M2, m2), and we have fαβ(M1, m1) = fM1m1, fαβ(M1, m2) = fM1m2, fαβ(M2, m1) = fM1m1,
fαβ(M2, m2) = fM2m2 shown in Fig. 3.2. Now we want to approximatefαβ(M, m) = fM m[10].
Assume thatfαβ(M1, ·), fαβ(M1, ·), fαβ(M1, ·) are linear functions, changes from fM1m1 tofM2m1
and fromfM1m2 tofM2m2 are also linear. Using 3.1, we have:
fL= (M2− M)fM1m2 + (M − M1)fM2m2 M2− M1 , fR = (M2− M)fM1m1 + (M − M1)fM2m1 M2− M1 , fM m = (m − m 1)fL+ (m2− m)fR m2− m1 (3.2) fM1m2 fM1m1 fM2m2 fM2m1 fR fL fM m fM m M2− M M − M1 m2− m m − m1 fαβ(M1, ·) fαβ(M, ·) fαβ(M2, ·)
Figure 3.2: Approximatefαβ(M, m) where (M,m) belongs to a rectangular cell.
If a point(M, m) belongs to a triangular cell formed by (M1, m1), (M2, m1), (M2, m2), and
want to approximatefαβ(M, m) = fM m[10]. Under the same assumption and a little modification, we have: fR = fM2m2 ,fL= (M2− M)fM1m2 + (M − M1)fM2m2 M2− M1 , fM m = (m − m1 )fL+ (m2 − m)fR m2− m1 (3.3) fM1m1 fM2m2 fM2m1 fL fR fM m fM m M2− M M − M1 m2− m m − m1
Figure 3.3: Approximatefαβ(M, m) where (M,m) belongs to a triangular cell.
3.2
Polynomial Approximation
It is convenient to investigate the measured data using an analyzable function. This is also called curve fitting and the polynomial function is commonly used [11]. Our goal in this section is to fit the samples of a finite number of FORC by polynomial functions.
3.2.1
Polynomial Curve Fitting
Consider the polynomial function given by
F (x) = c0+ c1x + . . . cnxnorF (x) =
n
X
j=0
cjxj (3.4)
to fitk pairs of data (xi, yi), i = 1, . . . k. Here we define the difference between the data and the fit
function at each point isri = F (xi) − yi, whereriis called the residual for the data point. We can
choose one normP |ri|, infinite norm max
alli |ri| or two normP r
2
i as the criterion for minimization.
For statistical and computational reasons [12], two norm is mostly used and it is given by
ρ = 1 k k X i=1 ri2 (3.5)
and √ρ is called the root mean square (RMS) error. In other word, our object is to find the
coefficientscj of a polynomialF (x) such that it can fit the data yiin a least squares sense.
Denote r= r1 r2 · · · rk T , and r = r1 r2 .. . rk = Pn j=0cjx j 1 − y1 Pn j=0cjx j 2 − y2 .. . Pn j=0cjxjk− yk = 1 x1 x21 · · · xn1 1 x2 x22 · · · xn2 .. . 1 xk x2k · · · xnk c0 c1 c2 .. . cn − y1 y2 .. . yk = Ac − y (3.6) where A= 1 x1 x21 · · · xn1 1 x2 x22 · · · xn2 .. . 1 xk x2k · · · xnk , c= c0 c1 c2 .. . cn , y= y1 y2 .. . yk . and now ρ = 1 k k X i=1 r2 i = 1 kkrk 2 2 = 1 kkAc − yk 2 2 (3.7) To find c, differentiatekrk22 = cTATAc
− 2yTAc+ yTy with respect to c and set it to zero.
2cTATA− 2yTA= 0 ⇒ ATAc= ATy⇒ c = (ATA)−1ATy (3.8)
Here we suppose thatxh 6= xlforh 6= l. Because Aa = 0 is only valid under the condition a = 0
ask ≥ n, A is full rank implies that ATA is nonsingular.
3.2.2
Polynomial Approximation
For an particular M, fαβ(M, m) is a monotonic decreasing function of single variable m
shown in Fig. 3.4, so using quadratic polynomial functions to fit fαβ(M, m) is a suitable choice.
That is,
fαβ(M, m) = c0(M) + c1(M)m + c2(M)m2 (3.9)
where ci(M) is a notation to indicate the coefficient corresponding with M. However, the
coef-ficients of curve fitting method varies with M, so to know all of them is also impossible. Hence
Mk < M < Mk+1, ci(M) = Mk+1ci(Mk) + Mkci(Mk+1) Mk+1 − Mk , wherei = 0, 1, 2. (3.10) 0 0.5 1 1.5 2 2.5 3 3.5 4 0 5 10 15 20 25 30 35 input voltage m (V) displacement ( µ m) measured data fitting curve
Figure 3.4: Using a function to fitfαβ(M, m).
According to (2.13) and (2.14), the Preisach model output can be determined by a storage part
and a update part. Now we substitute for The update part by using (3.9). If the input u(t) is
ascending,
f (t) = f (tn)
| {z }
storage part
+ [fαβ(u(t), u(t)) − fαβ(u(t), mn)]
| {z }
update part
(3.11)
= f (tn) + [c2(u(t))u2(t) + c1(u(t))u(t) + c0(u(t))] − [c2(u(t))m2n+ c1(u(t))mn+ c0(u(t))]
= f (tn) + c2(u(t))(u2(t) − m2n) + c1(u(t))(u(t) − mn) (3.12)
Defineci(p∆) = csi[p], i = 0, 1, 2, for u(t) = p∆ + r, 0 ≤ r < ∆,
ci(u(t)) = (1 − r ∆)csi[p] + r ∆csi[p + 1] (3.13) then (3.12) becomes f (t) = f (tn) + θTx1 (3.14)
where θ = cs2[p] cs2[p + 1] cs1[p] cs1[p + 1] , andx1 = (u(t) − mn) (1 − r ∆)(u(t) + mn) r ∆(u(t) + mn) (1 − r ∆) r ∆ (3.15)
If the inputu(t) is descending,
f (t) = f (Tn) | {z } storage part + [fαβ(Mn, u(t)) − fαβ(Mn+, Mn)] | {z } update part (3.16) = f (Tn) + c2(Mn)(u2(t) − Mn2) + c1(Mn)(u(t) − Mn) (3.17) ForMn = p∆ + R, 0 ≤ R < ∆, ci(Mn) = (1 − R ∆)csi[p] + R ∆csi[p + 1] (3.18) then (3.17) becomes f (t) = f (Tn) + θTx2 (3.19) where θ = cs2[p] cs2[p + 1] cs1[p] cs1[p + 1] , andx2 = (u(t) − Mn) (1 −R ∆)(u(t) + Mn) r ∆(u(t) + Mn) (1 − R ∆) R ∆ (3.20)
For higher resolution, both the table method and the polynomial approximation needs more memory; but the memory growth of the polynomial approximation is less than the table method.
Consider a division of the input region into N parts, the number of data should be storage is N(N +1)2
for table method and2N for polynomial approximation.
3.3
Adaptive Polynomial Approximation
The polynomial approximation is easier with respect to modification than the the table method. In this section we develop an adaptive algorithm for identification of polynomial coefficients which can improve model accuracy [14]. Base on polynomial approximation method discussed in
pre-vious section, we obtain a set of coefficients c that can be use to describe fαβ(M, m), however,
these coefficients may be inaccurate. A system model with inaccurate parameters can be regulated by update its parameters via adaptive filters such as LMS or RLS, shown in Fig. 3.5, the adaptive model output is compared with the plant output, and the error is used to update the adaptive system parameters. So now we adopt adaptive filters to modify polynomial approximation Preisach model for accuracy improvement.
f Input data u Plant output fm Model output Error data e Plant Model Adaptive algorithm
Figure 3.5: Block diagram of adaptive system identification.
3.3.1
Adaptive Polynomial Approximation with LMS
Steepest descent algorithm is the most common and easiest one of adaptive algorithms [16].
Assume the adaptive model outputfm(θ, u) is a function of parameters θ and the input u, and the
plant outputf (u) is a function of the input u. Define the error as
e = f (u) − fm(θ, u) (3.21)
and a cost function
J = e2 = (f (u) − fm(θ, u))2 (3.22)
The steepest descent algorithm is:
θk+1 = θk− µ∇θJ = θk− 2eµ∇θ(f (u) − fm(θ, u)) ⇒ θk+1 = θk+ 2µe∇θfm(θ, u) (3.23)
for an FIR, fm(θ, u) = θTu, then ∇θfm(θ, u) = u, so the steepest descent algorithm can be
simplified as
θk+1 = θk+ 2µeu (3.24)
Using (3.14), (3.19), (3.24), we obtain following update formulations:
If the inputu(t) is ascending,
θk+1= θk+ 2µex1 = cs2,k[p] cs2,k[p + 1] cs1,k[p] cs1,k[p + 1] + 2µe(u(t) − mn) (1 − r ∆)(u(t) + mn) r ∆(u(t) + mn) (1 − r ∆) r ∆ (3.25)
If the inputu(t) is descending,
θk+1 = θk+ 2µex2 = cs2,k[p] cs2,k[p + 1] cs1,k[p] cs1,k[p + 1] + 2µe(u(t) − Mn) (1 − R ∆)(u(t) + Mn) R ∆(u(t) + Mn) (1 − R ∆) R ∆ (3.26)
3.4
Inverse Preisach Model
We have seen that how does the Preisach model character the system with hysteresis and introduced three methods to realize it. However, to establish the inverse Preisach model is also an important problem in control aspect. As the Preisach model is realized by the polynomial approximation, the system with hysteresis can be described by a series of polynomial functions. This section establishes the inverse Preisach model by solving polynomial functions.
3.4.1
Cubic Equation
Before combining the inverse Preisach model with solving polynomial functions, we prepare formula for cubic equation which is needed in solving polynomial functions later [15].
A cubic equation,u3+ pu2+ qu + r = 0 can be reduced to the simpler form,
x3+ ax + b = 0 (3.27)
by substituting foru the value, x − q/3. Where
a = 1 3(3q − p 2) , and b = 1 27(2p 3 − 9pq + 27r) (3.28)
For solution, let
A = 3 s −2b + r b2 4 + a3 27 , andB = 3 s −2b − r b2 4 + a3 27 (3.29)
then the values ofx will be given by
x = A + B , − A + B2 +A − B
2 √
−3 , −A + B2 − A − B2 √−3 (3.30)
However, only real roots are needed in our problem. For this reason, next we introduce trigono-metric solutions of the cubic equation. Consider the trigonotrigono-metric identity
4 cos3θ − 3 cos θ − cos 3θ ≡ 0 (3.31)
Letx = m cos θ, then
x3+ ax + b ≡ m3cos3θ + am cos θ + b ≡ 4 cos3θ − 3 cos θ − cos 3θ ≡ 0 (3.32)
Hence 4 m3 = − 3 am = − cos 3θ b ⇒ m = 2 r −a3 , cos 3θ = 3b am = D (3.33)
Ifa < 0 and |D| ≤ 1, then
cos 3θ = D ⇒ θ = 13arccos D , θ = 1
3arccos(D ±
2π
3 ) (3.34)
and the cubic equation has three real roots:
x = m cos θ = m cos(1 3arccos D) , m cos( 1 3arccos(D ± 2π 3 )) (3.35) Ifa < 0 and D > 1,
cos 3θ = D = cosh(j(3θ + 2kπ)) ⇒ θ = −13j cosh−1D , − 1
3j cosh
−1
(D ±2π3 ) (3.36)
and the cubic equation has one real root:
x = m cos θ = m cosh jθ = m cosh(1
3cosh −1D) (3.37) Ifa < 0 and D < −1, cos(3θ − π) = −D = cosh(j(3θ + (2k − 1)π)) ⇒ θ = −1 3j cosh −1 D , − 1 3j cosh −1 (D ±π 3) (3.38) and the cubic equation has one real root:
x = m cos θ = m cosh jθ = m cosh(1
3cosh −1 (−D)) (3.39) Ifa > 0, then m = 2 r −a3 = j2r a 3 = jm ′ ,D = 3b am = −j 3b am′ = −jD ′ (3.40) cos 3θ = D ⇒ j cos 3θ = D′ = j sin((2k + 1 2π) − 3θ) = sinh(j[(2k + 1 2π) − 3θ]) ⇒ θ = π6 + 1 3j sinh −1D′ , 5π 6 + 1 3j sinh −1D′ , 3π 2 + 1 3j sinh −1D′ (3.41) and the cubic equation has one real root:
x = m cos θ = m cos(3π 2 + 1 3j sinh −1D′ ) = m cos(3π 2 ) cos( 1 3j sinh −1D′ ) − m sin(3π 2 ) sin( 1 3j sinh −1D′ ) = m sin(1 3j sinh −1D′ ) = jm sinh(1 3sinh −1D′ ) = m′ sinh(1 3sinh −1D′ ) (3.42)
3.4.2
Solving Polynomial of the Preisach Model
Because of the existence of formula for solving polynomial function, the input can be obtain by solving (3.12) or (3.17) for the desired output.
As the desired outputfdis ascending, we have to solve (3.12). Forp∆ ≤ u < (p + 1)∆,
ci(u) = (u − p∆)c
si[p + 1] + ((p + 1)∆ − u)csi[p]
∆
= csi[p + 1] − csi[p]
∆ u + ((p + 1)csi[p] − pcsi[p + 1]) , i = 1, 2 (3.43)
Replaceci(u(t)) in the (3.12) with (3.43)
u3 + ((p + 1)cs2[p] − pcs2[p + 1])∆ + cs1[p + 1] − cs1[p] cs2[p + 1] − cs2[p] u2 + ((p + 1)cs1[p] − pcs1[p + 1])∆ + (cs2[p] − cs2[p + 1])m 2 n+ (cs1[p] − cs2[p + 1])mn cs2[p + 1] − cs2[p] u + (pcs2[p + 1] − (p + 1)cs2[p])m 2 n+ (pcs1[p + 1] − (p + 1)cs1[p])mn+ f (tn) − fd cs2[p + 1] − cs2[p] ∆ = u3+ pu2+ qu + r = 0 (3.44)
This is a cubic equation, and can be solved by formula introduced last section.
As the desired outputfdis descending, we have to solve (3.17). ForM = p∆+R, 0 ≤ R < ∆,
c2(Mn)u2+ c1(Mn)u + f (Tn) − fd− c2(Mn)Mn2− c1(Mn)Mn = au2+ bu + c = 0 (3.45)
whereci(Mn) can be computed using (3.18). (3.45) is a quadratic equation, and
u = −b ±
√
b2− 4ac
Chapter 4
Applications for Piezoelectric Actuators
Piezoelectric actuators have advantages of high stiffness, fast frequency response and high precision [16] have been widely used in precise industrial applications [17, 18]. The displacement of a micro piezoelectric actuator varies with actuator voltage and the relation between actuator voltage and displacement is a kind of hysteresis phenomenon. In this chapter, we use different realizations of Preisach model to predict the displacement of a micro piezoelectric actuator and the inverse Preisach model introduced in this thesis to control it. In these applications, the input of the Preisach model is actuator voltage and the output of the Preisach model is displacement.
4.1
Experimental Setup
The micro piezoelectric actuator studied in these experiments is model PMT 150/40
Transla-tion Stage PMT with Central Piezoelement, which has a maximum applied voltage of 150V and
a maximum displacement of 45µm. A dSPACE system (Model no. DS1104 PPC) is used for
acquisition digital data and real-time control. A power amplifier (Model no. PosiCon 150-3) with a gain of 30 is used to drive the micro piezoelectric actuator, so the actuator voltage given from
PC-Based system is0 ∼ 5V . The block diagram of the entire experimental system set-up is shown
in Fig. 4.2 and the real system is shown in Fig. 4.1.
To implement the Preisach model, several first order reversal curves measured for micro
piezo-electric actuator first. Divide 0 ∼ 5V into 10 segments, the input voltage of micro piezoelectric
actuator is increased from 0V to some divided point and then decreased to 0V . After the input
begins decreasing from some divided point to 0, we sample the displacements of micro
piezoelec-tric actuator every 0.45V change of input voltage and store these sampling data in Table 4.1. The
first column denotes the value ofM, and the first row denotes the value of m. In the predication
experiment we use a sinusoidal and input voltage that with a frequency of 1 rad/sec, an amplitude
Figure 4.1: Real system for experiment.
the predication experiment shown in Fig. 4.3, and the phase transition shown in Fig. 4.4. In the tracking control experiment we use a sinusoidal desired output that with a frequency of 1 rad/sec,
an amplitude of15µm and a bias of 16µm. The sampling time of the experiment is 0.01 sec., the
total time is 50 sec..
u piezoelectric f
actuator sensor A/D
D/A
PC with dSPACE Plug−in Board
power amplifier
Figure 4.2: The block diagram of the entire experimental system.
4.2
Table Method
The displacement of micro piezoelectric actuator and the output of Preisach model realized via table method are plotted together in Fig. 4.5(a). The difference between these two outputs, defined as the model error, are shown in Fig. 4.5(b). The phase transition of model and of plant
are plotted together in Fig. 4.6. The RMS of model error is1.0098µm, and the required number of
M/m 0 0.45 0.90 1.35 1.80 2.25 2.70 3.15 3.60 4.05 4.50 0 0 0.45 0.02 3.43 0.90 0.09 4.07 7.73 1.35 0.11 5.28 9.38 12.76 1.80 0.14 5.36 10.13 14.11 18.15 2.25 0.23 5.68 10.79 15.32 19.34 23.32 2.70 0.28 5.94 11.31 16.22 20.58 24.48 28.03 3.15 0.36 6.17 11.72 16.68 21.27 25.41 29.01 32.65 3.60 0.45 6.35 12.09 17.29 21.96 26.27 30.13 33.51 36.49 4.05 0.48 6.61 12.35 17.58 22.40 26.81 30.74 34.35 37.47 40.15 4.50 0.67 6.87 12.87 18.07 22.98 27.33 31.38 34.93 38.13 40.96 43.18
Table 4.1: Sampling data of FORC.
0 5 10 15 20 25 30 35 40 45 50 0 1 2 3 4 time (sec) input voltage (V) 0 5 10 15 20 25 30 35 40 45 50 0 10 20 30 40 time (sec) output displacement ( µ m)
Figure 4.3: The input and output of the micro piezoelectric actuator.
4.3
Polynomial Approximation
We apply the least-square approximation to sampling numbers of Table 4.1 to identify the polynomial coefficients directly in this experiment. The polynomial coefficients are stored in
Ta-ble 4.2. Base on (3.12) and (3.17), the polynomial coefficients c0 is not needed for the Preisach
0 0.5 1 1.5 2 2.5 3 3.5 4 4.5 0 5 10 15 20 25 30 35 40 45 input voltage (V) output displacement ( µ m) phase transition
Figure 4.4: The phase transition of micro piezoelectric actuator.
0 5 10 15 20 25 30 35 40 45 50 0 10 20 30 40 time (sec) displacement ( µ m) 0 5 10 15 20 25 30 35 40 45 50 −2 −1 0 1 2 3 time (sec) error ( µ m) model plant (a) (b)
Figure 4.5: Prediction of micro piezoelectric actuator using Preisach model realization via table method.
0 0.5 1 1.5 2 2.5 3 3.5 4 4.5 −5 0 5 10 15 20 25 30 35 40 45 input voltage (V) output displacement ( µ m) model plant
Figure 4.6: Phase transition of micro piezoelectric actuator using Preisach model realization via table method.
Preisach model which is realized via table method shown in the Fig. 4.7(a), and the model error shown in the Fig. 4.7(b). The phase transition of model and of plant are plotted together in Fig. 4.8.
The RMS of model error is1.4817µm, and the number of memory for polynomial coefficients is
22. M c2 c1 c0 0 0 0 0 0.45 0 7.5778 0.02 0.90 -0.79012 9.2 0.09 1.35 -2.2099 12.328 0.1275 1.80 -1.1111 11.949 0.174 2.25 -1.0326 12.543 0.2525 2.70 -1.117 13.302 0.24619 3.15 -1.0256 13.439 0.37625 3.60 -1.091 13.954 0.39576 4.05 -1.058 14.075 0.49745 4.50 -1.0752 14.269 0.73091
0 5 10 15 20 25 30 35 40 45 50 0 10 20 30 40 time (sec) displacement ( µ m) 0 5 10 15 20 25 30 35 40 45 50 −4 −3 −2 −1 0 1 time (sec) error ( µ m) model plant (a) (b)
Figure 4.7: Prediction of micro piezoelectric actuator using Preisach model realization via poly-nomial approximation. 0 0.5 1 1.5 2 2.5 3 3.5 4 4.5 0 5 10 15 20 25 30 35 40 45 input voltage (V) output displacement ( µ m) model plant
Figure 4.8: Phase transition of micro piezoelectric actuator using Preisach model realization via polynomial approximation.
4.4
Polynomial Approximation with Adaptive Identification
In this experiment we develop the least mean square (LMS) adaptive algorithm to obtain ac-curacy polynomial coefficients and the data in Table 4.2 are used as the initial parameters. We test
the step sizeµ = 10−4,µ = 10−3andµ = 0.0005
1+u(t) and the process of each is shown in Fig. 4.9. The
polynomial coefficients obtained by LMS adaptive algorithm with step size µ = 1+u(t)0.0005 are stored
in Table 4.3. 0 50 100 150 200 250 300 350 400 450 500 0 0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8 time (sec)
mean square error
µ = 0.0005/(1+u(t))
µ = 10−3
µ = 10−4
Figure 4.9: Mean square error of the prediction of micro piezoelectric actuator for varying step-size
µ using adaptive polynomial approximation.
M c2 c1 0 0.5718 3.8438 0.45 -2.7323 7.3383 0.90 -2.3325 9.2038 1.35 -2.1626 10.8419 1.80 -1.7700 11.8186 2.25 -1.4588 12.3549 2.70 -1.1540 12.4685 3.15 -1.2125 13.2395 3.60 -1.2565 13.9011 4.05 -1.1797 14.0324 4.50 -1.0938 13.9965
Now we use new polynomial coefficients stored in Table 4.2 to predict the displacement of micro piezoelectric actuator. The displacement of micro piezoelectric actuator is compared with output of Preisach model which is realized shown in the Fig. 4.10(a), and the model error shown in the Fig. 4.10(b). The phase transition of model and of plant are plotted together in Fig. 4.11.
The RMS of model error is1.5898µm.
0 5 10 15 20 25 30 35 40 45 50 0 10 20 30 40 time (sec) displacement ( µ m) 0 5 10 15 20 25 30 35 40 45 50 −1 −0.5 0 0.5 1 1.5 2 2.5 time (sec) error ( µ m) model plant (a) (b)
Figure 4.10: Prediction of micro piezoelectric actuator using Preisach model realization via adap-tive polynomial approximation.
Obviously a bias is exists in the Fig. 4.10(b). If we survey Fig. 4.4 again, then we can observe that there is another branch in the Fig. 4.4. This fact means the Preisach model should be modified in the beginning of the micro piezoelectric actuator work. The Fig. 4.12 shows the phase transition of the model error and the input in the beginning ascending branch. The relation between input
u(t) and error e(t) is roughly like a linear function u(t) = 0.83e(t) for 0 < u(t) < 2.7; for u(t) > 2.7, e(t) is roughly like a constant 2.2414. Then (3.12) in the beginning ascending branch
is modified :
f (t) =
f (tn) + c2(u(t))(u2(t) − m2n) + c1(u(t))(u(t) − mn) + 0.83u(t) if 0 < u(t) < 2.7,
f (tn) + c2(u(t))(u2(t) − m2n) + c1(u(t))(u(t) − mn) + 2.2414 ifu(t) > 2.7
(4.1) With this modification, we predict the displacement of micro piezoelectric actuator again. The displacement of micro piezoelectric actuator is compared with output of Preisach model which is
0 0.5 1 1.5 2 2.5 3 3.5 4 4.5 −5 0 5 10 15 20 25 30 35 40 45 input voltage (V) output displacement ( µ m) model plant
Figure 4.11: Phase transition of micro piezoelectric actuator using Preisach model realization via adaptive polynomial approximation.
0 0.5 1 1.5 2 2.5 3 3.5 4 4.5 −0.5 0 0.5 1 1.5 2 2.5 3 input voltage (V) error ( µ m)
Figure 4.12: Phase transition of the model error and the input.
transi-tion of model and of plant are plotted together in Fig. 4.14. The RMS of model error is0.4539µm.
The model error of the Preisach model realization via the table method and the modified polyno-mial approximation are plotted together in Fig. 4.15 for comparing.
0 5 10 15 20 25 30 35 40 45 50 0 10 20 30 40 time (sec) displacement ( µ m) 0 5 10 15 20 25 30 35 40 45 50 −1.5 −1 −0.5 0 0.5 time (sec) error ( µ m) model plant (a) (b)
Figure 4.13: Prediction of micro piezoelectric actuator using Preisach model realization via modi-fied polynomial approximation.
4.5
Tracking Control
This section presents experiments of tracking control of the micro piezoelectric actuator using the inverse Preisach model introduced in chapter 3. The block diagram of the tracking control system is shown in Fig. 4.16. The desired output compared with the displacement of the micro piezoelectric of this experiment are plotted together in Fig. 4.17(a). The deference between these two outputs, defined as the tracking error, are shown in Fig. 4.17(a). The RMS of tracking error
is 0.7079µm. To reduce tracking error, we use combine the inverse Preisach model with a PID
controller shown in Fig. 4.18 [9, 13]. The parameters of the PID controller are Kp = 0.0452,
Ki = 2.9828 and Kd = 0.0024. The desired output compared with the displacement of the micro
piezoelectric of this experiment are plotted together in Fig. 4.19(a). The tracking error are shown
0 0.5 1 1.5 2 2.5 3 3.5 4 4.5 0 5 10 15 20 25 30 35 40 45 input voltage (V) output displacement( µ m) model plant
Figure 4.14: Phase transition of micro piezoelectric actuator using Preisach model realization via modified polynomial approximation.
0 5 10 15 20 25 30 35 40 45 50 −2 −1.5 −1 −0.5 0 0.5 1 1.5 2 2.5 3 time (sec) error ( µ m) polynomial approximation table method
Figure 4.15: Model error of the Preisach model realization via the table method and the modified polynomial approximation.
inverse
Preisach model
piezoelectric
actuator
desired output fd real output f control signal uFigure 4.16: Block diagram of the tracking control system using inverse Preisach model.
0 5 10 15 20 25 30 35 40 45 50 0 5 10 15 20 25 30 time (sec) displacement ( µ m) 0 5 10 15 20 25 30 35 40 45 50 −2 −1 0 1 2 time (sec) displacement ( µ m) real output desired output (a) (b)
Figure 4.17: Tracking of micro piezoelectric actuator using inverse Preisach model.
inverse Preisach model piezoelectric actuator desired output fd real output f control signal u PID controller
Figure 4.18: Block diagram of the tracking control system using inverse Preisach model with PID controller.
0 5 10 15 20 25 30 35 40 45 50 0 5 10 15 20 25 30 time (sec) displacement ( µ m) 0 5 10 15 20 25 30 35 40 45 50 −1 −0.5 0 0.5 1 time (sec) error ( µ m) real output desired output (a) (b)
Figure 4.19: Tracking of micro piezoelectric actuator using inverse Preisach model with PID con-troller.
Chapter 5
Conclusions
A system with hysteresis nonlinearities is often difficult to describe accurately and may result in unstable behaviors if not controlled appropriately. In this thesis, we adopt the Preisach model to characterize hysteresis behaviors and propose to use polynomial approximation to realize the Preisach model; this approach is then used to model the displacement of a real micro piezoelectric actuator. The proposed approach can overcome the drawbacks of conventional table method; it requires less memory space and enables the parameter tracking of hysteresis elements. We have successfully obtained the polynomial coefficients to model the displacement of a micro piezoelec-tric actuator; the obtained model compared with that via the table method not only requires less
memory size but also yields a smaller modeling RMS error from1.0098µm via the table method to
0.4539µm. We also establish the inverse Preisach model base on the polynomial approximation;
this model is combined with the PID controller used for tracking control and yields small tracking
References
[1] D. Hughes and J. T.Wen, “Preisach modeling of piezoceramic and shape memory alloy
hys-teresis, ”Smart Mater. Struct., vol. 6, pp. 287-300, 1997.
[2] X. Zhou, J. Zhao, G. Song, and J. De Abreu-Garcia, “Preisach modeling of hysteresis and
tracking control of a Thunder actuator system,”Proc. SPIE, vol. 5049, pp. 112-125, 2003.
[3] B. N. Agrawal, M. A. Elshafei, and G. Song, “Adaptive antenna shape control using
piezo-electric actuators,”Act. Astronaut., vol. 40, no. 11, pp. 821-826, 1997.
[4] D. Croft, G. Shed, and S. Devasia, “Creep, hysteresis, and vibration compensation for
piezoactuators: Atomic force microscopy application,”ASME Journal of Dynamic Systems,
Measurement, and Control, vol. 123, no. 1, pp. 35-43, 2001.
[5] P. Mayhan, K. Srinivasan, S.Watechagit, and G.Washington, “Dynamic modeling and
con-troller design for a piezoelectric actuation system used for machine tool control,” J. Intell.
Mater. Syst. Struct., vol. 11, pp. 771-780, 2000.
[6] Goldfarb, M. and Celanovic, N., “A lumped parameter electromechanical model for
describ-ing the nonlinear behavior of piezoelectric actuators,”ASME Journal of Dynamic Systems,
Measurement, and Control, vol. 119, pp. 478-485, 1997.
[7] I. D. Mayergoyz,Mathematical Models of Hysteresis, New York: Springer-Verlag, 1991.
[8] L. Dupre, R. van Keer, and J. A. A. Melkebeek, “Identification of the relation between the
material parameters in the Preisach model and in the Jiles-Atherton hysteresis model,” J.
Appl. Phys., vol. 85, pp. 4376-4378, 1999.
[9] G. Song, Jinqiang Zhao, Xiaoqin Zhou, and J.A. De Abreu-Garcia, “Tracking control of a piezoceramic actuator with hysteresis compensation using inverse Preisach model,”
[10] P. Ge and M. Jouaneh, “Modeling hysteresis in piezoceramic actuators,”Prec. Eng., vol. 17, pp. 211-221, 2005.
[11] Gerald W. Recktenwald, Numerical Methods with Matlab: Implementations and
Applica-tions, Prentice-Hall, Upper Saddle River, NJ., 2001.
[12] J. L. Melsa and D. L. Cohn,Decision and Estimation Theory, McGraw Hill, 1978.
[13] P. Ge and M. Jouaneh, “Tracking control of a piezoceramic actuator,”IEEE Trans. Control
Syst. Technol., vol. 4, no. 3, pp. 209-216, May 1996.
[14] S. S. Haykin,Adaptive Filter Theory, Englewood Cliffs, NJ: 4th Ed., 2002
[15] Selby, Samuel M.,CRC Standard Mathematical Tables, CRC Press., 18th Ed., 1975.
[16] Slocum and Alexander H., Precision machine design, Englewood Cliffs, NJ: Prentice-Hall,
1992.
[17] Okazaki Y., “A Micro-positioning tool post using a piezo electric actuator for diamond
turn-ing machines,”Prec. Eng., vol. 12, pp. 151-156, 1990.
[18] G. S. Choi, H. S. Kim, and G. H. Choi, “A study on position control of piezoelectric
