Hysteresis behaviour of different magnetorheological elastomer models: modelling and simulation
Thaer M. I. Syam^{1} , Asan G. A. Muthalif^{2}
^{1, 2}Departament of Mechanical and Industrial Engineering, College of Engineering, Qatar University, Doha, Qatar
^{1}Corresponding author
Vibroengineering PROCEDIA, Vol. 31, 2020, p. 714.
https://doi.org/10.21595/vp.2020.21408
Received 4 April 2020; accepted 14 April 2020; published 7 May 2020
JVE Conferences
Smart materials are materials by which their properties can be altered due to external stimuli such as temperature, pressure, and magnetic field. Magnetorheological Elastomer (MRE) is a type of smart composite material consisting of a polymer matrix embedded with ferromagnetic particles. In the presence of an external magnetic field, its mechanical properties such as stiffness change due to the interaction between the magnetic particles. Vibration isolation has been of interest to researchers for decades. In this work, simulation studies were done by utilizing different MRE models from the literature such Bignham, BoucWen, Modified BoucWen, Dahl and Hysteresis models on a base motion isolation (base excitation) system in order to understand the characteristics of MRE. Results showed that as the magnetic field increases, the stiffness of MRE increases significantly. Shifting in the natural frequency of the system from the transmissibility curve was observed for all models which lead to vibration isolation.
 The main goal of this analytical study was to test the performance of different MRE models on a system subject to base motion (excitation) isolation.
 Results showed that as the magnetic field increases, the stiffness of MRE increases significantly.
 It can be concluded that the MRE stiffness properties vary more significantly than damping.
 Shifting in the natural frequency of the system from the transmissible curve was observed for all models which lead to vibration isolation.
 The results given here for the base excitation system will be useful to predict the performance of MRE in vibration isolation applications such engine mounting system.
Keywords: mathematical models, magnetorheological elastomer (MRE), stiffness and damping, transmissibility curve, and Simulink.
1. Introduction
The main goal in vibration control is to eliminate and reduce the vibrations in mechanical systems which also mean that vibration should be isolated from the systems. There are two types of interest namely; force isolation and motion isolation which can be identified from a system with either rotational or linear motion [1]. Base motion isolation is being studied analytically using different models of magnetorheological elastomer (MRE). MRE is a smart composite material that consists of a nonmagnetizable polymer matrix embedded with micro/nanosized ferromagnetic particles. These materials can alter their mechanical properties such as elasticity and stiffness when exposed to an external magnetic field [1, 2]. This property makes them attractive to many engineering applications such as dampers and vibration isolators [3, 4]. The internal magnetic particles which form like a chain in the direction of the magnetic field are the main factor in determining the mechanical properties of these elastomers [4]. Many phenomenological models have been developed to deal with the viscoelastic behaviour of MR dampers described by Hooke’s and Newton’s laws [5]. They are called viscouselastic materials. Table 1 shows the schematic representations and mathematical models of these models. One of the most used model to describe the viscouselastic behaviour of MR damper is the Bingham model. It assumes the material to be rigid before yielding [6]. [7] developed a model that shows the hysteretic behaviour of MR damper. Similar to BoucWen model, [7] adopted a 14parameter model; the modified BoucWen model. This model considers the high nonlinear behaviour of MR damper. Dahl model is one of the models developed to describe the behaviour of MR dampers [8]. [9] proposed a new model that utilizes a hyperbolic sine function to represent the hysteresis in order to describe the stiffness and damping characteristics. [10] used MRE as a semiactive vibration isolator to suppress torsional mode of vibration. Results showed that a shift of 3.9 Hz of the natural frequency and slight changes in the damping properties were observed. [11] developed the mathematical and experimental setup for dynamic tuning of torsional vibration; it is found that a shift in the curve was observed in the transmissibility curve (${\theta}_{2}/{\theta}_{1}$) as the magnetic field increases.
Table 1. Schematic and mathematical representations of different MRE models [710]
Schematic representation of models

Mathematical expression by the damping force,$F$

Bingham model

$F={f}_{c}sgn\left(\dot{x}\right)+{c}_{0}+{f}_{0}$
where ${f}_{c}$ is the friction force due to yield stress is, ${c}_{0}$ is the damping coefficient, ${f}_{0}$ is the stored force due to the accumulator

BoucWen model

$F=\alpha z+{c}_{0}\dot{x}+{k}_{0}(x{x}_{0})$
$\dot{z}=\gamma \left\dot{x}\rightz\left({\leftz\right}^{n1}\right)\beta \dot{x}{\leftx\right}^{n}+A\dot{x}$
where $\alpha $ is the scaling BoucWen model parameter related to yield stress of MR damper, ${x}_{0}$ is the initial displacement due to effect of the accumulator. ${c}_{0}$ and ${k}_{0}$ are the damping coefficient and spring stiffness, respectively. $z$ is the deformation of the model due to the hysteretic effect, $\gamma $, $\beta $, $A$ and $n$ are Boucwen model parameters

Modified BoucWen model

$F={c}_{1}\dot{y}+{k}_{1}\left(x{x}_{0}\right)$
$\dot{y}=\frac{1}{{c}_{0}+{c}_{1}}[\alpha z+{c}_{0}\dot{x}+{k}_{0}\left(xy\right)]$
$\dot{z}=\gamma \left\dot{x}\dot{y}\right{\leftz\right}^{n1}z\beta \left(\dot{x}\dot{y}\right){\leftz\right}^{n}+A(\dot{x}\dot{y})$
where ${k}_{1}$ is the stiffness of the accumulator, ${x}_{0}$ is the initial displacement of spring ${k}_{1}$ due to the accumulator effect, $y$ is the internal displacement of the damper, ${c}_{1}$ and ${c}_{0}$ are the viscous damping coefficients at low and high velocities, respectively. $z$ is an evolutionary variable. $\gamma $, $\beta $, $n$ and $A$ are modified Boucwen model parameters

Dahl model

$F={k}_{0}x+{c}_{0}\dot{x}+\delta z{f}_{0}$
$\dot{z}=\rho (\dot{x}\left\dot{x}\rightz)$
where ${k}_{0}$ is the stiffness of the spring, ${c}_{0}$ is the damping coefficient, $z$ is the intermediate variable and $\rho $ is the coefficient of stiffness

Hysteresis model

$F={c}_{0}x+{k}_{0}\dot{x}+\alpha z+{F}_{0}$
$z=sinh\left(\beta x\right)$
where ${F}_{0}$ is the offset of the isolator force, $\beta $ is the scale factor of the displacement of the isolator, $\alpha $ is a factor to scale the hysteresis, and ${c}_{0}$ and ${k}_{0}$ are the damping coefficient and spring stiffness, respectively

2. Analytical simulations of different MRE models using Matlab
2.1. Mathematical model development and simulation parameters
The purpose of this study is to use MRE in semiactive vibration isolation by which the parameters can be controlled. The system is presented as a simple base motion isolation (base excitation) system by which the excitation input to the base is to be isolated on the mass. This model is developed to represent the MRE mathematically. The mathematical model of the system is a description in terms of equations. Newton’s and conservation of energy are the physical laws as well as the basics in building the model. This system considers one degree of freedom base excitation system as shown in Fig. 1. The system is excited by an input displacement $y\left(t\right)$ applied to the base and an output displacement $x\left(t\right)$ on the mass. Five models are selected from the literature to be modelled and simulated on the base excitation system. MRE properties vary under a magnetic field. Each model has its mathematical representation and equations by which their parameters are a function of a magnetic field; this magnetic field is inputted as an electric current. The simulation is done for each model with and without MRE. When the system is without MRE, it means that the effect of stiffness and damping is not related to the current. For this simulation, the current is varied from 0 A to 4 A.
Fig. 1. Schematic representation of base isolation model utilizing MRE
The equation of motion is based on Newton law and it becomes:
Knowing that the force of the MRE is given by:
Then, Eq. (4) becomes:
This is the equation of motion presents the dynamics of the system modelled in terms of $m$, the mass, $y$, the input displacement, $x$, the relative output mass displacement, ${c}_{p}$, the passive damping coefficient of MRE when no magnetic field is applied, ${c}_{mre}$, the active damping coefficient of MRE when the magnetic field is applied, ${k}_{p}$, the passive stiffness of MRE when no magnetic field is applied, ${k}_{mre}$, the active damping coefficient of MRE when the magnetic field is applied. The passive term of the damping coefficient ${c}_{p}$ isn’t presented in Fig. 1 and this is due to the fact that its change under magnetic field is not considerable. However, it’s included in the total MRE damping force ${F}_{mre}$ as in Eq. (4) To simplify the equation, let $k={k}_{p}+{k}_{mre}$ and $c={c}_{p}+{c}_{mre}$, Eq. (5) becomes:
The system transfer function can be achieved by taking a Laplace transform to the differential equation by considering the Laplace transform ($a=A\left(s\right)$) on Eq. (6):
By using the conversion from ($s$domain) to ($\omega $ domain) through substituting ($s=j\omega $), the transmissibility factor is given by:
By taking the magnitude of Eq. (10) the amplitude of displacement transmissibility is:
The mass in this analysis is assumed to be 5 kg. The inputs signals are step input (step time: 0.5, initial value: 0, and final value: 1) and a sinusoidal wave (Amplitude: 1 m, frequency: 10 rad/s). Table 2 below briefly summarizes the five models that will be studied. Also, the mechanical properties of these models of MRE are changeable by applying a current. These relationships with their parameters are also presented which were taken from the literature [8, 1214]. Some of models parameters are constant and others are fielddependent; by which they depend on the magnetic field which related to the current as shown in Table 2.
Table 2. Fielddependency between models parameters and current [8, 1214]
MRE model

Fielddependent parameters equations

Bingham

${c}_{0}=147.5i+13.52$[N.s/m], ${f}_{c}=394.8i+38.14$[N], ${f}_{0}=663.56i+52.19$[N]

BoucWen

${k}_{0}={{k}_{0}}_{a}+{{k}_{0}}_{b}i$[N/m], ${c}_{0}={{c}_{0}}_{a}+{{c}_{0}}_{b}i$[N.s/m], $\alpha ={\alpha}_{a}+{\alpha}_{b}i$[N/m], $A$, $\beta $, $\gamma $, and$n$ are constants, where ${{k}_{0}}_{a}=$ 469 N/m, ${{k}_{0}}_{b}=$ 323.2 N/m.A, ${{c}_{0}}_{a}=$ 2100 N.s/m, ${{c}_{0}}_{b}=$ 350 N.s/m.A, ${\alpha}_{a}=$ 14000 N/m, ${\alpha}_{b}=$69500 N/m.Am, $A=301$, $\beta =$3630000 m^{2}, $\gamma =$190 s^{1} and $n=$ 2

Modified BoucWen

${k}_{0}={{k}_{0}}_{a}+{{k}_{0}}_{b}i$ [N/m], ${k}_{1}={{k}_{1}}_{a}+{{k}_{1}}_{b}i$[N/m], ${c}_{0}={{c}_{0}}_{a}+{{c}_{0}}_{b}i$ [N.s/m], ${c}_{1}={{c}_{1}}_{a}+{{c}_{1}}_{b}i$ [N.s/m], $\alpha ={\alpha}_{a}+{\alpha}_{b}i$ [N/m], $A$, $\beta $, $\gamma $, and $n$ are constants, where ${{k}_{0}}_{a}=$ 500 N/m$$, ${{k}_{0}}_{b}=$ 15.9 N/m.A, ${{k}_{1}}_{a}=$ 469 N/m , ${{k}_{1}}_{b}=$ 323.2 N/m.A, ${{c}_{0}}_{a}=$ 2100 N.s/m, ${{c}_{0}}_{b}=$ 350 N.s/m.A, ${{c}_{1}}_{a}=$ 28300 N.s/m, ${{c}_{1}}_{b}=$ 295 N.s/m.A, ${\alpha}_{a}=$ 14000 N/m, ${\alpha}_{b}=$ 69500 N/m.A, $\beta =$3630000 m^{2}, $A=301$, $\gamma =$190 s^{1} and $n=$ 2

Dahl

${k}_{0}={{k}_{0}}_{a}+{{k}_{0}}_{b}i$[N/m], ${c}_{0}={{c}_{0}}_{a}+{{c}_{0}}_{b}i$ [N.s/m], $\delta ={\delta}_{a}+{\delta}_{b}i$, $\rho $ is constant, where $\rho =$25 kg/m^{3}, ${f}_{0}=0.001$, ${{k}_{0}}_{a}=$ 800 N/m, ${{k}_{0}}_{b}=$ 250 N/m, ${{c}_{0}}_{a}=$ 24.64 N.s/m, ${{c}_{0}}_{b}=$ 40.41 N.s/m.A, ${\delta}_{a}=0.031$, and ${\delta}_{b}=1.12$

Hysteresis

${c}_{0}={{c}_{0}}_{a}+{{c}_{0}}_{b}i$ [N.s/m], ${k}_{0}={{k}_{0}}_{a}+{{k}_{0}}_{b}i$[N/m], $\alpha ={\alpha}_{a}+{\alpha}_{b}i$, $\beta $ is constant, where, ${{c}_{0}}_{a}=$ 315.9 N.s/m, ${{c}_{0}}_{b}=$ 2.74 N.s/m.A, ${{k}_{0}}_{a}=$ 246.7 N/m, ${{k}_{0}}_{b}=$ 31.2 N/m.A, ${\alpha}_{a}=$ 17653 N/m, ${F}_{0}=$ 0 N, ${\alpha}_{b}=$ 1.553 N/m.A, and $\beta =$ 8570000 m^{2}

The combination of the base isolation system and the damping force by these models is done by combining the equations from Table 1 and Eq. (2). For example, the mathematical expression for this methodology can be given as follows, reconsidering the Bingham model mathematical expression from Table 1, $F={f}_{c}sgn\left(\dot{x}\right)+{c}_{0}+{f}_{0}$, From Eq. (2) for base motion isolation system is:
By combining Eq. (1) and Eq. (12), and considering the relative motion between the mass output motion $x\left(t\right)$, and base input excitation $y\left(t\right)$, by which $x=\left(xy\right)$ and $\dot{x}=\left(\dot{x}\dot{y}\right)$, Bingham model damping force becomes:
The combination of Bingham model damping force and the base isolation system becomes:
The equation in this form is preferable since it’s a secondorder differential equation to build the block diagram on Simulink. All models followed the same procedure to be built on Simulink. An example of the combination of equations of these five models can be shown for Bingham model as in Eqs. (2224). In parallel, the block diagrams are related also to some MATLAB codes developed in order to be simulated. Firstly, the MATLAB code is run for each model which is consists of the fielddependent parameters equations as shown previously in Table 2. Then, the models will be automatically opened and simulated using ‘open’ and ‘sim’ commands in MATLAB. The simulations are done twice for a step input and sinusoidal input with the parameters shown previously. A manual switch is used to connect each input once at a time. Time domains for both inputs are taken as the output of these simulations. Finally, the system with and without MRE are shown in scope and extracted. Stiffness and damping coefficient of each model in the system was calculated from the 2nd order underdamped step responses. All of the five MRE models are summed up as subsystems so that simulations are done in at the same time.
3. Analytical results and discussion
Stiffness and damping coefficients values are plotted with respect to the applied current, as shown in Figs. 2, 3. It can be clearly seen from Fig. 2 that the stiffness of the MRE for all models increases significantly when the current is increased. For example, the stiffness increases from 150 N/m to 200 N/m at an applied current 0 A to 1 A, respectively. This shows a 33.3 % increase in the stiffness when the current is applied initially. In contrast, the damping coefficient increases from 13.5 N.s/m to 17.945 N.s/m and an applied current 0 A to 1 A, respectively as shown in Fig. 3. This shows a 26.82 % increase in the damping coefficient. For BW model, the stiffness increases by 78 % when 1 A current is applied, an increase of 42.6 % is obtained.
Fig. 2. Stiffness values for each MRE model with respect to the applied current
Fig. 3. Damping coefficient values for each MRE model with respect to the applied current
The stiffness shows a rise of 63.5 %, 78 % and 90 % for MBW, Dahl and Hysteresis models, respectively whereas the damping coefficient increased by a percent of 19.5 %, 42.3 % and 43.3 %, respectively. Also, the range of the stiffness values is more than the range of the damping coefficient values. For example, the minimum and maximum values of stiffness for MBW model are 212.32 N/m and 1322.46 N/m whereas damping coefficient values are 14.32 N.m/s and 21.32 N.m/s, respectively. This can be seen for all other models as well. Therefore, an increase in the stiffness and damping coefficient is noted when the current is increased. However, the change in the stiffness of the MRE compared to the damping coefficient is significant. This proves the theory mentioned in the literature by which the stiffness of the MRE varies significantly with the magnetic field, whereas the damping characteristics do not. Although the change in the damping in these simulations still considerable but it can be concluded that the MRE stiffness properties vary more significantly than damping. One reason for the significant change is that the equations relating the model parameters might have some magnetorheological fluids (MRF) characteristics. The damping behaviour of MRF is substantial with a magnetic field unlike MRE.
Fig. 4. Displacement transmissibility vs frequency for: a) Bingham model, b) BoucWem model, c) Modified BoucWen model, d) Dahl model, and e) Hysteresis model
After obtaining $c$ and $k$ values for each MRE model, the transmissibility factor is plotted with respect to the frequency for each model at different current. The range of the frequency for this simulation was set to be 0 H to 100 Hz. It can be clearly seen from Fig. 4 that the natural frequency of the system is increased and the transmissibility curve is shifted to the right. This clearly shows that the curve is being shifted towards the reduction region which means that the vibration is isolated. This shift in the transmissibility curve is evidence showing that the stiffness of the MRE is increased when the applied current increases. In addition, the amplitude is being reduced which means that the damping ratio is increased so that the system is damped better so the vibration is isolated faster. However, as concluded earlier that the change in the damping properties is slighter than the change in the stiffness of the MRE. The transmissibility curve for Bingham model is shifted in a small increment as shown in Fig. 4(a). The variation in the stiffness values is minimum, this due to the simplicity of this model as it has only three parameters related to the current. MBW (Fig. 4(c)) and Dahl (Fig. 4(d)) models show more shifting the transmissibility curve as the applied current increases. This due to the more accurate approximation of the fielddependent parameters that control the hysteresis loops. Altogether, the modelling approach presented in this work is adequate to predict the hysteresis behaviour of different MRE models obtained from the literature. These models were applied on a system subjected to base motion excitation and the transmissibility factor was studied. This study has vast application in vibration isolation. The results given here for the base excitation system will be useful to predict the performance of MRE in vibration isolation applications such engine mounting system.
4. Conclusions
The mathematical model was derived for a base motion isolation system and combined with the mathematical models of five different MREs for step and sinusoidal excitation inputs. It can be concluded from the analytical simulation study that different MRE models were able to reduce vibration in the base isolation system. Bingham, BoucWen, Modified BoucWen, Dahl, and Hysteresis loops were modelled using Simulink and MATLAB. Fielddependent parameters were identified by mathematical expressions. The stiffness and damping coefficients were approximated from the secondorder underdamped systems obtained from the step response. It can be clearly deduced that the change in the stiffness for these models is way more than the change in the damping properties. The change in the stiffness was found to be approximately 2 times higher than the change in the damping coefficient which confirms the literature survey. The transmissibility factor was calculated and plotted against a range of frequencies; results showed that there is a shift in the natural frequency of the system for all five models which is due to the change in the stiffness as the applied current is increased. This paper contributed to reviewing different MRE models and studied their performance on the base motion isolation system. The main goal of this study was to test the performance of different MRE models in terms of damping and stiffness properties under an external magnetic field. The transmissibility factor vs excitation frequency was plotted to visualize the efficiency of fielddependent properties in vibration isolation.
Acknowledgements
This work was partly supported by Qatar University Student Grants (QUST1CENG202018).
References
 Wen H., Guo J., Li Y., Liu Y., Zhang K. The transmissibility of a vibration isolation system with ballscrew inerter based on complex mass. Journal of Low Frequency Noise Vibration and Active Control, Vol. 37, Issue 4, 2018, p. 10971108. [Publisher]
 Jolly M. R., Carlson J. D., Muñoz B. C., Bullions T. A. The magnetoviscoelastic response of elastomer composites consisting of ferrous particles embedded in a polymer matrix. Journal of Intelligent Material Systems and Structures, Vol. 7, Issue 6, 1996, p. 613622. [Publisher]
 Bastola A. K., Li L. A new type of vibration isolator based on magnetorheological elastomer. Materials and Design, Vol. 157, 2018, p. 431436. [Publisher]
 Kumbhar S. B., Chavan S. P., Gawade S. S. Adaptive tuned vibration absorber based on magnetorheological elastomershape memory alloy composite. Mechanical Systems and Signal Processing, Vol. 100, 2018, p. 208223. [Publisher]
 Xu Z., Wang Q., Zhu K., Jiang S., Wu H., Yi L. Preparation and characterization of magnetorheological elastic polishing composites. Journal of Intelligent Material Systems and Structures, Vol. 30, Issue 10, 2019, p. 14811492. [Publisher]
 Kamath G. M., Wereley N. M. A nonlinear viscoelasticplastic model for electrorheological fluids. Smart Materials and Structures, Vol. 6, Issue 3, 1997, p. 351359. [Publisher]
 Spencer B. F., Dyke S. J., Sain M. K., Carlson J. D. Phenomenological model for magnetorheological dampers. Journal of Engineering Mechanics., Vol. 123, Issue 3, 1997, p. 230238. [Publisher]
 Dahl P. R. Solid friction damping of mechanical vibrations. AIAA Journal, Vol. 14, Issue 12, 1976, p. 16751682. [Publisher]
 Yu Y., Li Y., Li J. A new hysteretic model for magnetorheological elastomer base isolator and parameter identification based on modified artificial fish swarm algorithm. 31st International Symposium on Automation and Robotics in Construction and Mining, 2014, p. 176183. [CrossRef]
 Nor Alias F., Asan Muthalif G. A., Khairul Arpan A. M., Diyana Nordin N. H. Experimental investigation of static properties of magnetorheological elastomer. Iranian Journal of Science and Technology Transactions of Mechanical Engineering, Vol. 42, Issue 2, 2018, p. 185197. [Publisher]
 Hashi H. A., Muthalif A. G. A., Diyana Nordin N. H. Dynamic tuning of torsional transmissibility using magnetorheological elastomer: Modelling and experimental verification. Iranian Journal of Science and Technology Transactions of Mechanical Engineering, Vol. 40, Issue 3, 2016, p. 181187. [Publisher]
 Roussel N., Lanos C., Toutou Z. Identification of Bingham fluid flow parameters using a simple squeeze test. Journal of NonNewtonian Fluid Mechanics, Vol. 135, Issue 1, 2006, p. 17. [Publisher]
 Zhu X., Lu X. Parametric identification of BoucWen model and its application in mild steel damper modelling. Procedia Engineering, Vol. 14, 2011, p. 318324. [Publisher]
 Dasari P. K. R., et al. Adaptation of the modified BoucWen model to compensate for hysteresis in respiratory motion for the listmode binning of cardiac SPECT and PET acquisitions: Testing using MRI. Medical Physics, Vol. 41, Issue 11, 2014, p. 112508. [Publisher]