**3. Compensation algorithms**

The definition of suitable control algorithms to be embedded in actuators employing smart materials is sensibly affected by the memory phenomena involved. This because, roughly speaking, the controller has no memory of past history which, conversely, affects the plant behavior, and will be not able to provide suitable correction to track a desired trajectory. The system therefore will show poor performances. A simple way to overcome this problems, proposing quite a simple control strategy, is based on the following well known issues. As a preliminary step, let us assume that the plant is adequately modeled by the series connection of a hysteresis operator and a linear dynamic component. In Fig. 8-(a) it is shown a control scheme employing a generic nonlinear controller which gather both dynamics and *rate-independent* properties, even if specifying its details is non-trivial. Let us assume that the *rate independent* behavior of the plant can be modeled by the aid of a hysteresis operator *W*[*x*], as defined in section 2. If the inverse of the hysteresis operator partly modeling the plant were available, their series connection in the control scheme bears to a linear system, allowing the employment of all standard linear control strategies, as shown in Fig. 8-(b). In summary, in order to simplify the controller design, it is suitable to compensate nonlinearity and *rate-independent* memory of the actuator. This approach is generally referred to as *model-based* control strategy. The issue is now to address the definition of the compensator of a hysteresis

(a) The general controller should take into account memory effects (b) Model-based control scheme: a suitable modeling of the plant is performed which enable to provide a compensator

**Figure 8.** General control scheme (a) - Model-based control scheme (b)

operator and to define the conditions under which such an operator admits an inverse or *compensator*,[7, 27]. In particular:

**Definition 1.** *A hysteresis operator* <sup>W</sup>−<sup>1</sup> *with initial state <sup>φ</sup>*−<sup>1</sup> <sup>∈</sup> <sup>Ψ</sup>0*, is called a* compensator *(or inverse) of the operator* W*, with initial state ψ*−<sup>1</sup> *if, for any state φ* ∈ Ψ0*, there is a state ψ* ∈ Ψ0*, such that* <sup>W</sup>*<sup>ψ</sup>* ◦ W−<sup>1</sup> *<sup>φ</sup> <sup>x</sup>*(*t*) = <sup>W</sup>−<sup>1</sup> *<sup>φ</sup>* ◦ W*<sup>ψ</sup> x*(*t*) = *x*(*t*)*, for every input function x(t).*

Of course, the above definition lead us asking the conditions making this issue sound. This requires to specify the hysteresis operator with more detail. To this aim, we refer to the class of operators based on the Preisach deletion rule, stated in eqn. (5). We refer to this operators, to as P0-operators. In particular, we can define a piece-wise strictly monotone operator of this class, [7], according to the following

**Definition 2.** *The operator* W = ◦ <sup>Θ</sup> *with initial state <sup>ψ</sup>*−<sup>1</sup> *mapping the set of continous functions,* C[0, *T*]*, into itself, is piece–wise strictly increasing if*

$$\left[F\left(\Theta(\mathbf{x},\boldsymbol{\psi}) - F\left(\boldsymbol{\psi}\right)\right)\left[\mathbf{x} - \boldsymbol{\psi}(\mathbf{0})\right] > 0, \quad \mathbf{x} \neq \boldsymbol{\psi}(\mathbf{0}), \tag{13}$$

*where ψ*(0) *is the past input, while the symbol* Θ(*x*, *ψ*) *is the staircase attained by starting from ψ and applying the value of input x at time t.*

8 Will-be-set-by-IN-TECH

The definition of suitable control algorithms to be embedded in actuators employing smart materials is sensibly affected by the memory phenomena involved. This because, roughly speaking, the controller has no memory of past history which, conversely, affects the plant behavior, and will be not able to provide suitable correction to track a desired trajectory. The system therefore will show poor performances. A simple way to overcome this problems, proposing quite a simple control strategy, is based on the following well known issues. As a preliminary step, let us assume that the plant is adequately modeled by the series connection of a hysteresis operator and a linear dynamic component. In Fig. 8-(a) it is shown a control scheme employing a generic nonlinear controller which gather both dynamics and *rate-independent* properties, even if specifying its details is non-trivial. Let us assume that the *rate independent* behavior of the plant can be modeled by the aid of a hysteresis operator *W*[*x*], as defined in section 2. If the inverse of the hysteresis operator partly modeling the plant were available, their series connection in the control scheme bears to a linear system, allowing the employment of all standard linear control strategies, as shown in Fig. 8-(b). In summary, in order to simplify the controller design, it is suitable to compensate nonlinearity and *rate-independent* memory of the actuator. This approach is generally referred to as *model-based* control strategy. The issue is now to address the definition of the compensator of a hysteresis

> (b) Model-based control scheme: a suitable modeling of the plant is performed which enable to provide a

compensator

operator and to define the conditions under which such an operator admits an inverse or

**Definition 1.** *A hysteresis operator* <sup>W</sup>−<sup>1</sup> *with initial state <sup>φ</sup>*−<sup>1</sup> <sup>∈</sup> <sup>Ψ</sup>0*, is called a* compensator *(or inverse) of the operator* W*, with initial state ψ*−<sup>1</sup> *if, for any state φ* ∈ Ψ0*, there is a state ψ* ∈ Ψ0*, such*

Of course, the above definition lead us asking the conditions making this issue sound. This requires to specify the hysteresis operator with more detail. To this aim, we refer to the class of operators based on the Preisach deletion rule, stated in eqn. (5). We refer to this operators, to as P0-operators. In particular, we can define a piece-wise strictly monotone operator of this

**Definition 2.** *The operator* W = ◦ <sup>Θ</sup> *with initial state <sup>ψ</sup>*−<sup>1</sup> *mapping the set of continous functions,*

[(Θ(*x*, *ψ*) − (*ψ*)] [*x* − *ψ*(0)] > 0, *x* �= *ψ*(0), (13)

*<sup>φ</sup>* ◦ W*<sup>ψ</sup> x*(*t*) = *x*(*t*)*, for every input function x(t).*

**3. Compensation algorithms**

(a) The general controller should take into account

**Figure 8.** General control scheme (a) - Model-based control scheme (b)

memory effects

*that* <sup>W</sup>*<sup>ψ</sup>* ◦ W−<sup>1</sup>

*compensator*,[7, 27]. In particular:

*<sup>φ</sup> <sup>x</sup>*(*t*) = <sup>W</sup>−<sup>1</sup>

class, [7], according to the following

C[0, *T*]*, into itself, is piece–wise strictly increasing if*

Such a property guarantees the invertibility of P0-operators and is connected to the specific memoryless mapping functional . It turns to the classical definition of monotonicity for the (single valued) hysteresis branch when the state (i.e. the staircase in the Preisach plane) is specified. In the particular case of the classical Preisach operator, with mapping functional specified in eqn. (7), the following theorem holds [6]:

**Theorem 1.** *The Preisach operator,* H∈P0*, in eqn. 3 with distribution function μ, admits an inverse if and only if <sup>μ</sup>* <sup>≥</sup> <sup>0</sup>*, and E*(*u*, *<sup>v</sup>*) <sup>≥</sup> <sup>0</sup>*, for every* (*u*, *<sup>v</sup>*) <sup>∈</sup> **<sup>R</sup>**<sup>+</sup> <sup>×</sup> **<sup>R</sup>***.*

Before discussing an efficient algorithm to implement the inverse of Preisach operator, it is important asking why the inverse operator cannot be employed by just exchanging the role of input and output.1 In principle, this is not accurate, due to a specific property of the Preisach mapping functional, known as *congruency*. Let us refer to Fig. 9-(a), where it is sketched the state updating corresponding to the same *up and forth* input variation (from *x*<sup>−</sup> to *x*+) for two different states *h*<sup>1</sup> and *h*2. It is quite evident that two minor loops are generated, each corresponding to an initial state, and they are *congruent* by vertical translation, Fig. 9-(b). The inverse operator describe a characteristic which is placed symmetrically to the direct characteristic, with respect to the *y* = *x* line. As a consequence, any minor loop described by the model is transformed into one which is symmetrically placed with respect the same line.This implies that two congruent minor loops generated by the model, are transformed by the inverse operator into two minor loops which are congruent by a horizontal translation. Fig. 10-(a) illustrates such a situation. Conversely the operator which just exchanges the

(a) States corresponding to the same up and forth input variation

#### **Figure 9.** Illustration of the Congruency property of the Classical Preisach Problem

role of the variables that we call *pseudo-compensator* is still a Preisach operator and therefore generates vertically congruent minor loops, as shown in Fig. 10(b). It is so manifest that this operator cannot be strictly considered as the inverse of a Preisach operator, due to the limitations imposed by the congruency property. However this observation doesn't limit the possibility to adopt pseudo-compensator as a hysteresis model by itself which is able to

<sup>1</sup> In reality this is a procedure to define and approximation of the inverse, referred to as *pseudo-compensator* as detailed in [30]

**Figure 10.** Illustration of the Congruency property of the Classical Preisach Problem

guarantee good performances, as widely described in [20, 29]. In summary, the *inverse* of a classical Preisach operator *is not* a Preisach operator; the pseudo-compensator, conversely, *is* a Preisach operator and cannot rigorously represent the Preisach compensator.

It should be stressed that even if the inverse of the Preisach operator exists, it is not possible, generally speaking, to represent the compensator in a closed form, except than in very specific cases, [44], [26], [27], [43]. So it is mandatory the definition of well-defined and computationally efficient inversion algorithms in general working conditions. To this aim we will briefly sketch the issue, by referring to [12, 44] for details. When the conditions of invertibility of the Preisach operator are fulfilled the problem of finding its inverse has sense and the inversion procedure is based on the algorithm of state updating which assumed the output *y*(*t*) as known. In other words, the algorithm described hereafter, allows to provide the input field of the operator at time instant *t*, corresponding to a given and known output *y*(*t*), once the initial state is assigned. This procedure, since no numerical iterations are needed, allows to compute the inverse with the same computational efficiency and accuracy of the classical algorithm used to compute the Preisach operator. For this reason we refer to it as a "fast" inversion algorithm, [13] [14]. If the Preisach plane is divided into a Δ*x* × Δ*x* mesh, with Δ*x* = *xmax*−*xmin <sup>N</sup>* then the Everett function is uniformly sampled into a (*N* + 1) × *N*/2 mesh, as sketched in the Fig. 11. Therefore the solution of the above equation can be performed by a simple inspection of a look-up table, according to the rules described below. Moreover, let us assume that the Everett integral is limited to the triangle *ABC*. This is a quite reasonable hypothesis, since the function must vanish at infinity. In the same plane, it can be evidenced the initial state of the operator, *h*0(*r*, *t*). Let us further assume the output *y*(*t*) = *y*<sup>0</sup> + Δ*y*, is known. In the ordinary fashion of handling the Preisach operator, the output variation Δ*y* coincides with the integral of the distribution function over the area *S*, that is:

$$\int\_{S} \int\_{S} \mu(\mathbf{r}, \mathbf{s}) d\mathbf{r} \, d\mathbf{s} = \mathcal{E}(\mathbf{r}, \mathbf{s}) - \iint\_{S^{\*}} \mu(\mathbf{r}, \mathbf{s}) d\mathbf{r} \, d\mathbf{s}.\tag{14}$$

By exploiting the wiping out property, all the vertex of the staircase such that the constraint

$$
\Delta y \ge \mathcal{E}\_k(r, s) - \iint\_{S^\*} \mu(r, s) dr \, ds. \tag{15}
$$

**Figure 11.** The algorithm of inversion of the Preisach operator

10 Will-be-set-by-IN-TECH


(b) Congruent minor loops

**Figure 10.** Illustration of the Congruency property of the Classical Preisach Problem

Preisach operator and cannot rigorously represent the Preisach compensator.

guarantee good performances, as widely described in [20, 29]. In summary, the *inverse* of a classical Preisach operator *is not* a Preisach operator; the pseudo-compensator, conversely, *is* a

It should be stressed that even if the inverse of the Preisach operator exists, it is not possible, generally speaking, to represent the compensator in a closed form, except than in very specific cases, [44], [26], [27], [43]. So it is mandatory the definition of well-defined and computationally efficient inversion algorithms in general working conditions. To this aim we will briefly sketch the issue, by referring to [12, 44] for details. When the conditions of invertibility of the Preisach operator are fulfilled the problem of finding its inverse has sense and the inversion procedure is based on the algorithm of state updating which assumed the output *y*(*t*) as known. In other words, the algorithm described hereafter, allows to provide the input field of the operator at time instant *t*, corresponding to a given and known output *y*(*t*), once the initial state is assigned. This procedure, since no numerical iterations are needed, allows to compute the inverse with the same computational efficiency and accuracy of the classical algorithm used to compute the Preisach operator. For this reason we refer to it as a "fast" inversion algorithm, [13] [14]. If the Preisach plane is divided into a Δ*x* × Δ*x* mesh, with

*<sup>N</sup>* then the Everett function is uniformly sampled into a (*N* + 1) × *N*/2 mesh,

*<sup>S</sup>*<sup>∗</sup> *<sup>μ</sup>*(*r*,*s*)*dr ds*. (14)

*<sup>S</sup>*<sup>∗</sup> *<sup>μ</sup>*(*r*,*s*)*dr ds*. (15)

as sketched in the Fig. 11. Therefore the solution of the above equation can be performed by a simple inspection of a look-up table, according to the rules described below. Moreover, let us assume that the Everett integral is limited to the triangle *ABC*. This is a quite reasonable hypothesis, since the function must vanish at infinity. In the same plane, it can be evidenced the initial state of the operator, *h*0(*r*, *t*). Let us further assume the output *y*(*t*) = *y*<sup>0</sup> + Δ*y*, is known. In the ordinary fashion of handling the Preisach operator, the output variation Δ*y*

coincides with the integral of the distribution function over the area *S*, that is:

Δ*y* ≥ E*k*(*r*,*s*) −

*μ*(*r*,*s*)*dr ds* = E(*r*,*s*) −

By exploiting the wiping out property, all the vertex of the staircase such that the constraint

 *S* -

(a) States corresponding to the same up and

forth input variation

Δ*x* = *xmax*−*xmin*

is fulfilled, are cancelled. The suffix *k* of he Everett function evidence the vertex of the triangle. In the case reported in Fig.11-(a) the initial state *h*<sup>0</sup> represented by the sequence of dominant extrema {0, 1, 2, 3, 4, 5, 6, 7} is updated and the vertex {0, 1, 2, 3, 4} cancelled. So the state updating algorithm consists in applying the above condition and cancel all the vertex until it is verified. Once the last vertex has been cancelled, the following condition hold:

the grid points

$$
\Delta y = \mathcal{E}(\mathbf{x}, v\_0) - \iint\_{S^\*} \mu(r, s) dr \, ds. \tag{16}
$$

being *x* the unknown value of the input field to be determined. In Fig. 11-(b) it is shown how *x*(*t*) wipes out the vertex of the state function, and the operator attains the new state *h*(*r*, *t*) = {*x*, *P*, 7, *C*} evidenced in red in the same figure where it is also shown a zoom of the new state around the point *P*(*x*, *v*0) lying within a square determined by the points (*n* + 1, *k*), (*n* + 1, *k* + 1), (*n*, *k*) and (*n*, *k* + 1). As a first step, the Everett integrals at the points (*k* + 1, *v*0) and (*k*, *v*0) are found by linear interpolation:

$$\mathcal{E}(m, v\_0) = \mathcal{E}(m, n) - \left(\mathcal{E}(m, n) - \mathcal{E}(m, n + 1)\right) \frac{v\_0 - n\Delta x'}{\Delta x'} \tag{17}$$

being *m* = *k* + 1 and *m* = *k*. Finally, once these values have been determined, and exploiting eqn.(16) the final output value can be written as:

$$\mathbf{x}' = \frac{\Delta y + \iint \mu(r, s) dr \, ds - \mathcal{E}(k, v\_0)}{\mathcal{E}(k + 1, v\_0) - \mathcal{E}(k, v\_0)} \Delta \mathbf{x}' + k \Delta \mathbf{x}',\tag{18}$$

where *x*� = *x*/ (2) and Δ*x*� = Δ*x*/ (2). For this reason the equation does not change when the variables without index are concerned. Moreover, it should be stressed that the integral in the previous equation is a known quantity. Such an approach requires only three linear interpolations at most and does not need any lengthy iteration. The procedure only requires the knowledge of the matrix storing the samples of the Everett integrals and its precision is linked to the mesh adopted for the Preisach plane or, in other words, to the number of measured samples used for the identification. Moreover, the algorithm, as already mentioned, does not increase the computational weight of the algorithm which is equivalent to that of the direct operator.

## **4. Multi-variate systems with hysteresis**

Smart materials, in general application conditions, cannot be considered as SISO systems, since each output quantity is a function of at least two input functions. For sake of example, in a magnetostrictive material the strain depends both on the mechanical load (stress) and magnetic field. A similar behavior can be experienced for the magnetic *flux density*. The behavior is illustrated in Fig. 12-(a) where the strain is affected by both the stress and field variation. Moreover, not only the *H*-field affects the internal state of the system, which is also altered by the applied stress, as shown by the hysteresis branches sketched in Fig. 12-(b). The experimental evidence witnessed by Figs. 2 and 12, allows to conclude

**Figure 12.** Applied current and force vs. time to a magnetostrictive sample (a); displacement/force hysteresis branches (b).

that the modeling of these systems requires a generalization of all the classical approaches to hysteresis phenomena. This need started two decades ago [1], [46] and yielded to models able to couple together two different input variables (i.e. *ε* and *σ*). This effort become of fundamental importance whenever the coexistence of non stationary stress and fields cannot be overlooked.

The simplest way to define more general models is to keep the usual SISO structure of systems with hysteresis, as the ones discussed in the previous section. To this aim, let us consider the mapping:

$$\mathcal{L}: (\mu, \upsilon) \in \mathcal{C}^0[0, T] \times \mathcal{C}^0[0, T] \to \mathfrak{x} \in \mathcal{C}^0[0, T]. \tag{19}$$

where *ζ* is the operator relating the couple of continous functions (*u*, *v*), affecting the system behavior, to the input, *x*, of the memory operator, Θ. For the sake of example, in the case of a magnetostrictive material, *ζ* "links" the stress *σ* and magnetic field *H* to the input of the memory operator. Therefore

$$\mathcal{W} = \! \! \! \! \! / \, \circ \Theta \circ \! \! \! / \,, \tag{20}$$

represents a hysteresis operator with two inputs which generalizes the definition in [5] to a wider class of hysteresis operators. Before proceeding in the description of multi-variate systems with hysteresis in control applications, it is mandatory to provide some information concerning its identification. To this aim, let us now consider the following operator, specified for a magneto-elastic material

12 Will-be-set-by-IN-TECH

does not increase the computational weight of the algorithm which is equivalent to that of the

Smart materials, in general application conditions, cannot be considered as SISO systems, since each output quantity is a function of at least two input functions. For sake of example, in a magnetostrictive material the strain depends both on the mechanical load (stress) and magnetic field. A similar behavior can be experienced for the magnetic *flux density*. The behavior is illustrated in Fig. 12-(a) where the strain is affected by both the stress and field variation. Moreover, not only the *H*-field affects the internal state of the system, which is also altered by the applied stress, as shown by the hysteresis branches sketched in Fig. 12-(b). The experimental evidence witnessed by Figs. 2 and 12, allows to conclude


*<sup>ζ</sup>* : (*u*, *<sup>v</sup>*) ∈ C0[0, *<sup>T</sup>*] × C0[0, *<sup>T</sup>*] <sup>→</sup> *<sup>x</sup>* ∈ C0[0, *<sup>T</sup>*], (19)

W = ◦ Θ ◦ *ζ*, (20)

 

 

**Figure 12.** Applied current and force vs. time to a magnetostrictive sample (a); displacement/force

that the modeling of these systems requires a generalization of all the classical approaches to hysteresis phenomena. This need started two decades ago [1], [46] and yielded to models able to couple together two different input variables (i.e. *ε* and *σ*). This effort become of fundamental importance whenever the coexistence of non stationary stress and fields cannot

The simplest way to define more general models is to keep the usual SISO structure of systems with hysteresis, as the ones discussed in the previous section. To this aim, let us consider the

where *ζ* is the operator relating the couple of continous functions (*u*, *v*), affecting the system behavior, to the input, *x*, of the memory operator, Θ. For the sake of example, in the case of a magnetostrictive material, *ζ* "links" the stress *σ* and magnetic field *H* to the input of the

represents a hysteresis operator with two inputs which generalizes the definition in [5] to a wider class of hysteresis operators. Before proceeding in the description of multi-variate systems with hysteresis in control applications, it is mandatory to provide some information

 

(b)

 

direct operator.



hysteresis branches (b).

be overlooked.

memory operator. Therefore

mapping:

 

 

 

 

**4. Multi-variate systems with hysteresis**

 (a)

$$y = f \circ \Theta \circ \mathbb{Q}(\mathbf{x}\_1, \mathbf{x}\_2) + q(\mathbf{x}\_2), \tag{21}$$

where *y* = *ε*, *x*<sup>1</sup> = *H*, and *x*<sup>2</sup> = *σ*, being *ε*, *H*, and *σ* the strain, magnetic field and stress experienced by the material, respectively. Moreover, *q*(*σ*) is a pure elastic response to be identified with almost trivial mechanical measurements at zero magnetic field. Even it is not strictly necessary, the latter assumption allows to take into account separately the effects of vertical translation of the magnetostrictive response due to pure mechanical reasons. Let us preliminary assume that the set of *j* = 1, .., *M* anhysteretic2 curves, *ej* shown in Fig. 13-(b) is available. Moreover, assume that the interval of input variation (*x*<sup>2</sup> = *σ*) is [*σmin*, *σmax*] and that the magnetostrictive response experiences a monotonic decrease for increasing compressive stresses, that is *σmin* = *σ*∗. With reference to the figure, this corresponds to a 1.05 MPa applied mechanical stress. Finally, each curve is discretized into *N* samples such that *Hi*<sup>+</sup><sup>1</sup> − *Hi* = Δ*H*, with *i* = 1, ..*N*. Assuming further that *ζ*(*H*, *σ*∗) = *H* and *λ*(*σ*∗) = 0, if

**Figure 13.** Example of *Anhysteretic* curve (a) Experimental Anisteretic curves for a magnetostrictive material (b).

*σ* = *σ*∗, then the operator takes the following classical form

$$e(t) \equiv e(t) - q(\sigma) = \mathcal{F} \circ \Theta[H]\_{\prime} \tag{22}$$

which could be easily identified, according to the procedure sketched in sect. 2 and detailed in [28]. Once the Preisach distribution function is known, the procedure to identify the *ζ*(·, ·) function can start. Recalling now that in suitable conditions ◦ Θ admits an inverse and that an efficient inversion algorithm is available, as shown in section 3, the eq. (22) can easily be re-arranged as:

$$
\Gamma^{-1}[\varepsilon - q(\sigma)] = f(H, \sigma). \tag{23}
$$

By this procedure, a set of *M* × *N* measured samples is now available, which enable us finding the samples *ζ*(*Hi*, *σk*) of the unknown function.

<sup>2</sup> It is defined as *anhysteretic*, the curve in the input-output plane, as in Fig. 13-(a), each point of which corresponds to the output of the P<sup>0</sup> operator in the state represented by the horizontal curve *s* = *Hk* , that can be constructed by applying to the system the input *x*(*t*) slowly converging to *Hk*

**Figure 14.** Block scheme of the 2-variables operator (a); Example of reconstruction of the *ζ* function, named as "effective current" in the vertical axis of the plot (b).

In order to show the effectiveness of the described procedure, it is important now to show the performances of a multi-variable model as the one described in eqn. (23). In particular, in Fig. 14-(b) it is shown the reconstructed *ζ* function, while in Fig. 14-(a), the scheme of the model is shown. It undergoes the input variable *x*<sup>1</sup> and *x*2, according to eqn. (23), which in the specific case are the magnetic field and the stress, respectively. The output *y*∗(*t*) is compared to the actual measured strain function and the results are sketched in Fig. 15. In particular, in Fig. 15-(a) it is shown the time variation of the input variables, i.e. *x*1(*t*) = *H*(*t*) and *x*2(*t*) = *σ*(*t*), in the first two frames, while in the third, it is described the output of the *ζ*-function, to which the hysteresis operator W is subject. In Fig. 15-(b) it is shown the system response, compared to the time behavior of the measured strain. Here the performances of the model described so far are compared to those of a model proposed in [10], where the dependence of memory on the stress is not taken into account.

**Figure 15.** Comparison of the model's behavior to experiments for a prescribed history. Stress (load), magnetic field (curr) and *ζ* (corr-curr) vs. time (a); output displacement compared to the experimental data (b).

In conclusion of this section, it is worth to be further stressed that the procedure described so far is of general breath and can be exploited in modeling and compensation of any material/system with multi-variable hysteresis as most of multifunctional materials behave.
