Sequential neural network model for the identification of magnetorheological damper parameters

Magnetorheological (MR) dampers exhibit a complex nonlinear hysteresis which makes the modeling of their behavior with parametric or non-parametric models to be challenging. In case of parametric models, the generalization of the parameters identified for a particular excitation is difficult and requires high computation costs. On the other hand, non-parametric models are considered as black-box type with no association to physical phenomena. The objective of this study is to propose a new identification model combining the merits of a parametric model and neural network paradigm. The proposed model is a parametric type which exploits an algebraic model with a hyperbolic tangent hysteresis, while a series multilayer-perceptron (MLP) neural networks are used to identify the model parameters under different excitation conditions. This approach not only preserves the physical meanings of the model parameters but also prompts generalization to common excitation conditions. The data for training the MLP neural networks were generated from a test program on a RD-8041-1 MR damper covering a wide range of input conditions. Results show that the proposed sequential neural network model not only increases the accuracy of the predicted MR damper force but also exhibits higher robustness and better consistency under different excitation conditions than a conventional algebraic model.


Introduction
Magnetorheological (MR) dampers are categorized as semiactive control devices due to the unique rheological characteristics of the contained fluids.These adaptable dampers are widely utilized to mitigate the dynamic responses of Original content from this work may be used under the terms of the Creative Commons Attribution 4.0 licence.Any further distribution of this work must maintain attribution to the author(s) and the title of the work, journal citation and DOI.
structures and mechanical systems.The superior features of MR dampers, such as controllability, low input power requirement compared to active devices, fast response to the applied current, capacity to generate large damping force, high reliability and stability, and wide operating temperature range motivated researchers who conducted numerous studies on their usage in vibration control [1][2][3][4][5][6][7][8][9].Despite their qualities, MR dampers exhibit an inherent complicated nonlinear behavior which poses many challenges for accurately modeling of their response [10].Most existing models can typically be classified as parametric and non-parametric models [11].
Parametric models are conventionally used to describe the behavior of MR dampers.These models idealize a MR damper using combinations of springs, dashpots, and other physical mechanical components [12].In the past decades, many studies have been conducted to develop parametric models that accurately represent the hysteresis of MR dampers.Besides Bingham model [13], modified Bingham models have also been proposed, of which either a linear spring component was added in series with the Bingham model [14] or the damping coefficient and plastic force were modified to be dependent of current level [15].However, the MR damper behavior in the pre-yield region could not be accurately described by these models.Further, some model [14] had extremely stiff governing equations so the numerical solution to which required an extremely small time-step in the order of one microsecond [10].Snyder et al [16] proposed a nonlinear bi-viscous model to represent the hysteresis of MR dampers, but the model could not simulate the roll-off phenomenon when the damper velocity was near zero.To improve the modeling of strong nonlinear hysteresis of MR damper, Spencer et al [10] developed a simple Bouc-Wen model in which the Bouc-Wen hysteresis was exploited as a hysteresis component in parallel with a spring and a dashpot.Although this model could describe the post-yield region and the hysteresis behavior well, it failed to model the roll-off region in the nonlinear force-velocity response.The modified Bouc-Wen model [10] was proposed for better representation of the roll-off occurrence at low velocities to overcome this difficulty.First-degree and third-degree equations with respect to applied current were established respectively by Spencer et al [10] and Yang et al [17] to adjust the parameters of the modified Bouc-Wen model for small and large-scale MR dampers.Furthermore, Dominguez et al [18] developed a current-frequency-amplitude-dependent Bouc-Wen model.The resulting modified Bouc-Wen model [10] necessitates the identification of 14 constitutive parameters.The identification of such a large number of parameters under different excitation conditions is time-consuming.Therefore, efforts have been made by numerous researchers to either reduce the number of parameters in the models based on Bouc-Wen hysteresis [19] or replace the Bouc-Wen hysteresis component with another evolutionary hysteresis such as Dahl [20] and LuGre [21].The inclusion of evolutionary hysteresis (Bouc-Wen, Dahl, and LuGre) makes it cumbersome to identify the parameters as the governing differential equation of the hysteresis variable must be integrated at each step.To address this issue, the feasibility of using algebraic functions, such as the hyperbolic tangent function [22], the arctangent function [11,23,24] and the sigmoid function [25,26], has been explored to replace evolutionary hysteresis variables as hysteresis components in the MR damper parametric model.S ¸ahin et al [11] carried out a comparison between the accuracy of models based on ordinary nonlinear differential equations and those who employed algebraic functions.The results showed that the MR damper response predicted by algebraic models agreed well with the experimental data [11].Compared to differential parametric models, algebraic models were computationally efficient in identifying the parameters and thus have a decisive advantage to be implemented in control-related problems [11].Kwok et al [22] proposed a simple algebraic model in which the hyperbolic tangent function was used as a hysteresis variable in parallel with a viscous dashpot and a spring.The model could accurately describe the forcevelocity response of the MR damper.The authors generalized the proposed model to include the fluctuation of current levels.This generalization reduced the complexity associated with the Bouc-Wen hysteresis and lowered the number of unknown parameters to six, with each of them representing a specific physical behavior in the hysteresis [22].Zhang et al [27] performed a sensitivity analysis on the parameters of the algebraic model proposed by Kwok et al [22], and further extended the form of the model parameters to be dependent on the current and excitation characteristics.
Non-parametric models treat the physical characteristics of a MR damper as a black box [28].Due to excellent nonlinear regression and learning capability, numerous studies have utilized neural networks to simulate the behavior of MR dampers in non-parametric models.Chang and Roschke [29], Ekkachai et al [30], Tudón-Martínez et al [31] and Han et al [32] used the multilayer perceptron (MLP) as a feed-forward neural network (FNN) approach to model the dynamic behavior of a MR damper.To improve the accuracy of the FNN model, Wei et al [33] applied a Hilbert transformation to develop the input layer of the network by extracting the instantaneous variables such as amplitude and frequency and by adding them to the piston displacement, velocity and current.Khalid et al [34] identified the dynamic behavior of a small-scale MR damper using a recurrent neural network (RNN) model and compared its performance with the FNN model by Wang and Liao [35].Results showed that although the accuracy of the FNN model in predicting the MR damper force improved, the RNN model did not require monitoring the online force and providing feedback to the input layer.Du et al [36] suggested that the radial basis function could be more efficient in modeling the characteristic of MR dampers than the MLP networks.Boada et al [28] proposed a recursive lazy learning method based on a neural network.The model could be adapted to new datasets because of its recursive feature [28].Although non-parametric models demonstrated the flexibility and capability of readily learning, the parameters used in those models did not have clear physical meanings.Further, when used in control problems, a large dataset will be needed in the training phase.
As discussed earlier, the variables in parametric models can be associated with the physical behavior of a MR damper, but generalizing a parametric model under various excitation conditions necessitates a high computational cost.Maintaining an acceptable accuracy of the generalized parametric models under various excitation inputs remains a challenging task.On the other hand, although non-parametric models can accurately predict the MR damper force in the domain of the dataset used for neural network training, a large set of collected data is needed for training and their parameters do not have specific physical meaning.Nevertheless, models based on neural networks can be adaptive with the addition of new data.
The current study aims at developing a new identification methodology for MR damper parameters by combining the merits of a parametric model and neural network strategy.The proposed model is based on the hyperbolic tangent model while the identification of the model parameters are performed through a series of connected MLP neural network models referred to as sequential neural networks (SNNs) model in this study.The proposed approach would not only allow the model to clearly reflect the physical behavior of a MR damper but also to identify the model parameters with manageable computational cost than training a neural network model to directly describe the nonlinear hysteresis behavior of a MR damper.In addition, with the superior capability of neural networks in handling nonlinear regression, the damper parameters can be identified more accurately, and the model may be generalized for a wide range of excitation conditions.To collect the experimental dataset for training the neural networks for model parameter identification, the MR damper was tested under 162 different conditions.

MR damper specimen
MR dampers can be considered as controllable semi-active actuators, the generated resisting force may be adjusted with minimum power requirements.Although the configuration of MR dampers is similar to typical fluid-based dampers in which the damper fluid flows through a narrow gap during operation, the fluid used in MR dampers contains micron-sized magnetizable iron particles dispersed in the carrier oil.This fluid has a unique feature; when an external low-powered magnetic flux is introduced, the dispersed particles in the MR damper fluid would be aligned in a chain-like form parallel to the field.This changes the damper's fluid rheology from a free-flowing state to a semi-solid state within a few milliseconds and thus alters the mechanical properties of the fluid.This phase transition is reversible and tunable.Two main operational modes could occur: i.e. the flow mode and the direct shear mode.A combination of these two modes is also possible.
The MR damper utilized in the current study is the RD-8041-1 model manufactured by the LORD Corporation.The stroke length of the damper is 74 mm.The maximum durable current level under 30 s of continuous excitation is 1 A. Under an input voltage of 12 V DC, the electrical resistance of the coil varies between 5 Ω at ambient temperature to 7 Ω at an operating temperature of 71 • C. A schematic illustration of the MR damper model used in this study can be found in figure 1(a).

Experimental setup
To obtain the characteristics of the MR damper hysteresis behavior, a series of cyclic loading tests were performed using an MTS machine (Model No. 831.50S) under displacementcontrol loading mode, as shown in figure 1(b).The force and velocity capacity of the MTS machine are ±10 kN and ±200 mm s −1 , respectively.The Multipurpose Testware software was used to acquire the displacement and force signals.
A harmonic excitation is applied to the damper piston under varying current levels.The amplitude and frequency of the excitation and the current level used in the test are listed in table 1.The sampling rate is fixed at 100 Hz.Considering the maximum allowable velocity of the MR damper piston is 200 mm s −1 , 28 cases are designed for each current level based on different excitation amplitude and frequency combinations.This yields a total of 162 cases covering all six current levels.

Experimental results
The main experimental results, in the form of the hysteresis behavior of the MR damper, are presented in this section.As illustrated in figure 2, the force-velocity hysteresis loop follows a counterclockwise path, of which the upper and lower parts of the loop correspond respectively to a negative and positive acceleration of the damper piston motion.In general, the hysteresis of a MR damper can be divided into two distinct regions of the pre-yield and post-yield zones by the rolloff effect.Roll-off occurs when the damper piston velocity is close to zero and the piston starts to reverse its motion direction.Under such a condition, the fluid between the piston and cylinder would be blown-by or bled, resulting in a considerable drop of the generated damper force [10].The pre-yield region demonstrates a strong hysteresis characteristic which stems from the MR fluid viscoelastic property, while a plastic behavior can be observed in the post-yield zone [12].The preyield and post-yield zones are clearly discernable by examining the force-displacement (F-D) hysteresis.The two vertical lines of F-D hysteresis describe the damper performance in the pre-yield area, and the top and bottom parts of the curve reflect the plastic behavior of the MR damper in the post-yield region.Note that the F-D hysteresis loop of the MR damper progresses along a clockwise path.

Effect of the excitation frequency.
This section investigates the MR damper performance under varying loading frequency.Figure 3 shows the hysteresis response of the MR damper under a loading corresponding to a stroke amplitude of 15 mm and a current level of either 0.4 A or 0.8 A when the excitation frequency varies from 0.5 Hz to 2 Hz.As can be seen from figures 3(a) and (c), an increase in the excitation frequency would expand the pre-yield zone and widen the velocity range associated with roll-off, leading to a larger force-velocity hysteresis loop.The post-yield region would extend to a larger velocity but at a decreased slope.The generated damper force increases with the excitation frequency, and this phenomenon is more pronounced at the high current level of 0.8 A. In addition, it is noticed that while the preyield and post-yield regions in all four studied excitation frequency cases can be easily distinguished at the current level of 0.4 A when the current level increases to 0.8 A, the post-yield zone becomes very limited at the frequency of 0.5 Hz and thus distinct from the roll-off region.This implies that a currentdependent threshold velocity exists such that when exceeded, the damper operation would cover both the pre-yield and the post-yield regions.respectively.The observation indicates that a higher excitation amplitude would result in a larger hysteresis, and such an influence is more considerable at a lower current.Further, when the excitation amplitude is small, the maximum damper force is found to be more sensitive to the amplitude effect.The pattern of the force-velocity curves in figures 4(a) and (c) shows that the slope of the force-velocity hysteresis in the post-yield region gradually decreases as the excitation amplitude becomes larger.This phenomenon is more visible at the current level of 0.6 A. However, compared to the frequency effect shown in figures 3(a) and (c), the width of the pre-yield region is much less sensitive to the change in the excitation amplitude.

Effect of the current level.
The influence of the current level on the MR damper behavior is depicted in figure 5 for three different excitation conditions: (i) A = 5 mm, f = 2.5 Hz; (ii) A = 10 mm, f = 1.5 Hz; and (iii) A = 25 mm, f = 1.0 Hz.Results show that while the resistive force generated by the MR damper increases with the current level, it gradually approaches saturation at a high current level.For instance, figure 5(d) depicts that when the current increases from 0 A to 1 A with an increment of 0.2 A, the difference between the forces generated at two successive current levels are 294, 314, 311, 159, and 111 N, respectively.Moreover, compared to the effect of the frequency and amplitude, the current level is observed to have a more significant impact on the amount of resistive force generated by the MR damper.As indicated by the force-velocity hysteresis loops shown in figure 5, with the increase of the applied current level, the slope of the postyield and pre-yield regions increases whereas the occurrence of roll-off is delayed to a slightly lower velocity.The postyield region becomes wider at high current levels, which is more obvious when the excitation amplitude is smaller.This phenomenon is induced by the presence of an accumulator at the bottom of the damper cylinder as shown in figure 1(a).A larger excitation amplitude would cause a more sizable change in the volume between the damper piston head and the diaphragm which corresponds to a larger difference in pressure between the top and bottom of the piston head.This results in the accumulator compensating more variation of pressure into the volume to prevent the possible occurrence of cavitation in the cylinder.Meanwhile, when the excitation amplitude remains constant and low with the increase of current level, the gap where the MR fluid flows through (as shown in figure 1(a)) would be narrower; thus, the accumulator would be required to compensate much more pressure difference.

Tangent hyperbolic model
The experimental results presented in the previous section suggest that the resistive force generated by a MR damper depends on the current level, the excitation frequency, and the excitation amplitude.Therefore, the hysteretic characteristics of the MR damper can be expressed in terms of these three variables as where F (t) is the damping force, I is the applied current, x and ẋ represent the displacement and velocity of the damper piston, respectively.The establishment of a model for a MR damper entails modeling the hysteresis of the damper force.Parametric models can be categorized according to the type of hysteresis component.Kwok et al [22] proposed a model based on a tangent hyperbolic function to estimate the force generated by a MR damper, as portrayed in figure 6.The model consists of three components placed in parallel: (i) a linear spring defined by its stiffness describing the accumulator effect; (ii) a linear viscous dashpot representing the viscoelastic response of the MR damper in the post-yield region; and (iii) a tangent hyperbolic function to model the hysteretic behavior of the damper.The resisting force can thus be expressed as: where c and k are respectively the viscous and stiffness coefficients, α is the hysteresis scale factor, F 0 is the damper force offset, and z is a hysteresis variable defined by equation ( 3), as follows: where parameters β and δ control the hysteresis loop shape.Figure 6 illustrates that the term cẋ represents the slope of the post-yield region with a larger c value for a steeper slope.The force generated by the spring, kx, creates a horizontal ellipse in the force-velocity hysteresis, resulting in an opening at the ends.A larger k value corresponds to a larger opening at the ends.The hysteresis variable z defined in equation ( 3) by a tangent hyperbolic function is always bounded between −1 and 1.Thus, the parameter α is used to scale the nondimensional variable z to the MR damper resistive force range.Furthermore, β and δ can determine the slope and width of the force-velocity hysteresis, respectively.As shown in

Identification of model parameters
In the algebraic model defined by equation ( 2), a tangent hyperbolic function is used to represent the hysteresis behavior of the damper.Thus, the associated parameters can be determined by applying curve fitting to the experimental data rather than using optimization approaches.In addition, applying a curve fitting to the hysteresis based on tangent hyperbolic function is computationally less demanding than other conventional hysteresis models, such as the Bouc-Wen hysteresis operator with an evolutionary variable [22].The dataset from the 162 experimental cases is used to calculate the optimal values of the model parameters (c, k, α, β, δ, and F 0 ) in equations ( 2) and (3), for each hysteresis curve, by employing the nonlinear least-squares method in the MATLAB Curvefitting Toolbox.The curve-fitting procedure yields 162 optimal values for each of the six parameters from the experimental data.Figure 8(a) reveals that the lowest and the average R 2 values of fitted curves to the experimental data are within the interval of 0.9963 and 0.9999, respectively, meaning that the damper parameters were estimated accurately, especially for cases with higher frequencies.regression is used to fit the identified 162 optimal values for different current levels, the trend of parameter variation cannot be accurately approximated.As shown in table 2, with the exception to the parameter α, the fitted curves for the rest of the parameters do not show acceptable accuracy, especially in the case of hysteresis parameters, β and δ.
Since it is challenging to measure the vibration amplitude and frequency of the damper in real time, Zhang et al [27] proposed to use the maximum damper velocity to represent the excitation characteristic to address this issue.Moreover, a direct calculation of the maximum damper velocity within a cycle is difficult in the case of random excitation.Thus, Wang et al [25] proposed to express the maximum velocity of the MR damper piston as a function of the instantaneous excitation displacement, velocity, and acceleration as: where x, ẋ and ẍ are the displacement, velocity, and acceleration of the MR damper piston motion, respectively, and v m is its maximum velocity.
Zhang et al [27] observed that the parameter α could have a linear relation with the current level and be represented by a quadratic polynomial equation of the maximum velocity to yield a good fit.Thus, parameter α can be described as follow: Furthermore, based on the observation that the parameter β decreases nonlinearly when the velocity amplitude increases and the effect of the current level on parameter β can be neglected, expression (6) of parameter β was proposed [27] where β 1 and β 2 are constant coefficients.The parameters δ and c can be expressed as linear functions of the current level as where δ 1 , δ 2 , c 1 and c 2 are constant coefficients.
Zhang et al mentioned that since the sensitive values of parameters k and f 0 to the algebraic model were low [27], the mean values of these two parameters can be used in the algebraic model.Referring to table 2, it can be seen that Zhang et al's model [27] can only improve the estimation of parameter β, whereas the prediction of the other parameters remains an issue.

Sequential neural network model
As shown in figure 9, the identified parameters of the tangent hyperbolic-based algebraic model manifest a substantial nonlinear variability; analytical expressions could not provide accurate prediction.Alternatively, neural networks are recognized to be an efficient method to model problems with strong nonlinearity and conduct regression analysis.This study proposes an application of a neural network methodology to estimate the parameters of the MR damper algebraic model.In this regard, a SNN model containing six neural networks, corresponding to the six model parameters (c, k, α, β, δ, and F 0 ), is developed to predict the parameters of the hyperbolic algebraic model.Artificial neural networks have been widely used to model complex nonlinear problems due to their strong learning capability and robustness.Multilayer FNNs, such as MLP networks, are well-known approaches to estimating any continuous functions from a dataset [37].The MLP neural networks are developed to identify the parameters of the tangent hyperbolic model expressed by equations ( 2) and ( 3).The objective is to improve the accuracy of the MR damper model.The proposed approach would not only exploit the capabilities of neural network methodology in predicting nonlinear functions but also allows retaining the physical meaning of parameters in the algebraic model.

MLP model.
A MLP network is a type of backpropagation FNN that includes an input layer, one or more hidden layers, and an output layer.Each layer is a set of neurons, and the input signals propagate layer by layer through the network [38].Figure 10(a) depicts the architecture of a fully connected MLP neural network consisting of n neurons in the input layer, m hidden layers and one neuron for the output layer.The sigmoid function is used as the transfer activation function in the hidden layers to ensure a smooth variation requirement and model the inherent nonlinearity of the expected model.A linear activation function is considered for the output layer characterizing the attempted regression problem.
The Levenberg-Marquardt (LM) backpropagation algorithm is used to train the neural network and optimize the weights and biases values.The LM algorithm addresses the difficulties of the Gauss-Newton method and gradient descent algorithm [39].To estimate the initial weights, LM uses a gradient descent type algorithm and then switches to Gauss-Newton until the cost function approaches the minimum value.Once a minimum cost function is reached, the LM returns to the gradient descent algorithm to linearly boost the model's accuracy [40].

Feature selection for the input layer.
To develop a neural network model for predicting the algebraic model parameters, input layer variables need to be determined first.As discussed in section 2.3, it was observed that the excitation frequency, amplitude and current applied to the piston's coil could all affect the hysteresis behavior of a MR damper.Determining each variable's importance to the MR damper's resistive force is challenging.A feature selection approach was utilized to rank the importance of variables.Three algorithms are commonly employed to evaluate the importance of input parameters in regression problems, i.e. the F-Test method, the minimum redundancy and maximum relevance (MRMR) algorithm and the ReliefF algorithm.These methods are briefly explained as follows.
• F-Test algorithm: F-Test is a statistical test that determines the importance of a variable by calculating the f -score.The f -score of each feature is the ratio of 'between variables variance' and 'within variables variance'.The features with higher f -scores are considered more important and are assigned higher ranks [40].• Minimum redundancy maximum relevance (MRMR) algorithm: the mutual information (MI) is a measurement used to determine the level of similarity or correlation between variables.The redundancy measure utilizes the value of MI between two features (input variables) to ensure whether the two features have a correlation.A large MI indicates redundancies between two features.The relevance measure usually uses the MI between each input variable and the target output, where a large value of MI indicates a strong correlation between the considered feature and the target output.The MRMR algorithm uses both the relevance and redundancy measures to select a subset of bene-ficial features to achieve a lower MI between the set of all features and a large MI with respect to a target output [41].• RReliefF algorithm: Kira and Rendell [42] proposed the original formulation of the ReliefF algorithm to estimate the quality of features based on how distinctively their values are between close instances [35].The ReliefF algorithm method does not examine the correlation between the features, thus, it is not possible to remove redundant features [42].To utilize ReliefF in regression problems, Kononenko and Robnik-Sikonja [43] adapted it called RReliefF.Six variables were used in the input layer of the proposed MLP network.They include the current level, the excitation frequency, the maximum velocity, and different algebraic combinations of these three variables.Based on the proposed input layer variables provided by the three algorithms, the neural network models of the MR damper model parameters were trained 50 times.In each training attempt, the R 2 value is determined based on the whole data set, of which 70% is used for neural network training, 15% for validation, and another 15% for testing.The average R 2 values were used as an index to evaluate the improvement made by the selected input layer variables for the network.The feature selection procedures used for the variable assessment and the corresponding results yielded from each method are listed in table 3.In the scenario 'All' in table 3 the input layer contains six variables, which are I, f, v m , I × f, I × v m , and f × v m .
• Parameter c: Results reported in table 3 indicate that the variables selected based on MRMR and F-Test algorithms are the same but the RReliefF algorithm chooses different features for the input layer.Although the MRMR shows an enhancement in the R 2 average, the difference is not sizable when compared by using three independent variables, I, f and v m .Moreover, the maximum R 2 of the three independent variables case is close to the referred three first-ranked features by the MRMR algorithm.Therefore, the best accuracy for estimating parameter c can be achieved when the input layer variables are the current I, the excitation frequency f, and the maximum velocity v m of the damper piston.• Parameter k: For the estimation of parameter k, the use of the three independent variables, I, f and v m results in a better R 2 than the input layers suggested by the three algorithms.• Parameter α: For the estimation of parameter α, all three methods indicate that the two features, I and I × v m should form the input layer.• Parameter β: The MRMR algorithm suggests that to have a more accurate prediction of parameter β, the input layer should consist of the three features, f × v m , I and I × v m .• In the case of parameters δ and F 0 , the three feature selection algorithms could not propose a combination of features that is significantly better than the three independent variables, I, f and v m .Moreover, including all six variables in the input layer does not show a sizable improvement in accuracy compared to the three independent variables.

Architecture of trained neural networks.
To estimate the parameters associated with the MR damper algebraic model expressed in equation ( 2), the structure of the proposed Table 3.     10.For instance, the neural network used to estimate parameter c is formed by three neurons in the input layer, i.e.I, f and v m , three hidden layers with ten neurons each, and the identified parameter c as the output neuron layer.Splitting the dataset between the training and testing parts is an important step in designing a neural network.A parametric study is needed to find the optimal ratio of data used in training and testing to prevent overfitting of the network and generalize the model to new data.Four ratios 80:10, 70:20, 60:30 and 50:40 were considered.For each of these ratios, 10% of the dataset was also allocated for the validation process.Joseph [44] used the number of input variables, p, to estimate the optimum weight between the training, validation and testing data.An optimal splitting ratio of p : √ p : √ p + 1 for training, validation and testing was proposed [44].In the present work, an empirical statistical approach was used where the neural network of each parameter was trained 100 times to find the optimal ratios.The sizes of the training and testing data sets were randomly selected from the available experimental dataset each time.The R 2 ratios of the training and testing dataset were collected from the generated statistics.The boxplots of R 2 value for each parameter are shown in figure 11.The optimal splitting ratio for each parameter is chosen based on the stability and accuracy in the training and testing results.The optimal ratio corresponds to the highest median R 2 value, and the difference between the 25th and 75th percentiles in the results should be acceptable.For instance, figure 11(d) shows that the trained network for identifying parameter β can have a large median R 2 when the ratios of training/testing are 80/10 and 50/40.Meanwhile, better stability and accuracy under the splitting ratio of 80/10 in testing is observed.Thus, a splitting ratio of 80/10 appears more appropriate to train the neural network for parameter β.The same criteria were applied to select the optimal data splitting ratio in the neural network training for the other five parameters.Table 4 shows the optimal splitting ratios (training:validation:testing) for the proposed MLP neural networks for the six parameters.Once the splitting ratio of the parameters was determined, the trained networks, which yield the highest R 2 in testing, were selected as the final regression models.The R 2 values of the testing data and the whole dataset are listed in table 4. Except for parameter F 0 , the other parameters all have acceptable R 2 .It means that parameter F 0 has weak correlation with variables I, f and v m .The development of regression procedure for the six neural networks was carried out using the MATLAB neural network toolbox [45].

Comparison of sequential neural network and algebraic models
Comparative studies were conducted between the proposed SNN model and the algebraic model [27].Both models include the effect of the applied current, excitation frequency and maximum velocity of damper piston.Figure 12 shows a better match between the proposed SNN model and the experimental data than the algebraic model.It also shows that the SNN model can describe the force-velocity (F-V) hysteresis behavior of MR damper in both the pre-yield and the post-yield regions accurately.For instance, as shown in figures 12(a) and (c), the pre-yield regions in the MR damper hysteresis curves predicted by the algebraic model do not fit accurately with the experimental data at a current of 1.0 A. This suggests that equations (7)   Figures 13(c) and (d) show the R 2 values of all the testing cases for the two models.In most cases, the R 2 values of the SNN model are close to 1 while it is not the case for Zhang et al's model [27].The improvement shown by the SNN model is particularly noticeable at zero-current conditions.
To further evaluate the efficacy of the two models, a comparison of the average and standard deviation of R 2 is portrayed in figure 14.The average R 2 represents the accuracy   of the model whereas the standard deviation describes the stability of the prediction accuracy.Figures 14(a)-(c) show the variation of the average and standard deviation of R 2 values with respect to the current level, the excitation frequency, and the excitation amplitude, respectively.In all three sets of comparison, the proposed SNN model is observed to manifest a consistency in the average R 2 value, under varying excitations conditions which is nearly 0.99.As shown in figure 14(a), compared to the algebraic model by Zhang et al [27], the proposed SNN model offers a 31.3%improvement in the average R 2 at zero-current and an almost 1.6% improvement when the current is non-zero.Besides, the SNN model not only shows a stable accuracy by keeping the standard deviation of R 2 less than 0.002 under different current levels, but also reduces it in a range of 32.4% to 90.4% from Zhang et al's model [27].When frequency increases, the average R 2 by Zhang et al's model [27] decreases from 0.95 to 0.91 (figure 14(b)), while the proposed SNN model remains to show a consistency in accuracy with an average R 2 being over 0.99 for all studied frequencies.Further, while the standard deviation of R 2 by Zhang et al's model [27] increases gradually from low frequency to high frequency, that of the SNN model remains less than 0.007.Figure 14(c) shows that Zhang et al's model yields smaller average R 2 values than the SNN model, especially at lower excitation amplitudes.Moreover, the standard deviation of R 2 for the SNN model is always less than 0.017 and shows consistency over the studied amplitude range, whereas those by Zhang et al's model [27] are larger.The SNN model is found to improve the standard deviation of R 2 values by more than 95% at lower amplitudes.
Results of the above comparison clearly indicate that the proposed SNN model can accurately predict the behavior of MR dampers in terms of the force-velocity and forcedisplacement hysteresis curves.Moreover, compared to the existing algebraic models, the model demonstrates consistent and stable performance under different currents and excitation conditions.

Conclusion
A SNN model is proposed in the present study to model the nonlinear hysteresis behavior of MR dampers.The proposed SNN model constitutes six feedforward MLP neural networks, each of which is used to predict one of the six MR damper parameters in the tangent hyperbolic model.A RD-8041-1 MR damper manufactured by the LORD Corporation was tested under different harmonic sinusoidal excitation conditions and currents.The measured hysteresis results of 162 testing cases were fitted to the hyperbolic tangent model and the identified six damper parameters were collected to generate a database for training the SNN model.The optimal splitting ratio of the data used for training and testing was determined based on a statistical approach.The MR damper force predicted by the proposed SNN model was compared with the experimental data and the prediction by the existing algebraic model.The following are the conclusions from the current study: • The proposed SNN model retains the benefits of both parametric and non-parametric models.Not only does each of the model parameters represent a specific physical phenomenon associated with the hysteresis characteristic of a MR damper, but the capability of MLP neural networks in tackling problems with a nonlinear regression nature allows it to accurately predict the parameters in the tangent hyperbolic model.• Compared to the excitation frequency and amplitude, the current level is observed to have a more significant impact on the resistive force generated by the MR damper.However, the trend of the six identified parameters with respect to the current indicates that the current variable alone cannot be used to estimate the model parameters.The excitation conditions such as frequency and maximum velocity of MR damper piston, should also be selected as the main input variables for the proposed SNN model.• The manageable computational cost and the adaptability to different operational conditions make the proposed SNN model a competitive candidate to be implemented in realtime semi-active control strategy for MR dampers.

Data availability statement
The data cannot be made publicly available upon publication because they are not available in a format that is sufficiently accessible or reusable by other researchers.The data that support the findings of this study are available upon reasonable request from the authors.
Figure 4  portrays the effect of the excitation amplitude on the hysteresis behavior of the MR damper at current levels of 0.2 A and 0.6 A when the excitation frequency is 1 Hz.It can be seen from

Figure 1 .
Figure 1.MR damper specimen (a) schematic of MR damper structure; (b) experimental test set-up.

Table 1 .
figures 4(b) and (d) that for the applied currents 0.2 A and 0.6 A, when the excitation amplitude changes from 5 mm to 25 mm, the resistive force generated by the MR damper increases by 34% (343 N-455 N) and 23% (892 N to 1095 N), respectively.The observation indicates that a higher excitation amplitude would result in a larger hysteresis, and such an influence is more considerable at a lower current.Further, when the excitation amplitude is small, the maximum damper force is found to be more sensitive to the amplitude effect.The pattern of the force-velocity curves in figures 4(a) and (c) shows that the slope of the force-velocity hysteresis in the post-yield region gradually decreases as the excitation amplitude becomes larger.This phenomenon is more visible at the current level of 0.6 A. However, compared to the frequency effect shown in figures 3(a) and (c), the width of the pre-yield region is much less sensitive to the change in the excitation amplitude.

Figure 3 .
Figure 3.Effect of excitation frequency on the hysteresis behavior of a MR damper under amplitude of 15 mm and current level of: (a) and (b) I = 0.4 A; (c) and (d) I = 0.8 A.

Figure 8 (
b) shows a sample of a fitted force-velocity hysteresis curve of the MR damper by MATLAB curve fitting tool, of which the damper is subjected to an excitation amplitude and frequency of 28 mm and 0.5 Hz, respectively, and the supplied current is 0.6 A. The identified parameters are: c = 3.08 N • s mm −1 , k = 1.302N mm −1 , α = 729.3N, β = 0.10, δ = 0.509 and F 0 = 13.9N with an R 2 of 0.9997.The identified six parameters for the tangent hyperbolic model in equations (2) and (3) are portrayed in figure9.It is found that for each parameter when the first-order polynomial

Figure 4 .
Figure 4. Effect of excitation amplitude on the hysteresis behavior of a MR damper under frequency of 1.0 Hz and current level of: (a) and (b) I = 0.2 A; (c) and (d) I = 0.6 A.

Figure 6 .
Figure 6.Physical interpretation of the hysteresis model.

Figure 7 .
Figure 7. Effects of parameters β and δ on basic hysteretic loop.

Figure 8 .
Figure 8. Accuracy of fitted model to experimental data; (a) root-mean-squared error of fitted-curves (b) a sample of fitted-curve (A = 28 mm, f = 0.5 Hz and I = 0.6 A).

Figure 10 .
Figure 10.Architecture of MLP; (a) Schematic of a fully connected MLP feedforward neural network; (b) input layers of MLPs for the six parameters.
and second-ranked features change.

Figure 13 .
Figure 13.Comparison between the SNN model and the algebraic model by Zhang et al [27].

Figure 14 .
Figure 14.Comparison of the average (solid-line) and standard deviation (dashed-line) of R 2 values of the SNN and Zhang et al [27] models under different levels of (a) current levels; (b) frequencies; and (c) amplitudes.

Table 2 .
[27]of polynomial regression lines and Zhang et al[27]'s model for the identified parameters.

Table 4 .
Configuration and accuracy of trained neural networks corresponding to the model parameters.
• The proposed SNN model manifests a superior capability of modeling the hysteresis behavior of MR dampers.It cannot only accurately predict MR damper performance in the preyield and post-yield regions, but more importantly, within the most challenging roll-off region in the vicinity of zero velocity.• Compared to training a neural network model to directly describe the nonlinear hysteresis behavior of a MR damper, employing six feedforward MLP neural networks in the proposed SNN model, with each for identifying one specific model parameter, significantly reduces the computational cost in training.• The SNN model has an excellent generalization in modeling MR damper behavior under different excitation conditions and current levels, of which the predicted damper force exhibits consistent and stable accuracy.