A Review on Parametric Dynamic Models of Magnetorheological Dampers and Their Characterization Methods

Magnetorheological (MR) fluids are capable of manifesting a rheological behaviour change by means of a magnetic field application and can be employed in many complex systems in many technical fields. One successful example is their use in the development of dampers: magnetorheological dampers (MRDs) are widespread in vibration control systems, as well as civil engineering applications (i.e., earthquake or seismic protection), impact absorption and vibration isolation technology in industrial engineering, and advanced prosthetics in biomedical fields. In the past, many studies have been conducted on MRDs modeling and characterization, but they have usually been focused more on the theoretical models than on the experimental issues. In this work, an overview of both of them is proposed. In particular, after an introduction to the physics of the magnetorheological effect, a short review of the main mathematical models of MRDs is proposed. Finally, in the second part of this study an overview of the main issues that occur in MRDs experimental characterization is reported and discussed.


Introduction
Recently, the interest in the so-called "smart materials" has been considerably increased.One of the areas in which their possible applications are more fruitful, regards the damping of the mechanical vibrations.In this field, different materials have been proposed: piezoelectric materials ( [1][2][3][4][5][6]), shape-memory alloys ( [7][8][9]), electro-rheological fluids ( [10]), magnetorheological (MR) fluids, etc.The magnetorheological fluids are composed of a suspension of magnetically responsive particulate matters in a liquid phase.The particles arrange themselves to create very strong chains if a magnetic field is applied [11,12] since they are aligned and forced to stay along their respective flux lines.In this way they prevent the carrier fluid flow forming a barrier.The MR fluids can switch from a semi-solid to a liquid state in a few milliseconds.This change is attributable to a significant increase in its yielding shear stress τ [13][14][15] which is the lowest stress that the external forces have to overcome to begin the flow.
As mentioned above, an interesting and widespread application of MR fluids is the MR damper (MRD) that can be considered as a hydraulic cylinder containing the MR fluid.MRDs can be applied in active and semi-active vibration control devices since they are able to change the damping characteristic, fit it to real loads and show wide force ranges, low sensitivity to impurities or contaminants, good integration in automatic control systems, rapid response, mechanical simplicity, noiseless work, robustness and low power requirements [13,16] (In fact, they may reach yield stresses of 50 kPa or more under a magnetic field, which may result from a low power electromagnet (e.g., 2-24 V and 1-2 A) [13,17]).
MR dampers are widely used for vibration isolation devices, seismic protection, impact absorption and cable-stayed bridge systems (Table 1).Special MR dampers can also be used to mitigate vibrations in drill string as reported by [30] or in leg rehabilitation [31].Since the MRD is an intrinsically nonlinear device, even today a challenging activity is represented by the design and modeling of the proper control algorithms [32].In fact, during the past years, MRD systems design made substantial advances [13].Moreover, the availability of joints with variable or controllable damping could be used to adjust the behaviour of more complex systems [33] in E(3) [34] or in SE(3) [35] .
However, the use of MR fluids in dampers may be affected by some drawbacks (Table 2):  [38,39] If dampers remain unused for long time periods but an instantaneous action is required, settling may represent a critical issue, e.g., in civil engineering.This issue can be overcome by mixing the fluid after it has settled: some tests revealed that a MRD may recover the nominal value of the force after only one cycle of use, even if it has not been used for one year but a high quality fluid has to be employed [36,40].
The damper system may be subjected to hydraulic imbalance, limiting the range of damping forces that can be achieved.The phenomenon of imbalance takes place when the MRD is in the compression phase and the pressure drop across the base valve is lower than the pressure drop across the piston.Therefore pressure lags in the above piston chamber due to the fluid stream through the base valve may occur; to avoid that phenomenon damping forces should be tuned [40].Moreover, the time response of a MR damper could be limited by eddy currents in the coil's core, so an optimal design of the controller is often needed [41,42].
Since many of the above issues may be overcome at the design level of the devices, a more detailed and suitable modelling should be approached considering also the error sources that affect their experimental characterization and validation.However, very little attention has been paid to uncertainty evaluation on damper modelling and testing in the scientific literature [43][44][45][46].In this work, the authors reviewed the parametric dynamic models commonly used for predicting MR dampers behaviour and their errors.A particular contribution is given by describing the experimental setups and methods used for model validation and device assessment, as well as the common source of errors in device testing.

MR Fluid Characteristics
Usually MR fluids are fluids with two phase built by scattering in a non-magnetic liquid huge quantities of micron sized, highly magnetizable, solid particles (up to % 50 vol ).To this aim iron particles produced by the thermal decomposition of iron pentacarbonyl are utilized because of their high magnetization [47].The MR rheology depends on the particle shape distribution, particle concentration, additional additives, properties of the carrier fluid, temperature and the applied magnetic field [48].Polyesters, mineral and silicone oils, synthetic hydrocarbons, polyethers, and water are common carrier liquids [40,[49][50][51][52][53][54][55][56].In order to reduce aggregation and sedimentation, additives may be used.They also add additional lubricating properties [32,40]: they encompass surfactants, polymers, and thixotropic agents.The MR effect is usually associated to the magnetization field induced to the suspended particles.If a magnetic field is introduced, the suspended particles magnetize and align themselves across the field lines.As a result the fluid shows typically high values of yield stress, viscoelasticity and shear rate.The mechanical behaviour of a MR fluid relies on the diameter, composition and particle volume fraction [47,57].By changing the magnetic field magnitude and/or the magnetization saturation of the particles is possible to vary their mechanical properties.Table 3 summarizes the main characteristics of commercial MR fluids.The rheology of a MR fluid can be described in terms of two conditions (pre-yield and post-yield): where τ is the yield stress, γ e the deformation and G * the complex (or dynamic) shear modulus.
The MR fluids post-yield behaviour has been mathematically described and experimentally observed on many works [40,48], anyway the most common model is the particle magnetization model [40,47] that considers the action of the particles as magnetic multi-domains [40,47,59].In general, the solid particles dispersed in the fluid have a randomly aligned dipole moment if no magnetic field is applied.When a magnetic field is introduced they become ordered and aligned in one direction [48,60].The authors of [61] and [62] showed that, depending on the field intensity, other expressions should be adopted for more accurate results.

Parametric Dynamic Models of MR Dampers
To date, the dynamics of magnetorheological dampers have been outlined by developing many models.They differ in accuracy, simplicity and robustness and these characteristics are usually taken into account to choose the most suitable model for the application of interest.The existing models can be classified into quasi-static and dynamic models (Table 4), the latter can be further divided into parametric and non-parametric [63].Quasi-static use of the Bingham plastic model for MR fluids is suitable for the MR damper, but the non-linear behaviour under dynamic loading is not accurately described and the hysteresis phenomena is not considered.This issue is settled by dynamic models that can be sorted in parametric and non-parametric models [63].These modelling techniques describe the dynamic behaviour by a set of linear and/or non-linear springs and dampers.Bingham dynamic models, biviscous models, Bouc-Wen based models are the most commonly parametric dynamic models used.These require the individuation of a set of parameters by minimizing the error between the experimental results and the model simulated outputs.The process starts by assuming an initial values of the parameters, which are gradually updated by the simulations [63].Instead, the non-parametric models use analytical expressions derived from experimental data and MR damper physics analysis.The black-box model, neural network and fuzzy model are commonly reported in the literature.
The force F can be expressed as: where ẋ is the speed of the external solicitation, sgn is the sign function, f c is the frictional force, f 0 the force offset that takes into account the accumulator of a MR damper and c 0 is the damping coefficient.The model is characterized by a discontinuous jump in the damper force-velocity response (Figure 2).Although it can easily be implemented in a control system, it cannot describe the hysteresis behaviour, therefore some corrections should be done.In [74] the model errors have been evaluated by means of an experimental setup with actuators, displacement and force transducers.The following expressions have been introduced and used as a reference for quantitative evaluation of normalized errors: F is the simulated or predicted force, F exp represents the experimental or measured force, and Fexp is the measured force's average value along the T period.E t , E x and E v are normalized errors between model and experimental force in time, displacement and velocity domains respectively.Errors up to 15% for time, 4% for displacements and 13% for velocities have been reported for different current values.

Gamota and Filisko Model
In [43] an expansion of the Bingham model is proposed.The model is shown in Figure 3 and the force F is expressed by: where c 0 is the damping coefficient of the Bingham model and k 1 , k 2 , and c 1 are the springs and damper coefficients.The model parameters have to be evaluated for a specific set of frequency, amplitude and current excitation.The force-velocity behaviour is closer to the experimental one if compared to the Bingham model, but the numerical solution of ( 4) is more tricky and time consuming.In [75] the model error has been evaluated by means of displacement and force transducers: errors up to 20% for time, 7% for displacements and 30% for velocities model simulations have been calculated according to the expressions in (3).

Bi-Viscous Models for MR Dampers
These models are based on the hypothesis that the fluid assumes a plastic behaviour both in the pre yield and post yield area as introduced by [76,77,79,80].The force-velocity curve is divided into three sections depending on the velocity as described by the following equation where the post yield and pre yield damping coefficients, c ps , c pre respectively, and the force F yield represent the model parameters, also illustrated in Figure 4.The velocity v y ield can be calculated as [80] v yield = F yield c pre − c ps (6) An enhanced model that takes into account a further parameter in order to better represent the pre yield hysteresis is represented by the non-linear hysteretic bi-viscous model [76][77][78][79].The hysteresis loop is split in two sets of three equations.In Figure 5 are reported the four characteristic parameters (F yield ,c ps ,c pre ,v 0 ) where v 0 is the intercept of the force-velocity diagram with the horizontal axis.The authors of [80] proposed that the relationship between the four parameters and the current could be described with second or fourth order polynomials.This improved model cannot describe accurately the F(t) if the current value change, because the characteristic parameters must be evaluated for a particular MR damper excitation status.However, further works by [81,82] brought forward a lumped parameters bi-viscous model that describes the motion of the fluid considering: 1. an additional non linear function (quasi-steady velocity) 2. the control signal of the damper 3. the fluid compressibility and inertia 4. the piston's mass.This model proposed hyperbolic tangent functions to take into account the current effect on F yield , c ps , c pre parameters.The lumped parameters model accurately describes the MR damper behaviour and the effect of the current on certain parameters; on the other hand the number of characteristic parameters to be evaluated (i.e., 10) has more than doubled.
According to Figure 6 the damping force F is: where c 0 is the viscous coefficient, k 0 is the stiffness coefficient and z is the evolutionary variable related to the Bouc-Wen block: α, β, γ, and n need to simulate the hysteresis loop.
The eight model parameters (n, A, β, γ, x 0 , α, k 0 , c 0 ) must be estimated for every excitation frequency, amplitude and input current.The simple Bouc-Wen model [72] equations are easier to be numerically solved than the Gamota-Filisko models, anyway it can't simulate the roll-off effect of yield region shown in the experimental data.As in the Gamota and Filisko model, errors (3) are evaluated in [75] they range from 6% for displacements to 13% for velocities and 16% for time.

Modified Bouc-Wen Model for Large-Scale MR Dampers
This model takes into account also the shear thinning and inertial effects as well as the MR fluid stiction phenomenon [63,71,74].It has been observed that c and k decrease with excitation's amplitude but increase with input current.Instead, α increases with both x a and i c .The damper force is given by: where the z refers to the expression in (8), m takes into account the inertial and MR fluid stiction effect, k 0 is the stiffness of the accumulator, c( ẋ) is the post-yield damping coefficient that addresses the MR fluid shear thinning effect by means of a decreasing function it is assumed to be: where a 1 , a 2 and p are positive constants.This model is useful to simulate the large-scale MR dampers that are employed in practical engineering applications where the maximum damping force reach 200 kN.In [118] the modified Bouc-Wen hysteresis model has been experimentally tested through an electrohydraulic servo-controlled actuator, a load cell and a displacement transducer.It was found that the error E t ranges up to 18% but the errors E x and E v were not evaluated.

Spencer Dyke Model
In the model proposed by Spencer [11,[73][74][75]83,87] other parameters has been included to get more tidy results: a scheme is shown in Figure 7.The nonlinear roll-off effect, that can be seen in the force when the velocity is close to zero, should be taken into account by adding a dashpot in series with the Bouc-Wen model.On the other hand an additional spring is used to consider the accumulator's stiffness, if it is present in the MR damper.The equations are: The resulting model seems to better characterize a MR damper for a large range of operating conditions, as step voltage, random displacement at random voltage and random displacement at constant voltage.The characteristic parameters are: Errors E t , E x and E v of the model have been evaluated up to 4% for time, 2% for displacements and 4% for velocities.

Dominguez Model
In his work Dominguez [83] proposed a new Bouc-Wen model in order to better represent the hysteresis behaviour by considering the frequency, current and amplitude of excitation as variables.In fact the experimental results revealed that the hysteresis force depends on these quantities, as expressed in the following equations: Here, ω is the excitation's frequency, d 1 . . .d 4 are constants and z is the evolutionary variable.The characteristic parameters c 0 , α and F z0 vary exponentially for high current excitation and linearly for low current, therefore the following equations have been proposed: I c is the critical current that is peculiar for each MR damper.Then the characteristic model parameters are c 1 , c 2 , c 3 , c 4 , k 1 , k 2 , α 1 , α 2 , α 3 , α 4 , γ 1 , γ 2 , F z 01 , F z 02 , F z 03 , F z 04 .The evolutionary variable z could be updated as follows: The model can predict the hysteresis behaviour of the MR damper under any excitation conditions.In contrast, in the case of traditional Bouc-Wen models, the characteristic parameters have to be assessed for each set of amplitude, current and frequency of excitation.The model complexity is overcome by the possibility of using it in applications where the excitation conditions are variables.In this case, the model has errors of up to 22% for velocities (E v ), 10% for displacements (E x ) and 15% for time (E t ) [83].

Ali and Ramaswamy Model
In [44] the effect of displacement amplitude of the excitation is considered separately for c, k and α as a quadratic function of amplitude of the stroke (x a ): A nonlinear optimization is needed to obtain the values of the characteristic parameters with the following constraints c(x a , 0) ≥ 0, k(x a , 0) ≥ 0 and α(x a , 0) ≥ 0. This model is applicable when the excitation is known a priori so this may represent a limitation in many engineering applications.The normalized errors E t , E x , E v are not available in the literature.

Non-Symmetrical Bouc-Wen Model
In [45] a non-symmetric hysteretic response of Force/Velocity data has been reported, especially near the zero velocity (shifted Hysteresis).To take account of this phenomenon the expression of ż changed to ż = (A − (β + γsgn(z ẋ))|z| n ) ẋ (22) and: ẋ = ( ẋ − µsgn(x)) (23) So (22) became: where µ represent the scale factor of the velocity adjustment.The accuracy of the parametric models depends strictly on the characteristic parameters assessment.For this purpose, in the last decades several optimization approaches have been used as non-linear optimization [76,77,79], particle swarm optimization [21] and genetic algorithms [45].Since the characteristic parameters have a coupling effect on the hysteresis and are submitted to severe non-linearities, a Genetic Algorithm is more suitable for the parameters identification [45].Further models developed in recent decades are represented by the non-parametric damper models that combine the experimental data processing and MR damper physics to derive the analytical equations in order to predict the damper behaviour.Among them, the neural network models are the most commonly investigated.These models consist of a phase of training and prediction based on the comparison between a measured variable of interest (e.g., the damper force) with the output of the neural model.The training data of the neural model are generally based on a parametric model, as in the case of [116].The neural networks can also be useful for parameter identification, the inverse dynamic modeling and the MRD control.These models are not covered by the present work but can be further explored in [97,[107][108][109][110]116,117].
In Table 5 the most common MR models used for numerical implementation have been reported.Every model is characterized by advantages and drawbacks, depending on the application of interest.Applicable when the MR behaviour is not [45] symmetrical but a Genetic Algorythm approach is needed for the parameters evaluation In Table 6 a subjective evaluation, based on the scientific literature, concerning the reported models applicability fields is presented.For example, the Bingham model is not recommended for the behaviour prediction of MR dampers because it cannot describe the hysteresis behaviour.On the other hand, the Gamota and Filisko model, that is an extension of the Bingham, has been developed to enhance this gap.Bouc-Wen based models could be implemented for both civil and automotive applications because they have been tailored for these purposes.Moreover, the authors propose a preliminary evaluation about the application of MR damper models in microvibration control, wherein vibration amplitude are lower than 2 mm.To this regard the Bouc-Wen models seem to be the most promising.
In all the studies where a model is proposed, experimental measurements are usually employed for its validation.Nevertheless, it can be observed that a detailed uncertainty analysis of the results is often not conducted.However, it is clear that the accuracy of the measured data influences dramatically the model validation and performances.For this reason in the next section an overview of the main uncertainty sources in the experimental setup for magnetorheological dampers characterization is described and discussed.

Characterization Methods for MR Dampers
Performances of dampers are usually studied on the basis of: (1) the provided force vs. its displacement and (2) the provided force vs. speed [84].These above characteristics are commonly considered when more dampers have to be compared and in logical models evaluation [85].In particular, the force-velocity graphs come from acquired data during a damper testing by means of a "shock dyno" or "damper dynamometer": a device that extends and compresses a damper at known speeds, providing measurements of the produced damping forces.
In particular, sinusoidal testers or, more generally, a Fatigue Testing Systems (FTS) represent two experimental methods to characterize a MR damper [88][89][90].The FTS may be set in a displacement controlled mode for the performance assestment using harmonic excitations over a large spectrum of amplitude and frequency.In fact the sinusoidal motion of the FTS is usually achieved by a crank [85] or a shaker [86].The crank quality is also described in other applications [138,139].On the other hand, the FTS should be integrated with proper sensors and an acquisition system (Table 7).The main sensors of a sinusoidal tester are: (1) a load cell (2) a sensor for the displacement or the acceleration of the damper; (3) a temperature sensor [87].Figure 8 shows an example of setup for damper testing: it is made by a frame withstanding an electric motor with a drive belt and pulleys that revolves a slider-crank mechanism linked to the damper shaft by means of a linear bearing.The damper piston can move up and down with different stroke lengths (e.g., by means of different crank configurations) and variable rotation speeds.The damper force is measured by the load cell [84].The performances of the testing machine are affected by several issues, i.e., the actuator friction, the frame and specimen stiffness, the machine resonance, the moving mass and its deflection, the excitation frequency and its waveform shape.Some further consideration should be undertaken to estimate the uncertainty due to static calibration and static testing.In other fields, researchers have studied the effect of dynamic calibration to dynamic devices performances [140,141].Static calibration uncertainty can be estimated following different standards [93,94], it can be evaluated from the Root Sum Squared (RSS) which takes into account different contributions, among which the main are the repeatability, hysteresis, resolution and uncertainty of the calibration force.Due to static calibration and resolution, the static testing uncertainty can be evaluated from the error.On the other hand, the uncertainty estimation of the FTS during its use can be very complex and the number of contributions can be high and case-dependent: in Table 8, an example of the main contributions is reported.The stability refers to the hydraulic supply (e.g., an error occurs when flow or pressure in the system are not adequate to provide the forces required) and depends strongly on the test condition, i.e., type of testing and component being tested.If manufacturer's specifications are not satisfied, power fluctuations may also significantly affect the test results and their uncertainty is case-dependent, as is the specimen alignment and preparation.Typically, during an MRD performance try-out, a certain current value is set and maintained.Then for each current value, numerous force vs. velocity curves are acquired in order to analyse the whole damper force extent.However, the MR damper measured force vs speed curve looks like a magnetic material hysteresis cycle [85].The parameters in a dynamic model of a MR damper have to be determined by means of the measured data in order to make the theoretical model able to fit the hysteresis loop.The parameter estimation can be considered a multi-dimensional numerical optimisation problem where the parameters are the decision variable and the optimisation algorithm should be introduced to adjust the variable value in order to achieve the minimum of an objective function (performance criterion).The most common objective functions are the quadratic error J 1 and the root-mean-square error J 2 (RMSE) between the measured and simulated damping forces.In particular the quadratic error J 1 can be expressed as: where k = 1, 2, . . ., N is the number of experimental samples for the measured force F measured .On the other hand the RMSE J 2 is expressed as: The strategy is to estimate the parameters values that make J 1 (or J 2 ) minimum.

Conclusions
Recently, a growing interest has been observed in magnetorheological (MR) fluids, and MR dampers are widely used in many fields of engineering and science.Their models are not always adequate to fit the experimental data, despite their complexity (e.g., great numbers of parameters).Furthermore, the parametric models' capability to accurately outline the typical hysteresis shape in the force-velocity diagram depends significantly on the parameters identification method.Generally, models with more complex physics and several characteristic parameters are more accurate at predicting the behavior of the damper, as long as the parameters' assessment is properly conducted.Most of the reported models successfully describe the force displacement MRD behaviour.On the other hand, currently, there is no simple dynamic parametric model that can provide the MRD behaviour thoroughly enough.Moreover, from the present work, it emerges that the MRD models were mostly investigated and designed for civil engineering applications (seismic and building protection) and later adapted for vibration control in the automotive field.Some models (Bouc-Wen based) are more suitable in applications where the damper operating conditions varying over time and, to this regard, studies on MRD real-time control have been successfully carried out.Some models are oversimplified and the result is not suitable for real-time control purposes (Gamota-Filisko,Bingham).Although many models have been developed for MR dampers, currently, only a few works have explored the possibility of using a MRD in aeronautical industry and bio-engineering.The authors believe that another promising and interesting applications may be seen in the control of microvibrations, e.g., in the whole body vibration protection or actuation systems.Hopefully, future efforts and researches might fill this gap in the literature and improve the techniques for a better control and characterization of MR Dampers.

Figure 1 .
Figure 1.Scheme of the Bingham model.

Figure 2 .
Figure 2. A typical measured characteristic of a MR Damper Force-Velocity response.

Figure 4 .
Figure 4. Typical Force-Velocity diagram of a Bi-viscous model.

Figure 8 .
Figure 8. Scheme of a testing setup used for dampers characterization.

Table 2 .
Examples of drawbacks in MR damping use.

Table 5 .
Main parametric dynamic models of MR dampers.

Table 6 .
Major applications of parametric dynamic models for MR dampers.Key: R = Recommended, U = Unfeasible, M = Request further research, P = Promising.

Table 7 .
Measurement ranges in MR dampers characterization.

Table 8 .
[93]ples of the main contributions to uncertainty in testing machine during use[93].