Nonlinear Modeling and Experimental Characterization of Twin-Tube Hydraulic Adjustable Damper: Effects of Cavitation On Hysteretic Behavior

: The hydraulic adjustable damper has attracted wide attention due to its superiorities of low energy consumption, fast response, strong durability, high reliability and simple structure. However, there has been no published detailed analysis about the effects of cavitation on the hysteresis of the hydraulic damper damping output. Furthermore, the existing damper models with simplified assumptions for the cavitation have not been completely studied. Therefore, a nonlinear model of twin-tube hydraulic adjustable dampers (twin-tube HAD) is proposed with an emphasis on the cavitation properties. Polytropic change in the gas content, seal friction, oil viscosity and the gas-oil emulsion flowing through orifices or valves are taken into consideration in the model. The cavitation form of twin-tube HAD valve is studied in depth and the dynamic cavitation number of hydraulic oil is formulated as a function of the gas volume fraction, then the damping force is characterized by the gas volume fraction. The model proposed in this paper can be used for accurately analytical investigation and it is useful in reducing damage from cavitation in similar nonlinear equipment. The mathematical model is validated by comparison against experimental results carried out on HONDA-EG8-RH twin-tube HAD in damper test facilities.


Introduction
The design of the suspension system plays an important role in vehicle handling and comfort. In order to improve the performance of the vehicle, it is necessary to implement the detailed design for each part of the suspension system [1]. The dynamics of the damper and its characteristics exert great effects on the overall performance of the suspension system. The damper should be designed and adjusted to various standards for achieving acceptable comfort and handling [2][3][4]. Therefore, it is significant to understand the damper's force characteristics and its dependence on the damper member' geometric and material properties. According to different ways of generating damping, dampers are divided into pneumatic damper, hydraulic damper and Magnetorheological damper [5]. Hydraulic dampers are widely used in automobiles due to its good performance and low design cost. Hydraulic dampers are divided into mono-tube dampers and twin-tube dampers with their own advantages and disadvantages. Mono-tube hydraulic damper has only one piston acting on the active cylinder, and it is sensitive to the impact force. Twin-tube hydraulic damper uses various valve systems and gas storage chambers to configure the flow of oil, which makes it good adaptability. With the higher requirement for automobile comfort, twin-tube HAD is widely applied in suspension systems. However, cavitation makes twin-tube hydraulic damper inoperative and its heat dissipation ability less effective. It is necessary to conduct an in-depth study on the damping output mechanism of twin-tube HAD. From previous studies, there are still many problems in the working process of twin-tube HAD. In order to further optimize the working performance of twin-tube hydraulic damper, scholars have done massive researches on the causes of damping output hysteresis. Most scholars believe that the damping output hysteresis is caused by oil cavitation phenomenon. However, there has been no published detailed analysis about the effects of cavitation on the hysteresis of hydraulic damper damping output [6][7][8][9][10].
There have been many studies on different aspects of the characteristics and performance of hydraulic dampers [11][12][13][14]. Guan et al [15] investigated the dynamic performance of twin-tube shock absorber by experiment and simulation. The virtual prototyping technology is proposed to further investigate dynamic properties of shock absorber. Badri et al [16] described the fluid flow through the internal orifices between the compression and the rebound chambers by the CFD analysis. Results showed a viscosity increase of 70% when the magnetic field excitation current was elevated from 0 A to 5 A. Lee [17] studied the analytical model of mono-tube hydraulic dampers. The parameters be determined from the test data are the four flowing coefficients. Lee used the test data to verify his mathematical model. Moon and Lee [18] developed an analytical model for displacement sensitive shock absorber by using fluid flow distribution principles. Alireza et al. [19] established a mathematical model for mono-tube hydraulic dampers to study the geometric and material characteristics of gasket stack, bypass valve characteristics, overall damping and hydraulic oil characteristic parameters. The continuity equations of the compression chamber and the rebound chamber are given, as well as the state equation of the gas in the air chamber, the motion equation of the floating piston, and the deflection of the spacer set. By combining the component equations into a set of nonlinear ordinary differential equations, the problem can be solved numerically. Lang [20] was the first to establish a bubble based on evaporation, known as cavitation, to study the damping output hysteresis of the hydraulic damper. However, Segel and Lang concluded that the more likely assumption is that bubbles are not formed by steam but by nitrogen in the reserve chamber. This is known as frothing. The physical factor considered by Lang were the compressibility of oil, leakage around damper piston, deflection of damper cylinder, vaporization of oil, and the flow through valves. He used an analog computer to simulate, which can achieve good consistency between the mathematical model and test results. Duym et al. [21] worked on an analytical model of the hysteretic behavior in twin-tube dampers. They proposed that besides the oil compressibility, the root of the damping output hysteresis was the compressibility of variable gas phase in the shape of bubbles. Yuming et al. [22] introduced a nonlinear stiffness and damping characteristics of hydro-pneumatic suspension, taking into account the multiphase of gas, friction sealing, and gas-oil emulsion in the model of the damper. The model is formulated with the consideration of one and two bleed orifices configurations of the strut. Shuai Yue et al. [23] introduced a nonlinear mathematical model of the liquid spring damper, including entrapped air, the variable liquid bulk modulus, flow inertial effects and cavitation phenomena. A very good match can be established between simulation results of analytical models and experimental data in damping output characteristic. Alexander Lion et al. [24]studied the thermal and dynamic effects of the damper, and the model parameters were identified based on the test and the thermo-mechanical coupled dynamics model that was established to analyze the temperature of the thermal equilibrium state. To sum up, the above literatures investigate the effects of flow distribution characteristics, valve structure [17][18][19] and cavitation characteristics [20][21][22][23] on the performance of hydraulic damper.
The damping force output is the most important performance parameter of twin-tube HAD [25,26]. However, the damping output hysteresis in the damper is the main reason affecting the instability of its performance. Therefore, it is necessary to establish a accurately mathematical model describe the working process of the twin-tube HAD. However, cavitation makes twin-tube hydraulic damper inoperative, the generation of the cavitation in the twin-tube HAD is not only due to the pressure difference between the two sides of the valve port, but also the high frequency reciprocating motion of the damper, which causes the increase of entropy [24,27]. Further, the oil density and the bulk modulus are lowered in a state in which the bubble and the liquid coexist. When the damper is in operation, the damping force cannot output in ideal state. The heat-generating, heat-conducting, and temperature-increasing processes in the damper will greatly affect its damping output [28]. It is known from relevant literatures that researches on the hysteresis mechanism of dampers are mostly considered at the macro level from the instantaneous gas content of the damper and the gas compressibility. However, less attention is paid to the entropy increase process. In order to facilitate the in-depth study of the twin-tube HAD hysteresis mechanism, it is necessary to consider the entropy increase and accurately describe the cavitation phenomenon.
In this paper, twin-tube HAD of HONDA-EG8-RH front suspension is taken as the research object and the mathematical model of the twin-tube HAD will be discussed. The analysis focuses on the influence of the process of cavitation generation on the damping performance. Firstly, a combination of non-parametric method with parametric method is used to simplify complex valve systems. The mathematical model will take into consideration the solid-liquid coupling, seal friction, valve opening, effective piston force area, bubble volume fraction, the influence of temperature change and oil entropy increase on the gas volume fraction and oil viscosity during the operation of the damper. The model gives the continuity equation of the compression valve and the rebound valve, as well as the continuity equation of chambers, the bypass valve opening adjustment equation.

Regulating mechanism for twin-tube HADs
There are mainly two kinds of HADs available for damping force control. Typical HAD employs damping force adjustment rod which can be controlled to lift up and down via rotation angle to form various needle throttle, generating different pressure drops between the upper oil chamber and the lower oil chamber. This method is not only difficult to operate, but also requires big torque, resulting in high processing cost. The other one employs bypass valve to change the flow rate during the compression and rebound strokes. This kind of HAD can be more accurately controlled by simple structure. In this paper, a prototype of twin-tube HAD controlled by bypass valve will be investigated. The physical model of the established twin-tube HAD is shown in Fig.  1.
The piston moves up and down in the working chamber with high frequency when the relative displacement occurs between the sprung mass and the unsprung mass. It will cause the transfer of the hydraulic oil from the high pressure chamber to the low pressure chamber. Damping force is generated due to the friction of the liquid between molecules themselves and between the liquid and the orifice plate. Thus, the parameters of the liquid and the structure of the valve are the key point to the generation of the damping force. The flow valve, adjustable valve and bottom valve constitute the whole valve system of twin-tube HAD. Simplified structure of twin-tube HAD in Fig. 1 mainly contains the compression valve assembly, the rebound valve assembly, the adjustable valve assembly, the bottom valve assembly, the components of piston assembly, the bottom valve assembly with relief valve (active in compression stroke) and refill valve (active in rebound stroke). The generation of damping force has great correlation with the parameters of those valve components.
The working stroke of the twin-tube HAD composes of the compression stroke and rebound stroke. Fig. 2 presents the flowing situation of the hydraulic oil during the compression stroke and rebound stroke. In the compression stroke, as shown in Fig. 2(a), the pressure from the hydraulic oil loads on the spacer of the compression valve, causing the deformation of the spacer of the compression valve and generating a throttle orifice. Meanwhile, the hydraulic oil transfers from the compression chamber to the rebound chamber. When the pressure in the compression chamber reaches the critical pressure, the hydraulic oil in the compression chamber gets through the relief valve to the reservoir chamber. In the rebound stroke, as shown in Fig. 2(a), when the pressure in the rebound chamber is greater than that in the compression chamber, the hydraulic oil flows from the rebound chamber to the compression chamber. The pressure of the hydraulic oil in the reservoir chamber is larger than pressure in the compression chamber so that it will open the refill valve by compressing the adjustment spring. As the spring compresses, the spacer of the refill valve generates vertical displacement and the effective open area of the refill valve increases from zero. Then the hydraulic oil in the reservoir chamber gets through the refill valve to the compression chamber.

Analysis of hysteretic theories and laboratory measurements
The twin-tube HAD comprises of two chambers connected via four one-way throttles and one throttle, as illustrated in Fig. 1. The compression chamber and the rebound chamber are filled with hydraulic oil. In order to guarantees the steady output of the damping force, the storage chamber is shared by both hydraulic oil and nitrogen with a certain initial pressure. In the running process, the flow state of oil in the twin-tube HAD is shown in Fig. 2. However, after the twin-tube HAD has been in operation for some time, the bubble dissolved in the oil will be separated out (At higher temperatures and higher turbulence of the oil, the bubble will precipitate above its vapor saturation pressure). After several reciprocating motion, bubbles are evenly distributed throughout the damper, as shown in Fig. 3(t 1 -t 7 ). A mixture of gas and hydraulic oil is produced and the effective volume modulus of the mixture will be lower than that of the hydraulic oil. During the compression stroke, with the downward motion of the piston, as illustrated in Fig. 3(t 5 -t 7 -t 2 ), the instantaneous pressure of the oil is higher than its vapor saturation pressure and the gas will be dissolved in oil. The volume difference is Δv, as shown in Fig. 3(t 5 →t 6 ). At this time, if high frequent excitation occurs in the twin-tube HAD, then the hysteresis is manifested in twin-tube HAD damping output, which is also called idle. In addition to the leakage through the piston seal, a portion of the mixture in the compression chamber flows through the piston compression valve to the rebound chamber. Since the volume of the compression chamber is greater than that of the rebound chamber due to piston rod, the oil entering the rebound chamber through the compression valve is not enough to fill the volume difference between the compression chamber and the rebound chamber. As a result, a part of the mixture liquid enters into the gas storage chamber through the bottom valve and the adjustable valve. During the rebound stroke, as shown in Fig. 3(t 2 →t 5 ), a part of the gas will be dissolved in oil because pressure exceeds its vapor saturation pressure. The volume difference is Δv, as shown in Fig 3.(t 2 →t 3 ). With the upward motion of the piston, the mixture of the rebound chamber flows into the compression chamber through the rebound valve. However, the mixture is insufficient to fill the compression chamber, and then parts of the mixture flows into the compression chamber through the bottom valve and the adjustable valve. Under diverse excitation frequencies and excitation times, the time of the bubble generation in the oil is different, which mainly depends on pressures of the two chambers, oil temperature and oil disorder degree. The variation of the gas content in the oil also directly affects the damping output characteristics of the twin-tube HAD. In order to verify the correctness of theories, it is necessary to design an experimental scheme and conduct experiments. First, the damper of the HONDA-EG8-RH front suspension will be used to set up test bench to carry out related experiments. The feasibility, necessity and experimental equipment of the experiment are illustrated in the following part.
In order to investigate the stiffness and damping properties of the twin-tube HAD with gas-oil mixture, laboratory experiments were conducted on dynamometer to acquire its response characteristics under pseudo-static and dynamic excitations. The twin-tube HAD was installed vertically between fixed cross-beam via force transducer and lower fixture, as illustrated in Fig. 4. The force transducer was mounted between the fixed cross beam and the strut to acquire the total damping force. A displacement sensor mounted on the bottom measured the displacement of the reciprocating motion of the piston. The piston speed could be calculated by differentiating the piston displacement with respect to time in the LabVIEW program. The data from displacement sensors were acquired by using NI device. At the same time, the PC could deal with these signals via the script in LabVIEW.
Three series of experiments were conducted to study the cavitation in the hydraulic oil, so as to characterize damping output instability properties of the twin-tube HAD. The first series of experiments was conducted at different vibration frequencies and the same amplitude. The initial position of the piston was in the middle of the inner cavity of the damper, and an amplitude of A=50mm was applied. The movement of the piston was given by the sinusoidal motion trajectory, and the excitation frequency of the damper was adjusted with f=0.417Hz, 0.834Hz, 1.669Hz, 2.503Hz and 3.338Hz. The second series of experiments was conducted at different amplitude and the same excitation frequency. The piston was positioned in the middle of the damper cavity and a harmonic excitation was applied to the piston with the amplitude A=10mm, 25 mm, 50mm, and the frequency f= 2.503Hz was applied. During the first and second series of experiments, the signals from the pressure transducers, force transducer and temperature transducer of the thermocouple were continuously acquired in the data acquisition system at a sampling frequency of 5Hz.
In the third series, the pressure-displacement, force-displacement and force-velocity characteristics of the twin-tube damper were measured under different harmonic excitations. The cylinder temperature was continuously monitored and kept within the range of 25±2 ℃, 75±2 ℃ through a large size fan to study the damping output characteristics of dampers at different temperatures. The measurement was mainly for complex road conditions. Then, various performance parameters of the damper were measured at different excitation frequencies of 0.417Hz, 0.834Hz, 1.669Hz, 2.503Hz and at different amplitudes of 10mm, 25mm, and 50 mm. After a series of experiments, its working time for the calculation of the normal working will be used as life of the damper. The wear of the damper during the experiment is unavoidable, and it is necessary to retest the performance of the tested damper to study the change of its performance, including different amplitudes, different excitation frequencies, and different maximum speeds. The above experimental processes ensure that the adjustable valve of the damper is in the same opening state. The third series of experiments was repeated by adjusting the opening of the adjustable valve in order to investigate its effects on the damper stiffness and damping properties.

Nonlinear model of twin-tube HAD under cavitation
The twin-tube HAD has complex working conditions and high nonlinearities. Therefore, this section focuses on the establishment of the nonlinear twin-tube HAD model which introduces the influences of cavitation on oil compressibility, bulk modulus and viscosity. Firstly, this chapter establishes a classic damper model based on force and velocity. The model describes in details the effect of the volume of the oil gas on compressibility, bulk modulus and viscosity. The damper cavitation process is analyzed to study the gas production. Compression stroke and rebound stroke are established by integrating the dynamic processes of the chambers. The compressibility of the fluid, entrapped air, inertial effects of the fluid flow and cavitation phenomena are included in this mathematical model with the following assumptions: 1. The deformation of the damper's outer cylinder is ignored during the large impact force.
2. The entrapped air and vapor are uniformly distributed in the liquid. Pressure, fluid density, effective bulk modulus and viscosity are defined as lumped parameters for each working chamber.
3. The vapor saturation pressure of the gas is not a criterion for the formation and collapse of the gas phase.
4. The pressure, oil temperature, oil viscosity and oil density in the same working chamber are the same. The parameters of the oil at the valve or porosity are neglected.
5. The bubbles generated by cavitation are spherical, and bubble and liquid do not transfer heat.

Continuity equation of chambers and piston motion
The proposed mathematical model of the twin-tube HAD contains the equation of the piston motion, equation of the mass or volume balance for each fluid chamber and equation of the thermodynamic energy balance for each fluid chamber. The cavitation and thermal effects occur mainly in the compression chamber and the rebound chamber. Only two chambers are considered in establishing the correlation equation, twin-tube HAD Model parameters is shown in Fig. 5. The equations of the volume balance for the fluid filling two chambers are as follows, Where V 1 and V 2 are the liquid volume of the rebound chamber and the compression chamber, respectively; T, Re and ρ represent the oil temperature, oil Reynolds number and oil density, respectively; β (p, T, Re) is the compressibility coefficient and α (p, T, Re) is the thermal effect coefficient, which depend on the pressure, temperature and Reynolds number of the fluid. The terms at the right of each equation represent the volumetric inflow rate for each chamber. The equilibrium of the fluid volume can be replaced by the equilibrium of the fluid mass: where m 1 is mass of the fluid in rebound chamber; m 2 is mass of the fluid in compression chamber; and Q m is mass outflow rate from the compression chamber.
The relationship of the general variable fluid volume for each chamber can be denoted as follows: where V 01 and V 02 are the initial volumes of the rebound chamber and the compression chamber respectively; A P,1 and A P,2 are piston areas in rebound chamber and compression chamber, respectively.
The equation of the thermodynamic energy balance shows that the increased enthalpy turns into internal energy changes and work done by fluid. For the rebound chamber, the enthalpy is expressed by fluid parameters inside the chamber: where c p is a constant pressure heat capacity. Conversely, for the compression chamber, the enthalpy is represented by the gas parameters delivered through the valve: Under the assumption that the lossless valve flow and enthalpy are equal, the above equations can be attributed to the overall thermodynamic energy balance of the damper. the work done on the fluid is equal and the amount of thermodynamic energy change is equal.
According to the Bernoulli equation [29], the pressure drop of the piston pore consists of the flow pressure loss ΔP f and inertance pressure loss ΔP i . The relationship between the flow pressure loss model and the flow rate is described as follows: where l is the length of the orifice; S pg is the orifice area; C is the discharge coefficient; μ a is the dynamic viscosity of fluid in the orifice. Due to the discontinuity of the viscosity of the two chambers, μ a is the average viscosity for the two chambers. The pressure drop generated by the fluid flowing from one end of the valve to the other end under a certain pressure is called an inertance pressure drop: where I i is the inertance of fluid and C c is contraction coefficient of the orifice. The liquid flow equation for the two chambers is: The sealing friction between the piston and cylinder is characterized by the viscous slip, Strikbeck effect and hysteresis characteristics during the periodic motion. Therefore, the friction force can be expressed as a function of the piston velocity v. In addition, a linear force transition near zero velocity (|v| < v h ) is introduced to describe the hysteresis behavior of the friction force, as shown in Fig. 6. When the piston speed exceeds the hysteresis transition speed v h (|v|≥v h ), the friction can be expressed as: where F s is stiction force; F c is Coulomb friction force; μ v is viscous friction coefficient; k is Stribeck coefficient, and sgn represents the sign function.
The frictional force f in the linear transition zone (xv<0 & |v|≤v h ) can be expressed as, Where x is the displacement of the piston.
In the same velocity direction, when the piston's velocity is in the transition region, the friction force is obtained from Eq.(10). When xv <0 and |v|≤ v h , the friction force increases linearly with the increase of the piston's velocity.
In the compression stroke, the piston rod moves down and the compression force can be formulated as, Where S a is the annular area of the piston; S pr represents the area of the piston rod.
In the rebound stroke, the piston rod moves up and the rebound force can be formulated as, Values of pressure in chambers 1,2 and 3 will be inserted in Eq. (11,12) as functions of time and the damper force can be determined.

Model of density, effective bulk modulus and viscosity
The presence of the entrapped air bubbles and vapor in liquid can significantly influence the density, effective bulk modulus and viscosity of fluid [7,30]. Therefore, it is necessary to analyze air bubbles and vapor in liquid. In the working process of the damper, there must be three states, including the initial stage, the second stage and the third stage. In the initial stage, only the oil and the negligible gas exist. In the second stage, the oil and air bubbles exist. In the third stage, the liquid, air and vaporized liquid coexist.
From previous studies, When the damper is severely impacted, the local pressure inside the hydraulic damper will fall below to the vapor saturation pressure and the cavitation would occur [31]. As a result, the vapor volume V v is introduced. It is used to detect the cavitation and its evolution over time. For a given moment, the mass balance of the undissolved gas-liquid mixture, the liquid vapor characterized by the ideal gas law in the working volume is used to define the vapor volume: where ρ m (P v ) is the density of the mixture; m is the total mass of chambers; V is the total volume of chambers; M v and P v are the molar mass of the liquid vapor and the vapor pressure, respectively; R is the gas constant; and T is the chamber temperature. By differentiating the equation (V v ) with respect to time, the change of vapor volume with time can be obtained. However, in each working chamber, the change in pressure is instantaneous. The initial definition of the vapor generation is that the pressure of the chamber is lower than the vapor saturation pressure, that is, P i ＜P v . At this time, the vapor volume starts to generate: The above vapor generation is based on the vapor saturation pressure. However, during the operation of the damper, the vapor will also appear when the pressure of the oil is higher than the vapor saturation pressure, which will be elaborated in details in Chapter 3.3.
Apart from the cavitation, another factor which significantly influences the density, effective bulk modulus and viscosity of fluid is air bubbles. It is necessary to describe the bubble volume fraction in the working chamber. It is assumed that the mass of the oil per unit volume is V liq (P a ) = 1, where P a is the atmospheric pressure. The relationship between the volume and mass of the gas in the oil can be expressed as, where X a is air content; V air (P a ) is the volume of air at the atmospheric pressure; m fluid is the fluid mass; ρ liq (P a ) is the density of liquid at P a , and ρ air (P a ) is the density of air at P a . Therefore, the density of the fluid under a certain pressure ρ fluid (P) in the absence of vapor is, where V liq (P) and V air (P) are volumes of liquid and air under a certain pressure P . Based on the above assumptions and the law of conservation of mass, V liq (P) and V air (P) can be obtained as, and where ρ liq (P) and B liq (P) are the density and bulk modulus of liquid under a certain pressure P; and γ is the polytropic exponent. The bulk modulus of liquid is assumed to change linearly along with the pressure [32], where C B is the pressure coefficient. According to the definition of the bulk modulus of fluid, the density of the liquid ρ liq (P) can be obtained with the integration of B liq (P) under the pressure P [33]. Then ρ liq (P) can be expressed as, According to the above formula, the oil density and bulk modulus can be expressed under a certain pressure P when there is only air bubble in the oil. The dynamic viscosity of the oil can be represented by the average viscosity between the chambers. Thus, the density, effective bulk modulus and viscosity of fluid can be expressed as, 1 1 where μ and K μ are the viscosity of liquid and experimental coefficient by Hayward [34]. When the pressure of the fluid changes, cavitation occurs and a volume of vaporized liquid is produced. In this case, it is assumed that the liquid, air and vaporized liquid coexist and they are uniformly distributed in the fluid. Then, V liq (P) can be written as, where T 0 and V v are the ambient temperature and the vapor volume, respectively. Therefore, without considering the effects of oil impurities and oil temperature on the cavitation, mixture of vapors produced by cavitation and gases in the oil, the density, effective bulk modulus and viscosity of mixture can be expressed as,

Dynamic process of cavitation
During the reciprocating work process of the twin-tube HAD, the cavitation is not only affected by the pressure, but the pressure difference, oil flow rate (degree of oil vortex) and gas core size in the dynamic process affect the formation of the cavitation in the damper [35]. The higher the temperature of the hydraulic oil, the greater the tendency of the hydraulic oil cavitation. Therefore, the temperature changing process and the state change of the hydraulic oil should be considered so as to establish a more accurate mathematical model. The occurrence of the cavitation is caused by the small passage of the hydraulic oil through the valve. The difference is, not only the oil which passes through the valve port has a pressure difference, but also the hydraulic oil reciprocates through the small passage on the valve at high speed. The dynamic factor of cavitation formation is the flow of the fluid, and the violent reciprocating motion of the damper is also a necessary condition. When the damper speed is high enough, the cavitation nucleus is dispersed in the oil. During reciprocating motion of the damper, the cavitation nucleus becomes bubble which is dispersed in the oil due to pressure drop, and a portion of the gas-liquid mixture forms emulsion. Since the emulsion is compressible, the damper cannot output the damping force in an ideal state. In this study, the change of oil velocity is considered during the dynamic operation of the damper. The relationship between oil temperature rise and cavitation formation is established. In addition, the damper which has been used for a period contains different cavitation gas cores (such as oil impurities, emulsion), and the conditions of the cavitation under this situation will be studied. Cavitation usually occurs in high-speed flowing liquids. During the operation of the damper, the piston reciprocates at a high speed with the hydraulic oil passing through the valve port, thus generating a complex liquid flow that satisfies the condition of cavitation. In terms of previous researches, most models use the method of the fuzzy cavitation to study the process of the liquid cavitation from the view of macroscopic point. In this study, a detailed damper cavitation model is established based on liquid cavitation conditions (such as the Reynolds number of the liquid, the gas content, the turbulent pulsating pressure). The cavitation number is introduced to quantify the cavitation phenomenon and define the degree of the cavitation. Usually the definition of the cavitation number can be written as, where P ∞ and μ ∞ are the static pressure and velocity of the undisturbed reference section fluid, respectively; ρ and P v are the density of the flowing liquid and the saturated vapor pressure, respectively.
The existence of cavitation core in the liquid flow is the material basis of the cavitation. To study the number of the cavitation in a number of liquid, it is necessary to clarify the effect of the initial bubble size on the cavitation. The molecules with higher kinetic energy can overcome the attraction of molecules on the bubble surface to inflate the bubble. According to equilibrium theory, the bubble equilibrium state (As shown in Fig. 7(a)) is expressed as: where P v is vapor saturation pressure of the gas; P g is internal pressure of the bubble; P is pressure of the liquid; τ is the bubble surface tension coefficient; R is the radius of the bubble; R 0 is the initial radius of the bubble; P g0 is the initial pressure inside the bubble; and P v (T) is a constant in Eq. (25). From Eq. (25), it can also be seen that P v is related to the temperature of the liquid, a curve of P=P(R) can be obtained, as shown in Fig. 7(b), where R c is the critical radius. On the left side of R=R c , it is a stable balance for liquid internal bubbles. On the right side of R=R c , it is an unstable balance that the bubbles inside the liquid will grow up. Cavitation will occur when the external pressure of the bubble drops below the critical pressure of the bubble. Therefore, cavitation is a state of the sudden burst of the local liquid, not the entire liquid.
The liquid flow of the damper, different from the flow of the hydraulic pump, propeller and orifice constant flow, is in a form of the reciprocating motion. Therefore, the formation mechanism of the cavitation is basically the same and the cavitation collapse process is different. The generation of cavitation is a process which is affected by many factors. However, the cavitation collapse process is generally intense and instantaneous. According to the classification of the cavitation, the cavitation of the damper is the vortex cavitation in the shear flow closely related to the turbulent action. The cavitation mechanism is shown in Fig. 8. The formation and development of shear flows are described by mixing processes of fluids with different flow rates in the same direction. The uniform flow velocity at the upper part of the separation stage is μ 1 . The lower uniform flow velocity is μ 2 , and μ 2 <μ 1 . The shear layer δ w (x) is formed at the boundary between the two flows under the action of the fluid viscosity.
This cavitation generally occurs in the wake of blunt bodies, vortex structure of valve and jet of orifices. The damper is a typical vortex structure of valve. The dynamic cavitation number of the damper is defined from the dynamic pressure and fluid velocity of the damper. The volume of gas is described by the cavitation number considering the number of the original cavitation core, and its definition is quantified in chapter3.2. In the process of cavitation collapse, the generated bubble instantly collapses when the bubble is pressed. From previous studies, the shear cavitation refers to the process that the gas core enters into the shear layer and cavitation core is sucked by a coherent vortex structure because of low vortex core pressure and then grows into a cavitation bubble. Its shear flow, the size of the original gas core, the integration of the gas core, the diffusion of dissolved gases and the turbulent pressure pulsation play important roles in shear cavitation.
The damper piston reciprocates to make the liquid reciprocate through the compression valve and the rebound valve. As shown in Fig. 9(a), the liquid passes through the two valves in different directions. Since the valves have the same structure, the same type of cavitation is formed. Here, the cavitation process of the compression valve in the compression stroke is modeled and analyzed, as shown in Fig. 9(b). The evolution of the free shear layer can be roughly divided into three stages, corresponding to the evolution characteristics of three different vortices. In the steady state, as shown in Fig. 10, the gradient distribution of velocity in the valve inlet is shown in section A and the outlet velocity gradient is shown in section B. The velocity profile along the x direction in the shear layer exhibits a curve of continuous derivatives and has an inflection point. The velocity gradient in the y direction is larger than that in the x direction. According to the classic theory of non-viscous flow instability, the perturbation of the reciprocating reflow is amplified and the dominant mode propagates downstream to form a periodic structure, which is the two-dimensional exhibition of Kelvin-Helmholtz vortex. The scale in the y direction is proportional to Δμ=μ 1 -μ 2 , and the vortex moves and develops downstream at the speed of μ=(μ 1 +μ 2 )/2. With increases of Reynolds number, the convection between the downstream vortex and Kelvin-Helmholtz vortex forms a flow vortex and a complex span vortex.
The main feature of this stage is the matching and merging of vortexes. There are two-dimensional vortexes in the flow field, including the Benard-karman vortex and linear vortex for the flow direction. With the original two-dimensional vortex developing into a three-dimensional flow vortex, the strength parameters of the vortex increase. At this time, the intensity of the vortex reaches the maximum, as shown in Fig. 11. After development and growth of vortex, the interaction and dissipation of three-dimensional vortexes cause large-scale vortex structures to break. A turbulent flow consisting of a large number of tiny-scale vortexes is randomly distributed in time and space eventually. The shear cavitation is the same as the fixed vortex cavitation. The cavitation inception occurs mainly in the vortex structure of the flow field. The vortex coherent structure plays a decisive role in the generation and development of the cavitation. There are three important factors affecting the vortex structure during the cavitation process, including the average pressure of the vortex core, the turbulent pressure pulsation in the vortex flow and the gas content in the oil ( including dissolved gas and gas core spectrum). Therefore, the cavitation number of the shear cavitation inception can be qualitatively expressed as, where C pm is the influence coefficient of the minimum pressure of the vortex core; the second item of Eq. 26 is the influence of turbulent pulsating pressure; P' is the turbulent pulsating pressure; the third item of Eq. 26 is the influence of gas content; and P t is the tensile strength of the oil.
The lowest pressure at the center of the vortex is determined by the vortex structure. It has been verified that the linear vortex model established by Rankine et al. can accurately describe the influence coefficient of the minimum pressure. The main control parameters are the velocity of the incoming flow, the strength of the vortex and the radius of the vortex core. The influence coefficient of the minimum pressure can be expressed as,  (27) where P m is the vortex core pressure; P 0 is the vortex structure boundary pressure; V θ is the vortex velocity; Γ is the vortex intensity; and a is the distance from the vortex core. Using the fluid viscosity to diffuse the vortex, the velocity of the vortex is equivalent to, where r is the vortex radius. The gas content contains free gas core and dissolved gas. The fusion of the gas core and the diffusion of the dissolved gas have a significant influence on the growth of the gas core. Therefore, the gas core spectrum (the concentration of gas core, the size distribution of the gas core) plays an important role in the generation of cavitation. The influence of the gas core and the dissolved gas on shear cavitation is mainly reflected in time difference when different gas core size and gas core concentration are captured by vortex and developed into cavity. The qualitative analysis can be performed from the motion growth process of a single bulb in the vortex flow, as shown in Fig.  11(B).
Regardless of the viscosity of the liquid, compressibility, heat conduction and gas diffusion, where v is the velocity vector of the fluid particle; and F is the mass force. Since the bulb is a symmetrical motion, there is only radial motion, and R is the bubble radius; r is the distance from the barycenter to the center of vortex; Therefore, the above formula can be expressed as, where φ is the velocity potential; and f(t) is the integral constant related to time. The boundary conditions of the bubble wall be expressed as, where P r is the static pressure of the outer wall of the bubble; P g and P v are the partial pressures of the gas and vapor in the bubble, respectively; and Τ is the surface tension coefficient of the liquid. The growth equation of the bubble can be obtained by merging the above formula: By combining the growth equation of the bubble with the trajectory equation of the bubble, the growth process of the whole trajectory of the gas core growing into a bubble can be obtained. The motion equation of a bubble can be obtained by Newton's law: where F i is the resultant force of the bubble, which is composed of the following forces: where F 1 is the force by which the oil drives the air bubbles together; F 2 is the resistance caused by the relative movement of bubbles and oil; F 3 is the relative acceleration force of bubbles and oil; F 4 is the mass force of the bubble; u is the speed of the oil; v is the bubble centroid speed; and C D is the resistance coefficient related to the Reynolds number of the bubble. The equation of motion is, By combining Eq. 33, Eq. 34 and Eq. 35, the equation of the bubble motion trajectory can be obtained: By combining Eq. 32 and Eq. 36, the growth process of the gas core along with the trajectory can be obtained.
In order to study the effect of the gas core density on the vortex cavitation, the above-mentioned motion growth process is simplified. Ignoring the influence of dissolved gases, the time T for a gas core with the radius R to reach the center of the vortex with the vorticity intensity Γ from the distance vortex r 0 is: It is assumed that the density of the gas core is n, and it is evenly distributed in the incoming flow. The shear cavitation generally occurs at the tail of the vortex structure, and distance from the starting position to the vortex is c. It can be obtained that the maximum radial distance between the gas core formed by the cavitation and the center of vortex is: It can be seen that the capture time of the gas core is inversely proportional to the square of the gas core density under the same vortex strength and gas core size, which is Γ∝n -2 . For gas core in the axis of the vortex structure, the capture time is inversely proportional to the density of the gas core, which is T=1/nπa 2 u, where a is the radius of the vortex structure. The above equation does not consider the diffusion of dissolved gases. Vortex cavitation forms a gas-bearing cavitation due to the diffusion of dissolved gases. The effect of dissolved gases will result in the cavitation number being expressed as, where -C ps is the number of the primary cavitation in the pure vapor cavitation; K is the empirical correction coefficient; α is the dissolved gas content; and β is the Henry's law constant.
Not only the internal minimum pressure and gas core density of the vortex core, but also the turbulent pulsating pressure plays a decisive role in the generation of cavitation under a certain vortex structure. The more intense the local flow is, the rapider the gas core grows. The pulsating pressure of the turbulent can be expressed as, whereP is the average pressure; and P' is the turbulent pulsating pressure. Therefore, the cavitation number under the turbulent pulsating pressure is considered in the shear cavitation: The above research on the generation process of the shear vortex cavitation generated in the damper is carried out in detail. Main factors affecting shear cavitation include the minimum vortex core pressure, gas core density, dissolved gas content, and turbulent pulsating pressure. In summary, the cavitation number of the shear cavitation for complex flow dampers is: i ps Through the analysis of the above formula and the combination of the cavitation number with the volume of gas produced by unit flow of oil, V air can be obtained by: This section mainly study bubbles which are generated by cavitation and have a great influence on the bulk modulus and density of oil.

Results and discussions
The mathematical model is established by using the corresponding computer code in MATLAB as outlined in Sect. 3.6 and is validated against the experimental test data. A HONDA-EG8-RH twin-tube HAD is used as the damper for tests. The geometric and material parameters used in the model are based on the HONDA-EG8-RH twin-tube HAD. These parameter values are shown in Table 1. The accurately described excitation of the piston motion equation is critical to the correctness of the research model. From the statistical data of the research literature, the vibration mode of the damper on the vehicle suspension is similar to the sinusoidal mode. In this paper, the piston motion equations are assumed to be sinusoidal motion.

Numerical results and experimental data comparison
When conducting the static load test on the dynamometer, different frequencies (like 0.417Hz, 0.834Hz, 1.669Hz, 2.503Hz, and 3.338Hz) are excited with a constant amplitude of the piston displacement (like 50mm). By combining the amplitudes of the damping force with the piston speed under different frequencies, the static characteristics of the twin-tube HAD can be obtained. The experimental results of the force versus displacement (FVD) curve of the HAD can be obtained as shown in Fig. 12.
With the mathematical model of the HAD, ignoring the change of the gas content and oil viscosity, different frequencies (like 0.417Hz, 0.834Hz, 1.669Hz, 2.503Hz, and 3.338Hz) are excited with a constant amplitude of the piston displacement (like 50mm). The numerical results of the force versus displacement (FVD) curve of twin-tube HAD can be obtained as shown in Fig. 13.
Through a comparison of Fig. 12 with Fig. 13, the obvious disagreement between the numerical FVD and experimental FVD data can be found. Firstly, the experimental data is smaller than the simulation data in the compression stroke and the rebound stroke. This error is caused by the viscosity of the oil. The viscosity of the oil become smaller with the oil temperature rising in the experiment, so the output damping is small. Secondly, at the beginning of the compression and rebound strokes, the experimental data showed a damping output idle time which can be improved by considering the cavitation effect. The cavitation effect is quantified by the volume content of gas bubbles. Therefore, the volume fraction of the bubble has been studied to adjust the error between the numerical result and the experimental result.

Parametric study of cavitation
One main objective of this study is to investigate the effect of the cavitation on the characteristics of the damper, and this subsection will discuss this effect under the condition of ignoring the thermal effect. The effects of oil parameters such as gas core density and the dissolved gas content are studied, as well as liquid flow parameters, such as the turbulent pulsating pressure and the minimum pressure of the vortex core. In order to accurately describe the correctness of the model, the stroke and frequency are set as 50mm and 1.669Hz, 2.503Hz, respectively. Other damper parameters are the same as those shown in Table 1. The oscillation of the pressure in the chamber is various at different excitation frequencies, as shown in Fig. 14. The higher the excitation frequency of the piston is, the greater the pressure inside the damper chamber is. The smaller the excitation frequency is, the smaller the minimum pressure of the vortex core is.
Both the gas core and dissolved gas in the oil can form cavity as the cavitation core. The amount of the original gas core in the oil and dissolved gases has an important effect on the cavitation. Fig.  15 shows the average cavitation core number in per unit volume in the treated vacuum and untreated vacuum. Due to the "screening effect", the cavitation core radius of the liquid has upper and lower limits, and the cavitation only generates when the radius of the cavitation core is bigger than 2μm. By comparing two types of the hydraulic oil, the amount of the cavitation core with the treated vacuum is found to be much lower than the hydraulic oil with the untreated vacuum. For the twintube HAD in HONDA-EG8-RH, the vacuumed hydraulic oil is used.
According to the Johnson-Hsieh (J-H) equation and the Rayleigh-Plesset (R-P) equation, the growth of the gas core along with the motion trajectory is obtained. To simplify the expression, some parameters are expressed in a dimensionless time period t/T 0 =1. The evolution process of the gas core at the vibration frequency of 2.503Hz with different initial radii is shown in Fig. 16. It can be observed that R is nearly constant with the increase of t when R 0 is 2 μm, while R starts to increase at t/T 0 =0.05 when R 0 is more than 2 mm. Meanwhile, the maximum value of R is reached at t/T 0 =0.34 and the value of R decreases rapidly to 0 mm at t/T 0 =0.375 approximately. The whole process of the bubble expanding and breaking has a duration of 0.074s. The growth of the bubble lasts for 0.06s approximately, while the collapse process only takes 0.014s. The reason lies in that the growth of the bubble is a process, while the collapse of the bubble is instantaneous.
At different excitation frequencies, the dynamic cavitation number of two different hydraulic oils in one cycle is shown in Fig. 17. In the figure, σ crit is the critical cavitation number when the cavitation has just been generated or just disappeared. When σ>σ crit , the cavitation disappears; when σ≤σ crit , the cavitation occurs. As shown in Fig. 17(a), when σ=0.051<σ≤σ crit =0.524, the vacuum's untreated oil generates cavitation. It is observed from Fig. 17(b) that the initial cavitation number increases and the cavitation area becomes smaller after the vacuum treated. Under different excitation frequencies, a higher excitation frequency leads to a more sensitive vibration of the cavitation number and a larger cavitation interval. According to the difference of the cavitation number in the working process of the damper, the volume fraction of the gas in two chambers in one working cycle is shown in Fig. 18. Fig. 18(a) and Fig. 18(b) show the change of the gas volume fraction in two chambers of oil without vacuum untreated. It can be seen that the greater the vibration frequency of the piston, the larger the peak value of the bubble volume fraction. Fig. 18(c) and Fig.  18(d) show the change of the gas volume fraction in two chambers of the oil with vacuum treated. The peak of gas volume fraction is approximately twice as small as that of oil vacuum untreated.
Through the above analysis, the cavitation model of the oil is taken into account in the working process of the twin-tube HAD. The volume fraction of the bubble is substituted into the volume change of the two chambers. Regardless of the compressibility of the oil, the FVD curve of the damper at the excitation frequency of 2.503 Hz and 1.669 Hz is obtained and compared with the experimental FVD curve by using the same parameters as in Table 1. From Fig. 19 and Fig. 20, it can be observed that the analytical results are much closer to the experimental results at the excitation frequency of 2.503 Hz and 1.669 Hz when the generation and collapse process of the cavity are taken into consideration. Furthermore, the experimental results are larger than the analytical results in the compression stroke and the experimental results are almost similar to the simulation results in the rebound stroke. This illustrates that the FVD curve fitting is insufficient for the damping output hysteresis during compression and rebound processes at both excitation frequencies. The effects of other parameters on the twin-tube HAD damping output are discussed below.

Effects of other damper parameters
Using the mathematical model, the effects of other damper parameters on the performance of the damper can be studied. For example, initial viscosity, initial gas content and initial temperature of the oil can be studied for their effects on the damping output characteristics. Fig. 21 shows the effect of the oil's initial viscosity on damper behavior. As the initial viscosity of the oil increases, the damping output only increases a little, which has little effect on the damping output hysteresis zone. The effect of the oil initial gas content on the damper performance is shown in Fig. 22. As the initial gas content of the oil increases, the damping output slightly decreases while damping output hysteresis zone greatly decrease. When the initial gas content of the oil is 10% or 7%, the damping output hysteresis zone is large. For the vacuumed hydraulic oil with the gas content of 2%, the damping output hysteresis zone becomes small. Fig. 23 shows the effect of the initial temperature of the oil on the damper behavior. The oil temperature slightly changes the slope and hysteresis zone of the FVD curve. Fig. 24 shows the effect of the status of the bypass valve on the damper behavior. The state of the bypass valve has a great influence on the damping force of the damper and has little effect on the size of the hysteresis zone. As can be seen from the results, the initial parameters of the twin-tube HAD have a great influence on the damping output.

Summary and Conclusions
A detailed mathematical model for the twin-tube hydraulic damper is developed and validated with experimental data. Results of experimental data and mathematical model are compared for many different situations. And that the experimental data and model results satisfactorily match for different situations. This study focuses on modeling cavitation and its effects on damper characteristics. The conditions for the formation of the cavitation, the generation and collapse of the cavity were studied. Some other parameters such as the minimum pressure at the vortex center, vortex strength, and cavitation number were also investigated. A twin-tube HAD mathematical model was established by using MATLAB, greatly enhancing the practicability of the model. Some of the damper parameters were studied by using the established mathematical model, such as the initial gas content of the oil, the initial oil viscosity, and the initial working temperature. According to the simulation and experimental results, the following conclusions can be drawn: (1) The twin-tube HAD mathematical model with the consideration of cavitation can accurately predict the damping output performance of the damper.
(2) By reducing the initial gas content of the oil, the damping output hysteresis region will decrease and the damping output will increase to a small extent.
(3) The increase of the oil initial viscosity improves the output of the damping force, but it has little effect on the damping output hysteresis region. The initial working temperature will also affect the slope and hysteresis region of FVD curve in similar fashion. The diameters of piston: D

mm
The total numbers of throttles: n cvt 8 The inner diameter of piston cylinder: D i 30 mm The outer radius of the valve shim: r cv1

mm
The total thickness of valve shims t cv 4.5 mm Length of hydraulic damper cylinder: L 275 mm The initial volume of lower chamber: V 01 2.56×10 5 mm 5 The initial volume of upper chamber: V 02 3.54×10 3 mm 3 Length of the throttle and orifice length: l h 13.5 mm The width of the throttle orifice: l fv 0.6 mm Thickness of the throttle orifice: h fv 0.2 mm The outer radius of the valve shim: r fv1

mm
Radius of the orifice in the relief valve: r fv2 8.