Identification of Preisach hysteresis model parameters using genetic algorithms

1 Introduction

Magnetic materials are widely used in many engineering applications like shape memory alloy, piezoelectric, piezoceramic, magnetostrictive, and electromechanical actuators. Magnetic field H (A/m), magnetization M (A/m), magnetic induction B (T), the susceptibility χ (dimensionless) and permeability µ (H/m) are some of the primary magnetic parameters. The magnetic behavior and properties of a material can be learned by studying its hysteresis loop. A hysteresis loop shows the relationship between induced magnetic flux (B) and magnetic field (H).

In the past, several efforts have been made to develop different types of static hysteresis model. The classical Rayleigh model of scalar ferromagnetism (Lord, 1887) represents the H-M relation by a Prandtl-Ishlinksil model of play-type. This model is valid for small fields and it is fully determined by four experimental parameters (saturation magnetic field H_s, saturation flux density B_s, remanence flux B_r and frequency f). The Fröhlich model (Frolich, 1881) is more adapted for the low hysteresis loss materials and is fully determined by four experimental parameters (H_s, coercitivity H_c, B_s and B_r). Despite their simplicity, these models have limited applicability as they are specific to a particular variety of magnetic materials. The Langevin’s models are more accurate as they describe ferromagnetic material’s magnetization more accurately. The best known is the model of Jiles and Atherton (1986). It is a physical model of magnetic hysteresis and represents the magnetization M as a combination of a reversible component M_rev which results from the bending of the domain wall (Bloch walls) in the magnetization process and an irreversible one M_irr which corresponds to the domain wall displacement. The implementation of Jiles-Atherton uses five parameters (4 numericals and one experimental). The identification of these parameters is based on an iterative procedure (Jiles et al., 1992) which may introduce convergence problems. Moreover, it is very sensitive to initial values of parameters chosen as starting point for the optimization (Marion, 2008).

One of the most famous phenomenological hysteresis models is the Preisach model (Preisach, 1935). The popularity of this model is linked to the inventive technique of handling the contribution of infinite number of hysteresis operators with the Preisach triangle, the staircase line and the Everett function (Mayergoyz, 1991). It is much appropriate for numerical implementation. Currently, Preisach model includes similar scalar and vector hysteresis models, like the product model (Kádár, 1987), the moving model (Della Torre et al., 1994) or the Prandtl–Ishlinskii model (Hassani et al., 2014). Depending on the used distribution function, the classical Preisach model is completely determined by 5 parameters: two of them are experimentally deduced (H_s, M_s) and the remaining are numerically identified.

In order to generate the magnetization process for a given magnetic material, it is necessary to identify the Preisach parameters and evaluate its performances with regard to experimental data. There are many solving algorithms, which are mainly classified into two groups: Deterministic methods and stochastic ones. Deterministic methods are rarely used as they are based on the resolution of the gradient of the objective function. While, stochastic methods can be adapted to different forms of problems. They are based on a random evolution and require a lot of evaluation of the objective function to constantly end up guessing the optimum. Among these stochastic methods, the most widely used for magnetic domain are neural network (Zakerzadeh et al., 2011), genetic algorithm GA (Anh and Kha, 2008; Belkebir et al., 2009), particle swarm optimization (PSO) and nonlinear least squares method (Levenberg, 1944; Marquardt, 1963).

This paper describes an approach of Preisach parameter identification using GA. The proposed modeling method has achieved significant improvements in both accuracy and computation time, compared to the nonlinear least square identification method. The proposed approach can be applied to identify any type of hysteresis models. To demonstrate the efficiency of the proposed model, experimental results on hysteresis of a fully-process non oriented Fe-3 wt%Si steel sheet were provided and compared.

2 Classical Preisach model

The Preisach analysis (Preisach, 1935; Mayergoyz, 1991), mainly used for describing static magnetic behavior for ferroelectric and ferromagnetic materials, has been used to model iron-silicon films, piezoelectrical materials, super elastic response of Shape Memory Alloys for seismic dampers and arterial stents and also thermal effects on magnetic behavior (Quondam et al., 2016).

The magnetic state of the material is represented by magnetic entities (hysterons) having two possible saturation states M = 1 and M = −1, as shown in Fig. 1.

Close

Fig. 1

Preisach elemental operator.

Export to PPT

where, $\alpha$ and $\beta$ are the up and down switching values, respectively, and the arrows mark the authorized paths equivalent to changes in input.

The distribution of the elementary cycles defines the Preisach distribution function is expressed by (1).

The Preisach density is defined in the domain $S=\left\{(\alpha \text{,}\beta )/-H_S⩽\beta ⩽\alpha ⩽H_S\right\}$ as shown in Fig. 2.

Close

Fig. 2

Preisach triangle.

Export to PPT

Indeed, for each t time of magnetization, the Preisach plane {S} is divided into two subdomains S+ and S− separated by a line L(t). This boundary line, with its stairway shape, generates the magnetic hysteresis. Actually, the vertical portions result from increasing excitation field values and the horizontal portions result from decreasing values of this field. It is generally described by two particular points; the starting point (α = H_s, β = −H_s) and the ending point (α = β = H(t)) which is located on the line α = β.

Different methods for identifying the distribution function have been developed. They can be identified using numerical approaches (Mayergoyz, 1991; Rousseau et al., 1997; Everett, 1955; Ben Abou et al., 2003) or using analytical approaches (Preisach, 1935).

The Lorentz Modified Function LMF is set by the coercitivity H_c, a regulator coefficient k and two parameters a and b. The distribution function by LMF is then given by (2):

Then, the total magnetization M(t) is expressed by (3)

3 Genetic algorithms

Genetic algorithms are stochastic algorithms. They consist on evolutionary optimization based on genetic mechanism and evolution of species of Darwin in 1960 with concepts like ‘mutation’, ‘crossing over’, ‘survival of the best individual’, ‘natural selection’… (Holland, 1975; Goldberg, 1989). Generally, the GA attempts to converge to the optimal solution by creating new individuals from an initial population according to a random process and the population progresses by successive generations recombining the best solutions together. Particularly, the GA works with ‘parameters coding’. Once genes codified, the identification process begins and follows the schematic presented in Fig. 3; first, an initial set of individuals that represents a solution of the considered problem is generated. Each individual is evaluated on its performance with respect to the fitness function imposed. The individual follows a selection process where the best (fittest) survives and selected to carry on the regeneration and others disappear. Parents exchange their genetic information randomly to produce innovative child population by crossover and mutation selected operators. The parents are then replaced in the population by the children to keep the population size stable.

Close

Fig. 3

Flow diagram of the GA.

Export to PPT

This reproduction (Selection, crossover and mutation) is repeated and takes place with a probability of crossover (P_c) and probability of mutation (P_m), and it is again subjected to an evaluation of the relevance of the solution until the solution tends to the global optimum in a fixed iteration number.

In order to test the success of the model, several statistical criteria can be used. We can note ME (mean error), RMSE (root mean square error), R² (coefficient of determination) and MAPE (mean absolute percentage error) which are based on comparing the estimated variables with the original one. They are expressed as follows (Hosseini et al., 2016):

The major ending criteria is the satisfaction of the imposed fitness value (error between experimental data and pretended ones). Then, the best fit time evolution is the second condition. If the parameter remains considerably unchanged after five consecutive iterations, then it should be a stop criterion (Consolo et al., 2006). Finally, the third criterion characterizes the number of population or generations initially fixed.

The benefit of GA is that they are suitable for complex problems. They are applied to functions with great number of parameters. In addition, GA is one of the most practical methods as it does not depend on the choice of the initial population. However, it requires many evaluations of the objective function to achieve the global minimum as it does not converge quickly. Several researchers have used the GA in their optimization studies and identification of parameters of magnetic models (Belkebir et al., 2009; Quondam et al., 2016; Wilson et al., 2001; Salvini and Fulginei, 2002; Che-Hang and Guangjun, 2007). All these studies have shown a considerable reliability of this algorithm.

4 Magnetic measurements test apparatus

Various magnetic measurement systems have been developed for different applications, such as magnetic behavior under mechanical stress and/or thermal modeling of magnetic components. The used magnetic experimental measurements test apparatus (Fig. 4) was validated by several experimental studies (Matsubara et al., 1989; Iordache et al., 2003). It consists of two U ferrite cores (Philips U100/57/25-3C90, MnZn material) maintained in contact with the sample close the magnetic circuit. Primary windings are wound on the central limbs of the yokes and secondary winding surround the specimen. This double-yoke arrangement leads to better homogeneous distribution of the magnetic field in the measurement zone, moreover it minimizes the negative effects of the overhang and of the eddy currents on the measurements accuracy.

Close

Fig. 4

Block diagram of the magnetic measurement test apparatus.

Export to PPT

The magnetizing current waveform is controlled by a function generator and an operational amplifier. The voltage V(t) represents the direct output of the secondary coil and i(t) the input magnetizing current.

The magnetic field is estimated according to the IEC method (IEC, 1992) by:

The magnetic flux Φ_S detected by the secondary winding is the sum of the flux experienced by the ferromagnetic specimen Φ_Sample and the magnetic flux in the air between the secondary winding and the specimen Φ_air:

The Faraday’s law requires that the secondary voltage V(t) on the secondary coil is proportional to the variation of the magnetic flux as:

The magnetic induction B(t) is calculated by integrating the secondary voltage V(t) with respect to Faraday’s law (13):

The used sample is a fully-process non oriented Fe-3 wt%Si steel sheet of 0.35 mm thick. The specimen is strips 20 mm wide and 250 mm long cut out in the rolling direction and vacuum annealed at 720 °C for 2 h in order to eliminate the residual stresses which originate from the manufacturing process (Hubert, 1998). Fig. 5 Shows the experimental obtained hysteresis.

Close

Fig. 5

Experimental hysteresis of a Fe-3 wt%Si steel.

Export to PPT

5 Results and discussion

In this section, we investigate the ability of GA to identify the static Preisach model’s parameters that ensure reproducing the measured hysteresis curve of the Fe-3 wt%Si material. A comparative study is developed between GA and non-linear least square method (LSQ) in order to spot the more accurate. The identification procedures are established according to the following diagrams (Figs. 6 and 7) and are based on the minimization of two kinds of errors: the mean squared error and the percentage error described by Eqs. (14) and (15) respectively.

Close

Fig. 6

Flow diagram of the used GA method (using ε₁).

Export to PPT

Close

Fig. 7

Flow diagram of the LSQ method (using ε₁).

Export to PPT

Concerning the GA approach, the GA Matlab toolbox is used. The evaluation is carried out for an initial population of 20 individuals and 200 iterations. The crossover probability P_c and the mutation probability P_m are taken respectively equal to 0.8 and 0.2. Minimum and maximum expected values for the five Preisach parameters {a, b, k, H_s, H_c} are respectively {0, 0.3, 1000, 3500, 380} and {0.5, 1, 1500, 4500, 500}. Once the initial population is created, the Preisach model calculates the magnetic induction, which is then evaluated using the objective function (ε₁ or ε₂). The fitting results between measured and simulated values are provided at each run and compared to the major stop criterion ε_min (=0.05). If convergence is not succeeded, the algorithm creates a new population by applying genetic operators (crossover and mutation). This new population is also evaluated and the process is repeated until convergence is reached.

Concerning the LSQ method, it consists on an iterative improvement of parameter values in order to reduce the sum of the squares of the errors between the function to be calculated and the measured data points. For this purpose, it requires an initial guess for the parameters to be estimated. So, we chose randomly the values a = 0.1, b = 0.75, k = 1200, H_s = 4400 and H_c = 420 for the initial guess. If the LSQ error is ‘acceptable’, the program ends and returns the optimal parameters. Otherwise, the algorithm resumes calculation for another initial parameters vector until 200 iterations.

According to the symmetry of the hysteresis, we just use the experimental half upstroke of the major loop for the identification process. The results of the optimization for the two methods for the conditions already mentioned are summarized in the following Table 1. The obtained curves using these identified parameters are shown in Fig. 8.

Close

Fig. 8

Comparison of simulated cycles by GA and LSQ methods and experimental major cycle.

Export to PPT

GA error ε₁ evolution is shown in Fig. 9. The random choice of the initial population leads to a high value of the error. The creation of a new population by selecting individuals that have the better fitness value from the previous population improves drastically the error. For GA implementation, even if the fitness value ε_min is not reached, and in case of the best individual of the population remains unchanged for a given number of generations (30), the algorithm converges. The obtained final solution is the optimal one. Fig. 9. Shows that GA converges at the 78th generation with a final error ε₁ equal to 0.024 T.

Close

Fig. 9

Error evolution.

Export to PPT

GA is based on a random choice of the initial population. Thus, we have different optimal points for the same upper and lower boundaries of designed parameters. Table 2 shows the results of five consecutives identifications. Identified half upstroke of the major loop are plotted on Fig. 10.

Close

Fig. 10

Comparison of optimization results by GA and the experimental branch.

Export to PPT

According to the sensitivity tests, we note that all the tests give errors less than 0.03 for ε₁ and less than 1% for ε₂. The parameter values are closely identical and gives superimposed curves. We note also that a strong dependency of the five model parameters a, b, k, H_m and H_c is existing.

6 Conclusion

In this paper we have proven how the static Preisach model’s parameters have been efficiently determined from limited set of experimental data (major loop) by the genetic algorithm approach and the non-linear least square method based on the Levenberg Marquardt algorithm. The results reached from the GA method are in a good agreement with the measures within an error of one percent and it produces the nearest cycle to the experimental reference one. So, it is more efficacious than the LSQ approach. Nevertheless, to prevent premature convergence and guarantee solution accuracy, the GA requires a large population size. Thus, a large population size usually demands more generations and consequently more computation time, which is a major weakness of this method.