Next Article in Journal
A Structural Characterisation of the Mitogen-Activated Protein Kinase Network in Cancer
Next Article in Special Issue
Inversed Model-Based Disturbance Observer Base on Adaptive Fast Convergent Sliding Mode Control and Fixed-Time State Observer for Slotless Self-Bearing Motor
Previous Article in Journal
Subordination Involving Regular Coulomb Wave Functions
Previous Article in Special Issue
Adaptive Sliding Mode Control Anticipating Proportional Degradation of Actuator Torque in Uncertain Serial Industrial Robots
 
 
Font Type:
Arial Georgia Verdana
Font Size:
Aa Aa Aa
Line Spacing:
Column Width:
Background:
Article

Classical Preisach Model Based on Polynomial Approximation and Applied to Micro-Piezoelectric Actuators

1
Department of Electrical Engineering, National Formosa University, Huwei 632, Taiwan
2
Department of Electrical and Control Engineering, National Chiao Tung University, Hsinchu 300, Taiwan
*
Author to whom correspondence should be addressed.
Symmetry 2022, 14(5), 1008; https://doi.org/10.3390/sym14051008
Submission received: 17 April 2022 / Revised: 9 May 2022 / Accepted: 12 May 2022 / Published: 16 May 2022

Abstract

:
In engineering applications, where we demand more and more precision, the modeling of systems with hysteretic nonlinearity has received considerable attention. The classical Preisach model (CPM) is currently the most popular for characterizing systems with hysteresis, and this model can represent the hysteresis with an infinite but countable first-order inversion curve (FORC). The table method is a method used to realize CPM in practice. The data in the table corresponds to a limited number of FORC samples. There are two problems with this approach: First, in order to reflect the timing effects of elements with hysteresis, it needs to consume a lot of memory space to obtain accurate data table. Second, it is difficult to come up with an efficient way to modify the data table to reflect the timing effects of elements with hysteresis. To overcome these shortcomings, this paper proposes to use a set of polynomials instead of the table method to implement the CPM. The proposed method only needs to store a small number of polynomial coefficients, and thus it reduces the required memory usage. In addition, to obtain polynomial coefficients, we can use least squares approximation or adaptive identification algorithms, which can track hysteresis model parameters. We developed an adaptive algorithm for the identification of polynomial coefficients of micro-piezoelectric actuators by applying the least mean method, which not only reduces the required memory size compared to the table method implementation, but also achieves a significantly improved model accuracy, and the proposed method was successfully verified for displacement prediction and tracking control of micro-piezoelectric actuators.

1. Introduction

With the development of computer and numerical techniques, the use of hysteresis models has become more and more popular in finite element (FE) method (FEM)-based computer-aided design (CAD) simulations for calculating magnetization properties and iron losses [1]. Hysteresis appears in many engineering devices or smart materials, such as micro-piezoelectric actuators (PZA), piezoelectric ceramics, shape memory alloys (SMA), and ferromagnetic components [2,3,4,5,6], which is a nonlinear phenomenon. Systems with hysteresis are often difficult to describe accurately, and without feedback control, may lead to unstable behavior of the system. When a device exhibits nonlinearity, the system often exhibits inaccuracies or oscillations, even instability caused by the non-differentiable hysteresis and non-memoryless properties of magnetic hysteresis [7,8,9]. Therefore, accurate hysteresis models are crucial for developing suitable control algorithms or applications using these systems [10]. Several models have been developed to characterize systems with hysteresis; commonly used models include the Bouc–Wen, Preisach-like Krasnosel’skii–Pokrovskii (KP), Prandtl–Ishlinskii model (PIM), and the classical Preisach model (CPM) [11]. CPM has long been known as the most familiar model for characterizing hysteretic behavior. It uses double integrals with relay operators and weighting parameters to describe the system input/output relation, which can be further represented by an infinite but countable first-order inversion curve (FORC), which can be experimentally measured in many applications [12]. The traditional method to implement CPM for systems with hysteresis approximate by sampling a finite number of FORCs and storing the sample data in a table [13]. The model output is then evaluated from the table data and the given input by simple linear interpolation. However, this approach has two drawbacks: One is that when the number and sampling resolution of FORC should be increased, it requires a large amount of memory space to accurately predict the hysteresis behavior. The other is that it is difficult to find an efficient way to adjust the data table to reflect changes in parameters with aging or time. As a result, execution accuracy may degrade over time and even require rebuilding of the table.
To overcome these drawbacks, we recommend implementing CPM using a set of polynomials rather than data tables. It only needs to store a small number of polynomial coefficients to reduce the required storage space. Specifically, polynomial coefficients can be obtained using a least squares approximation or an adaptive identification algorithm [14], allowing tracking of hysteresis model parameters. The proposed method has been applied to computer simulations of CPM output predictions. We apply the least mean squares (LMS) adaptation algorithm to develop an adaptive algorithm for identifying the polynomial coefficients of FORC. Compared to the implementation of the table method, the results show that it has higher model accuracy than the table method and requires significantly less memory. Furthermore, the proposed method combined with an adaptive algorithm can achieve dynamic modification of CPM parameters instead of rebuilding the table. According to [15], the proposed method has been applied to simulate the displacement of micro-piezo actuators. The modeled RMS error changes from 1.0098 µm for the table method to 0.4539 µm.
The remainder of this paper is divided to five sections including conclusions. Section 2 reviews the hysteresis and the Preisach model. Section 3 introduces the basic definition and the numerical implementation on the CPM. Section 4 present the computer simulations to verify capabilities of different realizations of the CPM. Section 5 demonstrates the applications of the micro-piezoelectric actuators. The final section presents the conclusions.

2. Realization Numerical of CPM

The classical Preisach model (CPM) [16,17] represents the input and output relations of a system with hysteresis by a double integral formula as below:
f ( t ) = a β μ ( α , β ) γ α β [ u ( · ) ] ( t ) d α d β
where u ( t ) denotes the CPM input, f ( t ) is the CPM output, μ ( α , β ) is the Preisach density function, and γ α β [ u ( ) ] ( t ) is an operator with two distinct outputs. The idea of the CPM is to regard a hysteresis transducer as the superposition of infinitely weighted hysteresis operators as shown in Figure 1.
The hysteresis operator is the most important core in the Preisach model. The values α and β of the hysteresis operator correspond to “up” and “down” switching values of input, respectively. In this paper, we assume that α β . The operator is characterized by its switching values α and β , which are represented by a rectangular loop on its input and output phase diagrams, as shown in Figure 2. As shown in the figure, the operator output of γ α β [ u ( ) ] ( t ) remains −1 until the input is higher than the value of α , and the output remains + 1 until the input is lower than the value of β . A practical expectation for practical applications is the limit of the integration, which requires us to set the saturation value of the input, that is, u s e t α β u s e t .
A progressive evaluation of the CPM is demonstrated by the following example, depicted in Figure 3. In the example, we calculate a first order reversal curve (FORC), a curve is formed after the first reversal of the input, and a second order reversal curve. Any higher order reversal curve can be constructed by using the analogical method.
Assume the system is first at negative saturation; that is, the CPM input is u s a t and all hysteresis operator outputs are 1 . Then, the initial CPM output is
f ( t ) = a β μ ( α , β ) γ α β [ u ( · ) ] ( t ) d α d β
As the CPM input increases from u max to M 1 , the hysteresis operator output will be switched to the 1 if the switching value α is less than the current input value M 1 . Then, the CPM output is
f ( T 1 ) = f s a t + 2 M 1 a β u s a t μ ( α , β ) d α d β
When the CPM input changes direction and decreases from M 1 to a value m 1 , the hysteresis operator output will be switched to the 1 state for the switching value β greater than m 1 . Then, the FORC is
f ( t 1 ) = f ( T 1 ) 2 M 1 a β m 1 μ ( α , β ) d α d β
Evaluation of f ( T 2 ) is the same with f ( T 1 ) , that is,
f ( T 2 ) = f ( t 1 ) 2 M 2 a β m 1 μ ( α , β ) d α d β
The wiping-out property appears as the input u ( t ) reaching M 1 and then reaches M 2 * . Hysteresis operators with a switching value α less than M 2 * but greater than M 1 stay at the + 1 state after T 2 , so they have no influence on the current output. Additionally, hysteresis operators with α greater than M 1 but less than M 2 * are originally at the 1 state and would be transitioned to + 1 . Summarizing all of the above,
f ( T 2 * ) = f ( T 1 ) + 2 M 2 * a β u s a t μ ( α , β ) d α d β = f s a t + 2 M 2 * a β u s a t μ ( α , β ) d α d β
It is clear that ( M 1 , F 1 ) and ( m 1 , f 1 ) , which were wiped out, are not necessary to determine f ( T 2 * ) in (5).
In the CPM, hysteresis operator is a concern associated with the impact of the past inputs for such scheme. The wiping-out property, however, clarifies that to determine the current CPM the output does not have to memorize all of the history of CPM but only parts of the extremum of the past inputs and outputs [18]. The CPM can be rewritten in such a way to dichotomize a CPM output into two terms: the shifted term, depending on the past outputs, and the integrated term, depending on the past inputs and the current input. We can define alternating series, which are kept to determine the limits of integration and the storage part of a CPM output, by generalizing the above-mentioned example as below. Consider a particular input u ( t ) for a CPM in the time interval [ t 0 , t ] and u ( t 0 ) = u s a t . We use the notations T k for the time the global maximum M k of u ( t ) is reached and t k for the time the global minimum m k of u ( t ) is reached, where M k and m k satisfy that
M k = max [ t k 1 , t ] u ( t ) = u ( T k ) , m k = max [ T k , t ] u ( t k )
The alternating series of the input in [ t 0 , t ] is defined as a tuple s i = { u s a t , M 1 , m 1 , , u ( t ) } , and the alternating series of the output in [ t 0 , t ] is s o = { f s a t , F 1 , f 1 , , f ( t ) } , where F k = f ( T k ) and f k = f ( t k ) . Furthermore, for a given alternating series s i and s o , except f ( t ) at some point, the CPM output is obtained by
f ( t ) = f n + 2 u ( t ) a β m n μ ( α , β ) d α d β
If u ( t ) is increasing, and
f ( t ) = F n 2 M n a β u ( t ) μ ( α , β ) d α d β
If u ( t ) is decreasing. In the ensuing discussion, the alternating series except f ( t ) are assumed to be recorded.
So far, we developed forthright mathematical Formulas (8) and (9) to simplify the CPM with the known Preisach density function and given alternating series. However, the identification of μ ( α , β ) requires the differentiation of experimental data, which may amplify errors in the experimental data. There are some algorithms [19,20] published to remove sensitivity of identification of μ ( α , β ) , but the approach using μ ( α , β ) to evaluate CPM outputs is still unappealing. Although we could find an accurate μ ( α , β ) of CPM, the double integration would be very time consuming and mathematically intractable. The congruency property of CPM supplies another way of representing a CPM through use of the FORC to replace integration with experimental data [16]. We introduce the property and how to use FORC to numerically implement a CPM by examining (8) and (9). It is expected that the integrated terms in (8) or (9) for a particular CPM input will be congruent with other CPM input if their limits of integration are the same. Ground on this fact, we can store integrated terms with variational limits of integration beforehand. The FORC measuring experiment is used for this notion. In order to connect FORC and the algebraic calculation of the CPM output, a function f a ( m , M ) is defined as the output of a CPM with s i = { u s a t , m , M } and a function f d ( M , m ) is the output of a CPM with s i = { u s a t , M , m } . Hence, f a ( m , M ) denotes an ascending FORC and f d ( M , m ) denotes a descending FORC. Thus, (8) is equivalent to
f ( t ) = f n + f a ( m n , u ( t ) ) f a ( m n , m n )
Additionally, (9) is equivalent to
f ( t ) = F n + f d ( M n , u ( t ) ) f d ( M n , M n )
From the above discussion, the numerical implementation of the CPM can be obtained by combining the alternating series and the FORC, which can be used for evaluation of the CPM output; this approach avoids the trouble of differentiation and integration shown in Figure 4.

3. Realizations of the CPM

In the previous section, we introduced the basic definition and the numerical implementation on the CPM. We also developed f a ( m , M ) and f d ( M , m ) as an algebraic representation of a CPM. Because the f a ( m , M ) and f d ( M , m ) can be explicitly obtained by experiments of the FORC-measuring experiments, they play the central roles in the realization of a CPM. In practice, however, it is impossible to obtain all values of FORC which requires infinite and uncountable experiments. One realizable method to approximate values of FORC is to find finite data of f a ( m , M ) and f d ( M , m ) , store them as tables, and we can obtain values of f a ( m , M ) or f d ( M , m ) by the interpolation using its neighbors stored in tables. In this section, we take an introductory look at different realizations of a CPM using finite-FORC, including the table method [21,22], polynomial approximation, and adaptive polynomial approximation.
Linear interpolation is a very intuitive way to estimate values between known values and we adopt this method to approximate f a ( m , M ) or f d ( M , m ) . If a point ( M , m ) belongs to a rectangular cell formed by ( M 1 , m 1 ) , ( M 1 , m 2 ) , ( M 2 , m 1 ) , ( M 2 , m 2 ) , and we have f d ( M 1 , m 1 ) = f M 1 m 1 , f d ( M 1 , m 2 ) = f M 1 m 2 , f d ( M 2 , m 1 ) = f M 1 m 1 , f d ( M 2 , m 2 ) = f M 2 m 2 . Then, f d ( M , m ) is approximated by
f L = ( M 2 M ) f M 1 m 2 + ( M M 1 ) f M 2 m 2 M 2 M 1 , f R = ( M 2 M ) f M 1 m 1 + ( M M 1 ) f M 2 m 1 M 2 M 1 , f d ( M , m ) = ( m m 1 ) f L + ( m 2 m ) f R m 2 m 1 .
If a point ( M , m ) belongs to a triangular cell formed by ( M 1 , m 1 ) , ( M 2 , m 1 ) , ( M 2 , m 2 ) , then f d ( M , m ) is approximated by
f R = f M 2 m 2   , f L = ( M 2 M ) f M 1 m 2 + ( M M 1 ) f M 2 m 2 M 2 M 1 , f d ( M , m ) = ( m m 1 ) f L + ( m 2 m ) f R m 2 m 1 .
The Formulas (12) and (13) also exist for f a ( m , M ) .
It is convenient to investigate the measured data using an analyzable function called curve fitting and the polynomial function is commonly used. Our goal is to fit the samples of a finite number of FORC using polynomial functions. For instance, we can use a quadratic polynomial function of single variable M to fit f a ( m , M ) with a particular m . That is,
f a ( m , M ) = a 0 ( m ) + a 1 ( m ) M + a 1 ( m ) M 2
Additionally, f d ( M , m ) also can be approximated as
f d ( M , m ) = d 0 ( M ) + d 1 ( M ) m + d 2 ( M ) m 2
where a i ( m ) and d i ( m ) are notations to indicate the polynomial coefficients corresponding with m and M , respectively. The coefficients of the curve fitting method vary with m or M , so to know all of them is also impossible. Hence, we also use the interpolation to approximate coefficients that are not derived from experiments. For m k < m < m k + 1 ,
a i ( m ) = m k + 1 a i ( m k ) + m k a i ( m k + 1 ) m k + 1 m k ,   where   i = 0 , 1 , 2 .
For M k < M < M k + 1 ,
d i ( M ) = M k + 1 d i ( M k ) + M k d i ( M k + 1 ) M k + 1 M k
According to (10) and (11), a CPM output can be determined by a sifted term and a combination of FORC. Now, we substitute for FORC by using (14) and (15). If the input u ( t ) is ascending,
f ( t ) = f n + a 2 ( m n ) ( u 2 ( t ) m n 2 ) + a 1 ( m n ) ( u ( t ) m n )
Define a i ( p Δ ) = a s i [ p ] , i = 0 , 1 , 2 ; for m n = p Δ + r = u , 0 r < Δ , then
a i ( m n ) = ( 1 r Δ ) a s i [ p ] + r Δ a s i [ p + 1 ]
Thus, (18) becomes
f ( t ) = f n + θ 1 T x 1
where
θ 1 = [ a s 2 [ p ] a s 2 [ p + 1 ] a s 1 [ p ] a s 1 [ p + 1 ] ] , x 1 = ( u m n ) [ ( 1 r Δ ) ( u + m n ) r Δ ( u + m n ) ( 1 r Δ ) r Δ ]
Additionally, u indicates u ( t ) . Similarly, if the input u ( t ) is descending,
f ( t ) = F n + θ 2 T x 2
where
θ 2 = [ d s 2 [ p ] d s 2 [ p + 1 ] d s 1 [ p ] d s 1 [ p + 1 ] ] , x 2 = ( u M n ) [ ( 1 r Δ ) ( u + M n ) r Δ ( u + M n ) ( 1 r Δ ) r Δ ]
Additionally, M n = p Δ + r , 0 r < Δ . The polynomial coefficients represent the parameters of the CPM and are suitable for model identification because they are adjustable.
Based on the polynomial approximation, we obtain a set of coefficients that can be used to describe FORC. However, these coefficients may not be accurate. For a system model with inaccurate parameters, we can use adaptive algorithms such as LMS or RLS [18] to update the parameters. Compared with the adaptive model output and plant output, this error is used to update adaptive system parameters. Assume the adaptive model output f m = ( θ , x ) = θ T x is a function of parameters θ and the input x , and the plant output f p ( x ) is a function of the input x . Define the error as
e = f p ( x ) f m ( θ , x )
The steepest descent algorithm is:
θ ( k + 1 ) = θ ( k ) + 2 η e x
where η is a step size. Using (20), (22), (25), we obtain following updated formulations:
If the input u ( t ) is ascending,
θ 1 ( k + 1 ) = θ 1 ( k ) + 2 η e x 1 = [ a s 2 ( k ) [ p ] a s 2 ( k ) [ p + 1 ] a s 1 ( k ) [ p ] a s 1 ( k ) [ p + 1 ] ] + 2 η e ( u m n ) [ ( 1 r Δ ) ( u + m n ) r Δ ( u + m n ) ( 1 r Δ ) r Δ ]
If the input u ( t ) is descending,
θ 2 ( k + 1 ) = θ 2 ( k ) + 2 η e x 2 = [ d s 2 ( k ) [ p ] d s 2 ( k ) [ p + 1 ] d s 1 ( k ) [ p ] d s 1 ( k ) [ p + 1 ] ] + 2 η e ( u M n ) [ ( 1 r Δ ) ( u + M n ) r Δ ( u + M n ) ( 1 r Δ ) r Δ ]
Then, we can obtain CPM outputs by applying polynomial coefficients adapted by (26) or (27) to (20) or (22).

4. CPM Simulation and Discussion

In this section, we demonstrate the computer simulations to verify capabilities of different realizations of the CPM proposed in Section 3. Consider a special CPM with Preisach density function,
μ ( α , β ) = k e x p ( 1 2 { ( α α c σ 1 ) 2 + ( β β c σ 2 ) 2 } )
We chose k = 50 , α c = β c = 0 , σ 1 = 1.5 , σ 2 = 1.8 , and u s a t = 3 in our simulations. The exact CPM outputs are calculated by numerical integral of Matlab. The FORC are measured and sampled in the tables first; these sampled data are used to obtain polynomial coefficients. The CPM outputs are predicted via different realizations and the root mean square error (RMS) is adopted as the judgment of how close the prediction is to the CPM. We can also observe the tracking ability of the adaptive polynomial approximation by disturbing μ ( α , β ) slightly.
To implement the table method, the samples of the ascending FORC and descending FORC are measured first. By dividing 3 ~ + 3 into 6 segments, the CPM input is decreased from + 3 to some divided point and then increased to + 3 . After the input begins decreasing from some divided point to + 3 , we sample the CPM outputs every unit change in the CPM input and store these sampling data in Table 1 which presents the samples of ascending FORC. The first column denotes the value of m of f a ( m , M ) , and the first row denotes the value of M . The samples of descending FORC are stored in Table 2 using a similar procedure and both ascending and descending FORC are plotted in Figure 5.
We apply the least-square approximation to sampling data of Table 1 and Table 2 to identify the polynomial coefficients directly in this simulation. Because of the establishment of saturated boundaries, FORC are much more gradual and cannot be fit by polynomials adequately. Hence, f a ( m , 3 ) and f d ( M , 3 ) are rejected in the procedure of the least-square approximation and the table method is used if they are needed. The polynomial coefficients are stored in Table 3 and Table 4. Based on (20) and (22), the polynomial coefficients a 0 and d 0 are not needed for the polynomial approximation. The polynomial coefficients of f a ( m , M ) and f d ( M , m ) are plotted in Figure 6.
In our simulations, we use a sinusoidal CPM input that with a frequency of 1 Hz has a sample time of 0.01 s, and a varying amplitude. The waveform of CPM input and the input–output phase diagram are shown in Figure 7. The simulation results are shown in Figure 8, Figure 9, Figure 10 and Figure 11. The CPM output and CPM output realized via the table method are plotted together in Figure 8a. The difference between these two outputs, defined as the prediction error, are shown in Figure 8b. The RMS of the prediction error is 0.5256. The result of the simulation of the polynomial approximation is shown in Figure 9, and its RMS is 1.0717. We use a step size η = 0.005 and Table 3 and Table 4 as initial parameters of the adaptive polynomial approximation and the simulation result is shown in Figure 10. The RMS of the adaptive polynomial approximation is 0.5909. We also plotted these three prediction errors together in Figure 11 for comparison. We can find that the performance of the adaptive polynomial approximation and the table method are in the same range roughly and are better than the polynomial approximation. This is because polynomial approximation uses a single function to characterize a FORC but the table method uses a number of linear functions. The adaptive polynomial approximation, however, adjusts the coefficients of the polynomial so it can be viewed using several functions to fit FORC equivalently and has a better effect than the polynomial approximation.
We saw the performance of dissimilar realizations of the CPM with a fixed Preisach density function. That is, we assume that the quality of the CPM is abiding for all situations. This assumption is not cohere with systems in real applications. It is broad that the table method and the polynomial approximation are static models and hard to be modify real time. An advantage of the adaptive polynomial approximation is that the parameters of such scheme are obtained on-line; this means that even if the CPM changed, the adaptive polynomial approximation will trace the CPM by adjusting its parameters. To witness this, we disturb μ ( α , β ) such that σ 1 becomes from 1.5 to 1.3 and σ 2 becomes from 1.8 to 2.1. Table 1 and Table 2 are still used to implement the table method; Table 3 and Table 4 are still used to implement the polynomial approximation and are the initial parameters of the adaptive polynomial approximation. We apply the same CPM input of previous simulations and the new CPM output compared with the original CPM output is plotted in Figure 12. Then, we use the different realizations to predict the new CPM output and plot the prediction errors together in Figure 13. The RMS of each realization is 1.9707 for the table method, 1.0746 for the polynomial approximation, and 0.6349 for the adaptive polynomial approximation. It is clear that the adaptive polynomial approximation has the best performance in this simulation.

5. Applications for Piezoelectric Actuators

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.
To implement the Preisach model, several first order reversal curves are measured for the micro-piezoelectric actuator first. By dividing 0∼5 V into 10 segments, the input voltage of micro-piezoelectric actuator is increased from 0 V to some divided point and then decreased to 0 V. After the input begins decreasing from some divided point to 0, we sample the displacements of the micro-piezoelectric actuator every 0.45 V change in input voltage and store these sampling data in Table 5. The first column denotes the value of M, and the first row denotes the value of m. In the predication experiment, we use a sinusoidal and input voltage with a frequency of 1 rad/s, an amplitude of 2.25 V, and a bias of 2.25 V; the input and output waveforms of the micro-piezoelectric actuator of the predication experiment are shown in Figure 14, and the phase transition is shown in Figure 15. In the tracking control experiment, we use a sinusoidal desired output with with a frequency of 1 rad/s, an amplitude of 15 µm, and a bias of 16 µm. The sampling time of the experiment is 0.01 s and the total time is 50 s.
The displacement of the micro-piezoelectric actuator and the output of the Preisach model realized via the table method are plotted together in Figure 16a. The difference between these two outputs, defined as the model error, is shown in Figure 16b. The phase transition of model and of plant are plotted together in Figure 17. The RMS of the model error is 1.0098 µm, and the required amount of memory to store the table is 66.
We apply the least-square approximation to sampling numbers of Table 5 to identify the polynomial coefficients directly in this experiment. The polynomial coefficients are stored in Table 6. The polynomial coefficient c0 is not needed for the Preisach model realization. The displacement of the micro-piezoelectric actuator is compared with the output of the Preisach model which is realized via the table method as shown in Figure 18a, and the model error is shown in Figure 18b. The phase transition of the model and of the plant are plotted together in Figure 19. The RMS of the model error is 1.4817 µm, and the amount of memory for polynomial coefficients is 22.
We develop the least mean square (LMS) adaptive algorithm to obtain accuracy polynomial coefficients and the data in Table 6 are used as the initial parameters. We test the step sizes µ = 10 − 4, µ = 10 − 3, and µ = 0.00051 + u(t) and the process of each is shown in Figure 20. The polynomial coefficients obtained by the LMS adaptive algorithm with step size µ = 0.00051 + u(t) are stored in Table 7.
The displacement of the micro-piezoelectric actuator is compared with the output of the Preisach model, which is shown in Figure 21a, and the model error is shown in Figure 21b. The phase transitions of the model and the plant are plotted together in Figure 22. The RMS of the model error is 1.5898 µm.
Obviously, a bias exists in Figure 21b. This means that the Preisach model should be modified in the beginning of the micro-piezoelectric actuator work. Figure 23 shows the phase transition of the model error and the input in the beginning ascending branch.
We predict the displacement of the micro-piezoelectric actuator again. The displacement of the micro-piezoelectric actuator is compared with output of the Preisach model which is shown in Figure 21a, and the model error is shown in Figure 24b. The phase transition of the model and the plant are plotted together in Figure 25. The RMS of the model error is 0.4539 µm. The model error of the Preisach model realization via the table method and the modified polynomial approximation are plotted together in Figure 26 for comparison.
An experimental inverse Preisach model for tracking control of micro-piezoelectric actuators was developed. The desired output compared with the displacement of the micro-piezoelectric of this experiment is plotted together in Figure 27a. The difference between these two outputs, defined as the tracking error, is shown in Figure 27a. The RMS of the tracking error is 0.7079 µm. To reduce tracking error, we combine the inverse Preisach model with a PID controller. 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 actuator of this experiment are plotted together in Figure 28a. The tracking error is shown in Figure 28b. The RMS of tracking error is 0.2438 µm.

6. Conclusions

We implemented the CPM with a polynomial approximation to characterize the hysteresis. Compared with the traditional table method, this method requires less storage space and enables parameter tracking of the hysteresis element through experiments. We successfully obtained polynomial coefficients to model the CPM with a particular Preisach density function in Section 4, and the polynomial coefficients to model the displacement of a micro-piezoelectric actuator in Section 5. In Section 4, the obtained model compared with that via the table method not only requires less memory but also yields a smaller modeling RMS error from 1.0746 via the table method to 0.6349 as the CPM is changed. In Section 5, the modeling RMS error changed from 1.0098 µm via the table method to 0.4539 µm as CPM changes. We also established the inverse Preisach model based on the polynomial approximation; this model was combined with the PID controller used for tracking control and yielded small tracking RMS error of 0.2438 µm.

Author Contributions

All authors adopted the Preisach model to characterize hysteresis behaviors and proposed to use polynomial approximation to realize the Preisach model; this approach was then used to model the displacement of a real micro-piezoelectric actuator. The proposed approach can overcome the drawbacks of the conventional table method. Finally, all authors have given final approval of the version to be published and agree to be accountable for all aspects of the work in ensuring that questions related to the accuracy or integrity of any part of the work are appropriately investigated and resolved. All authors have read and agreed to the published version of the manuscript.

Funding

This research was sponsored by the Ministry of Science and Technology of Taiwan under the Grant MOST 110-2221-E-150-038.

Acknowledgments

This paper is in memory of Mu-Huo Cheng who was a professor in the Department of Electrical and Control Engineering, National Chiao Tung University, Hsinchu, Taiwan.

Conflicts of Interest

The authors declare no conflict of interest. Moreover, the funders had no role in the design of the study; in the collection, analyses, or interpretation of data; in the writing of the manuscript, or in the decision to publish the results.

Glossary

u(t)signal input
f(t)Preisach model output
μ(α,β)Preisach density function
γ α β [ u ( · ) ] ( t ) hysteresis operator
α(β)hysteresis operator corresponds to ‘up’ (‘down’) switching values of input
± u s a t divides the range of input data
u (τ)input signal
{ M 1 , M 2 , , M n 1 }peaks of the previous input signal
{ m 1 , m 2 , , m n 1 }valleys of the previous input signal
f α β (M,m)polynomial function
C i (M)indicates the coefficient corresponding with M
C 0 ( M ) ,   C 1 ( M ) ,   C 2 (M)polynomial coefficients
f n represents the effect of previous input
f p (t)polynomial approximation
e (t)output error
f(t)plant output
θ update algorithm for the vector
Ttranspose operation
X(t)function of the input
ηconstant step size
Nnormalization factor

References

  1. Hussain, S.; Lowther, D.A. An Efficient Implementation of the Classical Preisach Model. IEEE Trans. Magn. 2018, 54, 7300204. [Google Scholar] [CrossRef]
  2. MacKenzie, I.; Trumper, D.L. Real-Time Hysteresis Modeling of a Reluctance Actuator Using a Sheared-Hysteresis-Model Observer. IEEE/ASME Trans. Mechatron. 2016, 21, 4–16. [Google Scholar] [CrossRef]
  3. Liang, H.; Zeng, C.; Giraud-Audine, C. Modelling and Simulation of Shape Memory Alloys Micro-Actuator. In Proceedings of the 2009 Fourth International Conference on Innovative Computing, Information and Control (ICICIC), Kaohsiung, Taiwan, 7–9 December 2009; pp. 1405–1408. [Google Scholar]
  4. Jia, W.; Luo, X.; Chen, X. Modeling and analysis of hysteresis in piezoceramic actuators excited by dynamic driving frequency. In Proceedings of the 2009 Symposium on Piezoelectricity, Acoustic Waves, and Device Applications (SPAWDA 2009), Wuhan, China, 17–20 December 2009; p. 116. [Google Scholar]
  5. Agrawal, B.N.; Elshafei, M.A.; Song, G. Adaptive antenna shape control using piezoelectric actuators. Act. Astronaut. 1997, 40, 821–826. [Google Scholar] [CrossRef] [Green Version]
  6. Mayhan, P.; Srinivasan, K.; Watechagit, S.; Washington, G. Dynamic modeling and controller design for a piezoelectric actuation system used for machine tool control. J. Intell. Mater. Syst. Struct. 2000, 11, 771–780. [Google Scholar] [CrossRef]
  7. Chen, X.; Su, C.; Fukuda, T. Adaptive Control for the Systems Preceded by Hysteresis. IEEE Trans. Autom. Control 2008, 53, 1019–1025. [Google Scholar] [CrossRef]
  8. Tao, G.; Kokotovic, P.V. Adaptive control of plants with unknown hystereses. IEEE Trans. Autom. Control 1995, 40, 200–212. [Google Scholar] [CrossRef]
  9. Zhang, J.; Merced, E.; Sepúlveda, N.; Tan, X. Modeling and Inverse Compensation of Nonmonotonic Hysteresis in VO2-Coated Microactuators. IEEE/ASME Trans. Mechatron. 2014, 19, 579–588. [Google Scholar] [CrossRef]
  10. Li, Z.; Zhang, X.; Su, C.-Y.; Chai, T. Nonlinear Control of Systems Preceded by Preisach Hysteresis Description: A Prescribed Adaptive Control Approach. IEEE Trans. Control Syst. Technol. 2016, 24, 451–460. [Google Scholar] [CrossRef]
  11. Riccardi, L.; Naso, D.; Turchiano, B.; Janocha, H. Adaptive Control of Positioning Systems with Hysteresis Based on Magnetic Shape Memory Alloys. IEEE Trans. Control Syst. Technol. 2013, 21, 2011–2023. [Google Scholar] [CrossRef]
  12. Stoleriu, L.; Stancu, A. Using Experimental FORC Distribution as Input for a Preisach-Type Model. IEEE Trans. Magn. 2006, 42, 3159–3161. [Google Scholar] [CrossRef]
  13. Ge, P.; Jouaneh, M. Modeling hysteresis in piezoceramic actuators. Prec. Eng. 2005, 17, 211–221. [Google Scholar] [CrossRef]
  14. Liu, V.-T.; Lin, C.-L.; Wing, H.-Y. New realisation of Preisach model using adaptive polynomial approximation. Int. J. Syst. Sci. 2012, 43, 1642–1649. [Google Scholar] [CrossRef]
  15. Liu, V.-T.; Lin, C.-L.; Wing, H.-Y. Realization of Preisach model using adaptive polynomial approximation method. In Proceedings of the 2010 5th IEEE Conference on Industrial Electronics and Applications, Taichung, Taiwan, 15–17 June 2010; pp. 1258–1263. [Google Scholar]
  16. Mayergoyz, I.D. Mathematical Models of Hysteresis; Springer: New York, NY, USA, 1991. [Google Scholar]
  17. Mayergoyz, I.D. Hysteresis models from the mathematical and control theory points of view. J. Appl. Phys. 1985, 57, 3803–3805. [Google Scholar] [CrossRef]
  18. Choi, S.B.; Han, Y.M. Hysteretic behavior of a magnetorheological fluid: Experimental identification. Acta Mech. 2005, 180, 37–47. [Google Scholar] [CrossRef]
  19. Iyer, R.V.; Shirley, M.E. Hysteresis Parameter Identification with Limited Experimental Data. IEEE Trans. Magn. 2004, 40, 3227–3239. [Google Scholar] [CrossRef]
  20. Tan, X.; Baras, J.S. Adaptive Identification and Control of Hysteresis in Smart Materials. IEEE Trans. Autom. Control 2005, 50, 827–839. [Google Scholar]
  21. Ge, P.; Jouaneh, M. Tracking control of a piezoceramic actuator. IEEE Trans. Control Syst. Technol. 1996, 4, 209–216. [Google Scholar]
  22. Tang, Z.-H.; Lv, F.-U.; Xiang, Z.-H. Hysteresis model of magnetostrictive actuators and its numerical realization. J. Zhejiang Univ. Sci. A 2007, 8, 1059–1064. [Google Scholar] [CrossRef]
Figure 1. Architecture of the Preisach model.
Figure 1. Architecture of the Preisach model.
Symmetry 14 01008 g001
Figure 2. The hysteresis operator: the output follows trajectory A B C D when the input is increasing and follows E F G H when the input is decreasing.
Figure 2. The hysteresis operator: the output follows trajectory A B C D when the input is increasing and follows E F G H when the input is decreasing.
Symmetry 14 01008 g002
Figure 3. A progressive evaluation of the CPM: (a) the input diagram; (b) the input–output phase diagram.
Figure 3. A progressive evaluation of the CPM: (a) the input diagram; (b) the input–output phase diagram.
Symmetry 14 01008 g003
Figure 4. Block diagram of the numerical implementation of the CPM.
Figure 4. Block diagram of the numerical implementation of the CPM.
Symmetry 14 01008 g004
Figure 5. FORC of the CPM: (a) ascending FORC; (b) descending FORC.
Figure 5. FORC of the CPM: (a) ascending FORC; (b) descending FORC.
Symmetry 14 01008 g005
Figure 6. The coefficients of the PA: (a) of f a ( m , M ) ; (b) of f d ( M , m ) .
Figure 6. The coefficients of the PA: (a) of f a ( m , M ) ; (b) of f d ( M , m ) .
Symmetry 14 01008 g006
Figure 7. The CPM input and output of simulations: (a) the waveform of the CPM input; (b) the input–output phase diagram.
Figure 7. The CPM input and output of simulations: (a) the waveform of the CPM input; (b) the input–output phase diagram.
Symmetry 14 01008 g007
Figure 8. Prediction of the CPM output using the table method: (a) the CPM output compared with prediction via the table method; (b) the prediction error of the table method.
Figure 8. Prediction of the CPM output using the table method: (a) the CPM output compared with prediction via the table method; (b) the prediction error of the table method.
Symmetry 14 01008 g008
Figure 9. Prediction of the CPM output using PA: (a) the CPM output compared with prediction via PA; (b) the prediction error of PA.
Figure 9. Prediction of the CPM output using PA: (a) the CPM output compared with prediction via PA; (b) the prediction error of PA.
Symmetry 14 01008 g009
Figure 10. Prediction of the CPM output using the adaptive PA: (a) the CPM output compared with prediction; (b) the prediction error.
Figure 10. Prediction of the CPM output using the adaptive PA: (a) the CPM output compared with prediction; (b) the prediction error.
Symmetry 14 01008 g010
Figure 11. Prediction error of the CPM via different realizations.
Figure 11. Prediction error of the CPM via different realizations.
Symmetry 14 01008 g011
Figure 12. The CPM output of a changed Preisach density function: (a) the new CPM output compared with the original CPM output; (b) the difference between the new and original CPM output.
Figure 12. The CPM output of a changed Preisach density function: (a) the new CPM output compared with the original CPM output; (b) the difference between the new and original CPM output.
Symmetry 14 01008 g012
Figure 13. Prediction error of the CPM of a changed Preisach density function via different realizations.
Figure 13. Prediction error of the CPM of a changed Preisach density function via different realizations.
Symmetry 14 01008 g013
Figure 14. The input and output of the micro-piezoelectric actuator.
Figure 14. The input and output of the micro-piezoelectric actuator.
Symmetry 14 01008 g014
Figure 15. The phase transition of micro-piezoelectric actuator.
Figure 15. The phase transition of micro-piezoelectric actuator.
Symmetry 14 01008 g015
Figure 16. Prediction of the micro-piezoelectric actuator using the Preisach model realization via the table method: (a) the displacement of the micro-piezoelectric actuator and the output of the Preisach model; (b) the error between these two outputs.
Figure 16. Prediction of the micro-piezoelectric actuator using the Preisach model realization via the table method: (a) the displacement of the micro-piezoelectric actuator and the output of the Preisach model; (b) the error between these two outputs.
Symmetry 14 01008 g016
Figure 17. Phase transition of the micro-piezoelectric actuator using the Preisach model realization via the table method.
Figure 17. Phase transition of the micro-piezoelectric actuator using the Preisach model realization via the table method.
Symmetry 14 01008 g017
Figure 18. Prediction of the micro-piezoelectric actuator using the Preisach model realization via polynomial approximation: (a) the displacement of the micro-piezoelectric actuator and the output of the Preisach model; (b) the error between these two outputs.
Figure 18. Prediction of the micro-piezoelectric actuator using the Preisach model realization via polynomial approximation: (a) the displacement of the micro-piezoelectric actuator and the output of the Preisach model; (b) the error between these two outputs.
Symmetry 14 01008 g018
Figure 19. Phase transition of the micro-piezoelectric actuator using the Preisach model realization via polynomial approximation.
Figure 19. Phase transition of the micro-piezoelectric actuator using the Preisach model realization via polynomial approximation.
Symmetry 14 01008 g019
Figure 20. Mean square error of the prediction of the micro-piezoelectric actuator for varying step-size µ using the adaptive polynomial approximation.
Figure 20. Mean square error of the prediction of the micro-piezoelectric actuator for varying step-size µ using the adaptive polynomial approximation.
Symmetry 14 01008 g020
Figure 21. Prediction of the micro-piezoelectric actuator using the Preisach model realization via the adaptive polynomial approximation: (a) the displacement of the micro-piezoelectric actuator is compared with the output of the Preisach model; (b) the model error.
Figure 21. Prediction of the micro-piezoelectric actuator using the Preisach model realization via the adaptive polynomial approximation: (a) the displacement of the micro-piezoelectric actuator is compared with the output of the Preisach model; (b) the model error.
Symmetry 14 01008 g021
Figure 22. Phase transition of the micro-piezoelectric actuator using the Preisach model realization via the adaptive polynomial approximation.
Figure 22. Phase transition of the micro-piezoelectric actuator using the Preisach model realization via the adaptive polynomial approximation.
Symmetry 14 01008 g022
Figure 23. Phase transition of the model error and the input.
Figure 23. Phase transition of the model error and the input.
Symmetry 14 01008 g023
Figure 24. Prediction of the micro-piezoelectric actuator using the Preisach model realization via the modified polynomial approximation: (a) the displacement of the micro-piezoelectric actuator is compared with output of the Preisach model; (b) the model error.
Figure 24. Prediction of the micro-piezoelectric actuator using the Preisach model realization via the modified polynomial approximation: (a) the displacement of the micro-piezoelectric actuator is compared with output of the Preisach model; (b) the model error.
Symmetry 14 01008 g024
Figure 25. Phase transition of the micro-piezoelectric actuator using the Preisach model realization via the modified polynomial approximation.
Figure 25. Phase transition of the micro-piezoelectric actuator using the Preisach model realization via the modified polynomial approximation.
Symmetry 14 01008 g025
Figure 26. Model error of the Preisach model realization via the table method and the modified polynomial approximation.
Figure 26. Model error of the Preisach model realization via the table method and the modified polynomial approximation.
Symmetry 14 01008 g026
Figure 27. Tracking of the micro-piezoelectric actuator using the inverse Preisach model: (a) the desired output compared with the displacement of the micro-piezoelectric; (b) the difference between these two outputs.
Figure 27. Tracking of the micro-piezoelectric actuator using the inverse Preisach model: (a) the desired output compared with the displacement of the micro-piezoelectric; (b) the difference between these two outputs.
Symmetry 14 01008 g027
Figure 28. Tracking of the micro-piezoelectric actuator using the inverse Preisach model with a PID controller: (a) the desired output compared with the displacement of the micro-piezoelectric actuator; (b) the tracking error.
Figure 28. Tracking of the micro-piezoelectric actuator using the inverse Preisach model with a PID controller: (a) the desired output compared with the displacement of the micro-piezoelectric actuator; (b) the tracking error.
Symmetry 14 01008 g028
Table 1. Sampling data of ascending FORC.
Table 1. Sampling data of ascending FORC.
m / M −3−2−10123
−3−21.58−13.70−1.0911.6019.0221.3121.58
−2 −13.70−12.40−5.905.7516.2921.58
−1 −1.101.549.3417.3521.58
0 11.6014.1818.0821.58
1 19.0220.2421.58
2 21.3121.58
3 21.58
Table 2. Sampling data of descending FORC.
Table 2. Sampling data of descending FORC.
M / m −3−2−10123
−3−21.58
−2−21.58−21.28
−1−21.58−19.89−18.59
0−21.58−17.78−12.62−9.98
1−21.58−15.67−6.641.213.79
2−21.58−14.29−2.758.5014.4815.70
3−21.58−13.70−1.0911.6019.0221.3121.58
Table 3. Polynomial coefficients of f a ( m , M ) .
Table 3. Polynomial coefficients of f a ( m , M ) .
m a 2 a 1 a 0
−31.839.55−9.67
−21.697.81−5.37
−11.344.972.28
01.021.5611.60
10−1.2221.46
20021.31
Table 4. Polynomial coefficients of f d ( M , m ) .
Table 4. Polynomial coefficients of f d ( M , m ) .
M d 2 d 1 d 0
−200−21.28
−101.30−17.30
0−1.261.38−9.98
1−1.615.010.59
2−1.857.728.03
3−1.849.0211.12
Table 5. Sampling data of FORC.
Table 5. Sampling data of FORC.
M/m00.450.901.351.802.252.703.153.604.054.50
00
0.450.023.43
0.900.094.077.73
1.350.115.289.3812.76
1.800.145.3610.1314.1118.15
2.250.235.6870.7915.3219.3423.32
2.700.285.9411.3116.2220.5824.4828.03
3.150.366.1711.7216.6821.2725.4129.0132.65
3.600.456.3512.0917.2921.9626.2730.1333.5136.49
4.050.486.6112.3517.5822.4026.8130.7434.3537.7440.15
4.500.676.8712.8718.0722.9827.3331.3834.9338.1340.9643.18
Table 6. The polynomial coefficients.
Table 6. The polynomial coefficients.
M c 2 c 1 c 0
0000
0.4507.57780.02
0.90−0.790129.20.09
1.35−2.209912.3280.1275
1.80−1.111111.9490.174
2.25−1.032612.5430.2525
2.70−1.11713.3020.24619
3.15−1.025613.4390.37625
3.60−1.09113.9540.39576
4.05−1.05814.0750.49745
4.50−1.075214.2690.73091
Table 7. The polynomial coefficients obtained using the LMS adaptive algorithm.
Table 7. The polynomial coefficients obtained using the LMS adaptive algorithm.
M c 2 c 1
00.57183.8438
0.45−2.73237.3383
0.90−2.33259.2038
1.35−2.162610.8419
1.80−1.770011.8186
2.25−1.458812.3549
2.70−1.154012.4685
3.15−1.212513.2395
3.60−1.256513.9011
4.05−1.179714.0324
4.50−1.093813.9965
Publisher’s Note: MDPI stays neutral with regard to jurisdictional claims in published maps and institutional affiliations.

Share and Cite

MDPI and ACS Style

Liu, V.-T.; Wing, H.-Y. Classical Preisach Model Based on Polynomial Approximation and Applied to Micro-Piezoelectric Actuators. Symmetry 2022, 14, 1008. https://doi.org/10.3390/sym14051008

AMA Style

Liu V-T, Wing H-Y. Classical Preisach Model Based on Polynomial Approximation and Applied to Micro-Piezoelectric Actuators. Symmetry. 2022; 14(5):1008. https://doi.org/10.3390/sym14051008

Chicago/Turabian Style

Liu, Van-Tsai, and Home-Young Wing. 2022. "Classical Preisach Model Based on Polynomial Approximation and Applied to Micro-Piezoelectric Actuators" Symmetry 14, no. 5: 1008. https://doi.org/10.3390/sym14051008

APA Style

Liu, V. -T., & Wing, H. -Y. (2022). Classical Preisach Model Based on Polynomial Approximation and Applied to Micro-Piezoelectric Actuators. Symmetry, 14(5), 1008. https://doi.org/10.3390/sym14051008

Note that from the first issue of 2016, this journal uses article numbers instead of page numbers. See further details here.

Article Metrics

Back to TopTop