Modeling and simulation of nonlinear magnetic hysteresis

Internal magnetic fields are present in all materials but with different intensity. In some materials the effect of the magnetic field is weak, but in ferromagnetic materials the intensity of the magnetic field is well pronounced. The magnetisation of a material or substance is characterised by magnetic induction B [T] and magnetic field strength H [A/m] as given by equation (1).

\[B= \mu_0 \cdot (H+M) \tag{1}\]

where μ₀ is the permeability of the free space and M [A/m] is the magnetisation.

\[\mu_0 = 4 \cdot \pi \cdot 10^{-7} \frac{Wb}{A \cdot m}\]

Hysteresis is seen in ferromagnetic materials as a tendency to stay magnetised, after the applied field force has been removed. The general behaviour of a magnetic core, described in terms of B and H, is described as in the image below.

The most important parameters in a magnetic hysteresis B(H) curve are: the saturation level (B_sat), the magnetic remanence (also called retentivity) and the coercivity (also called coercive force). These parameters have the following physical meaning:

• Saturation level: is the maximum field density sustained by the magnetic material;
• Retentivity: is the level of magnetization which remain in the material after the magnetic field intensity is set to zero;
• Coercivity: represent the value of the magnetic field intensity required to return the magnetization of the material to zero.

The mathematical model proposed in this paper will simulate a magnetic material which has the B(H) curve similar to the one in the image above. This model can be considered in the context of control applications, which requires an accurate characterisation of the relation between input currents and magnetic field.

Jiles-Atherton Model

The phenomenon of ferromagnetic hysteresis may be represented by the Jiles-Atherton (JA) model. The JA model is using the magnetic field H as independent variable, which means that the input into the model is H and the output is B. The equations which describe the JA model are a set of differential and algebraic equations (2):

\[\left\{\begin{matrix}
H_{eff} &= H+ \alpha \cdot M\\
M_{anh} &= M_{sat} \cdot \text{tanh} \left ( \frac{H_{eff}}{\alpha} \right )\\
\frac{dM_{irr}}{dt} &= \frac{dH}{dt} \cdot \frac{M_{anh}-M_{irr}}{k \cdot \delta – \alpha (M_{anh}-M_{irr})} \\
M_{rev} &= c \cdot (M_{anh}-M_{irr})\\
M &= M_{rev}+M_{irr}
\end{matrix}\right.\]

The magnetisation (M) is decomposed into its reversible (M_rev) and its irreversible component (M_irr), corresponding respectively, to reversible and irreversible magnetic domain deformation. The relationship between these two components and the anhysteretic magnetisation (M_an) is based on physical considerations of the magnetisation process. The effective field strength (H_eff) takes into consideration both the magnetic field intensity and the magnetisation level (M).

The hysteresis is implemented using:

\[\delta = \text{sgn} \left ( \frac{dH}{dt} \right )\]

which takes the value of 1 or -1 depending on whether it’s the positive or negative half cycle of the hysteresis curve.

The advantage of the JA model lies in the small number of required parameters for the definition of the mathematical model. The parameters symbol, name, description and effects on the B(H) curve is described in the table below.

Symbol	Physical property	Effects on the B(H) characteristic
Msat [A/m]	Saturation magnetization	Higher value leads to bigger saturation values for the magnetization curve
a [A/m]	Shape parameter of Man	Higher value leads to a decrease of the Man slope
a [-]	Quantifies domain interaction	Higher value leads to steeper slopes for anhysteretic and magnetization curves
c [-]	Reversibility coefficient	Higher value leads to a narrower hysteresis curve
k [A/m]	Average energy required to break pinning sites	Higher value leads to wider hystresis curve

The parameters of the JA model can be identified using several methods like least square techniques, genetic algorithms or advanced optimisation procedures.

Magnetic hysteresis Xcos block diagram model

The set of equation (2) are implemented into Scilab/Xcos in order to simulate the magnetic hysteresis phenomenon. A continuous-time model is built and the simulations are performed using different integrations method. All solvers should output the same results, in order to prove the numerical stability of the model.

The Xcos implementation of the equation (1) is the following:

The magnetic field strength H is fed as an input in the model. The JA model calculates the magnetisation M, which is used to calculated the magnetic induction B. Both magnetic filed strength and magnetic induction are saved in the Scilab workspace as the structure variables Hsim and Bsim.

The magnetic field strength H is generated through a Signal Builder with the following setup:

Spline method (0 … 7)	1
x	x_time
y	y_H
Periodic signal (y/n)?	n
Launch graphic window (y/n)?	n

The values of the magnetic field strength H and time are defined in a Scilab script (*.sce) as:

x_time = [  0,1,2,3,4,5,6,7,8,9,10,11,12,13,14,15,..
            16,17,18,19,20,21,22,23,24,25,26,27,28,29,30];
y_H = [     0.0E+00,0.0E+00,5.0E+02,-5.0E+02,1.0E+03,-1.0E+03,..
            2.0E+03,-2.0E+03,3.0E+03,-3.0E+03,4.0E+03,-4.0E+03,..
            5.0E+03,-5.0E+03,6.0E+03,-6.0E+03,7.0E+03,-7.0E+03,..
            8.0E+03,-8.0E+03,9.0E+03,-9.0E+03,1.0E+04,-1.0E+04,..
            1.1E+04,-1.1E+04,1.2E+04,-1.2E+04,1.3E+04,-1.3E+04,0.0E+00];

We can use a plot() function to visualise the variation of the magnetic field strength in time.

Even if we define the magnetic field strength as an input, we save it with a To workspace block, in order to have the same number of points as the magnetic induction variable. The period of the Clock block is set at 0.01 s.

The details of JA block diagram model can be seen in the image below. This model is implementing the set of equations (2) defined above.

Xcos simulation results

The simulated magnetic characteristic calculates the initialisation (first magnetisation curve), starting from an unmagnetised state, followed by its saturated hysteresis loop. For the first simulations the following parameters are used: M_sat = 76500 A/m, a = 7012 A/m, α = 0.01, c = 0.18 and k = 3942.

The parameters can be defined in the Scilab console or in the Set context menu within Xcos.

Ms=76500;
a=7012;
alpha=0.01;
c=0.18;
k=3942;
miu0=1.2566e-5;

Several simulation are performed for different values of the shape parameter (α) and the average energy required to break pinning sites (k), in order to highlight theirs influence on the B(H) curve.

After running the Xcos model for 30 s, we can use the following Scilab instructions to plot the magnetic hysteresis curve:

plot(Hsim.values,Bsim.values)
xlabel('H [A/m]','FontSize',3)
ylabel('B [T]','FontSize',3)
title('x-engineer.org','FontSize',3)

The first run, with the above parameters, gives the following magnetisation curve:

The second run is with k = 2000 and α = 0.01.

The third run is with k = 3942 and α = 0.1.

The benefit of this model is that it can be used for simulation purposes, integrated in more complex models of actuators or motors, due to its simplicity and fast calculation time. In order to have a complete and validated model of our actuator or motor, the JA parameters must be identified by one of the methods mentioned above. The simulation results show that the JA model is flexible and usable for different values of the parameters (M_sat, a, α, c and k).

The model needs further improvements due to the fact that here is an unexpected increase in total magnetisation (M) just after the input current changes direction (near the loop tips). The cause was the inappropriate calculation of the total magnetisation between the point at which the current changed direction and the point where the magnetisation curve crossed the an-hysteric curve. It was decided that the
magnetisation should be completely reversible during this portion of curve. This is justified since the domain walls, once deformed in a particular direction by the applied field, must first relax to equilibrium and then be deformed in the opposing direction before any further irreversible change occurs. In order to accommodate this change to the equations, it is necessary to include an additional conditional clause to obtain the correct calculation of the total magnetisation.

References

[1] Tom Stroehla, OlivierRadler, Ralf Volkert, Veit Zoppig – “Hysteresis compensation of electromagnets“, FACTA Universitatis – Mechanical Engineering series, 2006
[2] J.V. Leite, S.L. Avila, N.J. Batistela, W.P. Carpes Jr., N. Sadowski, P. Kuo-Peng, P.A. Bastos – “Real Coded Genetic Algorithm for Jiles-Atherton Model Parameters Identification”
[3] Sina Valadkhan, Kirsten Morris, Amir Khajepour – “A review and comparison of hysteresis models for magnetostrictive materials“, journal of Intelligent Material Systems and Structures, University of Waterloo (Ontario, Canada), 2006
[4] Mark C. Williams, Ronald S. Vogelsong, Kenneth S. Kundert – “Simulation and Modeling of Nonlinear Magnetics“, The Designers Guide Community, 2006.

For any questions, observations and queries regarding this article, use the comment form below.

Don’t forget to Like, Share and Subscribe!