# I. Introduction uality requirements of rubber bushing, such as stability, vibration transmission, and reliability life are becoming more and more demanding in recently decades for design of vehicle suspension systems (1)(2)(3). Rubber elements are the cheapest and most recommended parts when a vibration isolator with certain stiffness is required (4,5). Furthermore, rubber bushing are highly acceptable for reducing noise transmission, accepting misalignment of axes and accommodating oscillatory motions, etc (6). Rubber bushing can be found in vehicles, airplanes, and tractors, such as chassis, suspension, engine mount, exhaust muffler support, vibration isolator, etc. A modern transportation vehicle may have more than 30 pieces of elastomeric bushings integrated in its structure. Due to the nonlinear force-deflection relationship, the credible prediction of the mechanical behavior and math model of the rubber is a great challenge. The rubber bushing models of the suspension and the whole vehicle system are generally simple, which are often presented using a linear spring and a viscous dashpot connection in parallel or in series. Nevertheless, in comparison with other materials, the mechanical properties of rubber are very complicated. The forces-displacement and momentrotation relationships are characterized as nonlinear and time dependent (7), thus, the non-linear hysteresis behaviors of rubber bushing and the frequencydependent loss modulus should be included. Due to its nonlinearities, the cumulative hysteretic energy accumulation inside the elastomeric bushing can cause its mechanical properties to change significantly at various temperatures (8). The changes of carbon black concentration, filler's structure, processing and curing affect the dynamic properties, such as modulus and damping response of rubber (9), this makes the test procedure and instruments used for measuring important in determining the corresponding material parameters (10). Especially, the characteristics of timedependent force-deflection relationship of rubber bushing should take into account preload, excitation amplitude and excitation frequency in order to present the superposition of elastic force, viscoelastic force and friction force (11). In most of the published researches about rubber's behaviors or models of suspension components, the thermo-dynamical effects upon dynamic excitation received few attention, especially the thermal contribution from the hysteresis damping (12).The primary energy loss in rubber bushing upon harmonic excitation is attributed to the viscoelastic-frequency dependent behaviors commonly represented by spring and dashpot elements in classical models. The second energy loss resulted from the hysteresis damping which is attributed to the nearly rate-independent internal material friction between molecular chains and carbon black particles (13). The composition of rubber and properties of fillers determined the ratio of two types of damping and its rough quantification is out of the discussion of this paper (1). The effect of hysteresis damping on the energy dissipation can be fully described using the non-ellipse loop obtained from rubber components when viscoelastic behavior is negligible, usually quasi-static harmonic excitation. Plenty of comprehensive researches have been done to identify the behavior of friction-resulted hysteresis damping in the rubber components. Hysteresis damping is slightly mentioned in early classic models characterizing the damping effect of viscoelastic materials since just the frequencydependent dynamic excitation has been investigated. Another widely accepted damping in structural analysis is Rayleigh damping coefficient, which expresses the damping as proportional to mass and stiffness of component, that expression has experimentally validated especially under higher frequency condition (5). Different viscoelastic models are mentioned trying to characterize the dynamic behavior of rubber components and balance the burden of experiment measurement as well as computational complexities. In sum, a good model should be concise enough to match the complicated mechanical properties without bringing in redundant parameters (14). The main purpose of this paper is to predict the structure response and temperature distribution of rubber bushing through combining the heat accumulation caused by viscoelastic and hysteresis damping. Considering the balance of computation complexity of adding more elements, practical applicability and reliability, Berg's model is employed in this paper. The principle of proper selection of elements to represent mechanical and physical properties of rubber without bringing in excessive calculation is discussed. The proposed hysteresis damping upon static loading is modeled with three springs in series connection and then in parallel to another three in series connected springs. Parameter identification is realized by employing the multivariable nonlinear optimization analysis to guarantee the minimum errors of predicted curves against experiment measurements. Quasi-static and dynamic harmonic excitations with ramp frequencies and amplitudes are applied on the rubber components to represent amplitude dependency of stiffness and damping. With the identified parameters, the dynamic behavior of rubber bushing under harmonic excitation is simulated in COMSOL Multi-physics, which tallies with the expected distribution. Furthermore, the estimation of temperature distribution is presented with great satisfaction. Many models are propose dtrying to accurately reproduce viscoelastic materials' mechanical behaviors using lumped or discrete springs and dashpots (15,16). # II. Principle of Model Construction and Constitutive Relationship In the commercial software like Abaqus, the generalized Maxwell model using in-series connection of spring and dashpot is recommended to properly represent frequency-dependent dynamic behaviors of rubber. Berg's three branches model (Fig. 1) incorporating the comprehensive effect of elasticity, viscosity and especially friction effect is a great progress characterizing the amplitude-related hysteresis damping of rubber component upon static loading conditions. The total force of applied uniaxial load is shared by the three branches while each branch displays the individual strain. Spring elements k 2 is selected to represent the elastic behavior of rubber busing, especially the modulus of viscoelastic material after infinite relaxation time (branch 1). Spring elements k 1 and dashpot in series forms Maxwell model (branch 2) which describes the attributes of the amplitude and frequency dependency of viscous damping. Maxwell module can properly predict frequency-dependent stiffness and damping properties of rubber components but not applicable for higher frequencies excitation and constant stress conditions. The third branch is a friction element which describes static friction (hysteresis damping) F f caused by the relative movement between carbon black and molecular chain. This construction leads to a nonlinear reproduction of the parameter dependency of hysteresis damping and dynamic stiffness in a large range of frequency. a) Quasi-Static Loading(Hysteresis Damping) In Berg's model, a pure algebraic relationship between displacement and friction force is set after configuring the turning point of hysteresis branches (17). The stiffness increases from the lower turning point to the upper turning point representing the auxiliary stiffness induced by the commonly hysteresis harshness. Although the algebraic approach describes the amplitude dependency of stiffness well, it's short of comprehensive coverage of the stiffness variation at every stage of loading and difficult to formulate the random excitation, the latter case is also commonly used in reality. To avoid this drawback of missing particular stiffness evolution, multiple spring elements with special stiffness value for each stage shown in Fig. 2 are suggested to present the typical static characteristic, hysteresis loop. Similar to the bilinear model describing the mechanical properties of elasticplastic materials, the nonlinear stiffness in the static force deflection stage can be simulated step by step using tri-linear model. ) ? (? b ) ? (? b ) ? (? b ) ? (? b 2 ? K 2 ) ? (? K 2 ) ? (? K 2 ) ? (? K 2 ) ? (? K 2 ? K E 5 4 5 4 3 4 2 3 3 1 2 2 2 5 6 2 5 2 4 5 2 4 2 3 4 2 2 2 2 3 2 1 2 2 2 2 1 1 f ? ? ? ? ? + ? + ? ? ? ? ? ? + ? + = (1) The typical S-shaped force-deflection curve also presents larger elastic elongation accompanying the stiffness reduction again, which is named as Mullins effect, it explains the micro-level damage mechanisms of rubber subjected to large strain (18). The rateindependent response upon static excitation is the crucial phenomenon in modeling internal material friction force and exploring the corresponding hysteresis damping. Under static or quasi-static conditions, time independent and amplitude dependent hysteresis loops are reported in plenty of literatures, even a change of magnitude of deformation rate couldn't lead to significant variation of the hysteresis loops. This interesting conclusion indicated that the widely accepted model of viscoelastic material is insufficient to characterize the quasi-static force-deflection response of rubber bushing. Therefore, the proposed tri-linear model of springs seems to coherently represent the response of rubber component due to the fact that the stiffness relaxation, hardening and relaxation occur with the ramping of strain. Furthermore, the tri-linear model can be adjusted according to the turning point of stiffness if more complicated hysteresis loop is produced in reality. The mechanism of stiffness hardening and relaxation received great attention and the commonly accepted statement can be explained from the relative motion between molecular chains and carbon black particles. It is well known that tiny carbon black particles are frequently added into polymer to enhance materials' damping, stiffness and abrasion resistance. Agglomerates of carbon black particles are linked together by those molecular chains of polymer of different lengths. Irreversible slipping process between polymer molecule chains and fillers can happen when the critical deformation of chain is reached. From the force-deflection test, the reduction of stiffness happened over a large range of strain, which depends on characters of different chain lengths. For molecular chain with different lengths, slipping process starts at different critical stretch. Furthermore, with the additional increase of strain, the relative slipping among carbon black agglomerates is stimulated when the elastic limit of those filler particles is exceeded (14,19). Another logical explanation focuses on recoverable straightening and rupture of molecule chains subjected to recycle loading. Molecular chains cluster are partially trenched straight while remaining molecular chains are still in relaxation because of the different lengths and different elastic limits of the molecular chains. It is assumed that the movement of carbon black particle is instantaneous and molecular chains are more delayed. During initial loading or under small strain, molecular chains cluster are mainly slack without any contribution of rubber hardness. When carbon black particle starts move, the stiffness of rubber increases obviously since more and more molecular chains stretched and starting to undertake load. Further increase of loading or lager strain push the carbon black particle move more, as a response, molecular chains cluster will start to break In Berg's model, maximum force and turning point of amplitude is needed to derive the half of the maximum force which is useful parameter to describe the frequency-independent hysteresis stiffness. Similarly, Fig 2 also displays four turning points which split the loading and unloading curves into three pieces. The auxiliary stiffness of hysteresis damping upon static friction effect are represented with spring K 2 and K 5 which are located between turning points and much more hardened than the remaining springs. Due to the different stiffness of each piece of curve, a tri-linear model with multiply springs K i is suggested, where i=1, 2, 3, 4, 5, and 6. Stresses of different pieces of hysteresis loop should be written as: 3 3 b K + = ? ? . In order to balance the model complexity and computation efficiency, the six springs stiffness can be simplified by setting K 2 =K 5 ; K 1 =K 4 and K 3 =K 6 according to the acknowledged profile of hysteresis curve. The total energy lost (E f ) per cycle caused by hysteretic damping can be produced by calculation of the loop area, once the elastic limit is reached. Due to different length of molecular chains, the hardening and reduction of stiffness happens gradually over a certain range of strain, which brings certain difficulty to detect the strain of turning points. The straightening and breakage of molecular chains are partially reversible during unloading. As a result, the non-overlapping forcedeflection curves upon loading and unloading caused noticeable static hysteresis damping. ? ? ? = = (2) Where ?, ? are total stress, ? 1 , ? 1 are stress and strain of Maxwell elements in series with mass, and ? 2 , ? 2 are stress and strain of spring mass elements. The constitutive equation of the standard liner solid is expressed as Equation (3), dt d? ? ? ? ?k ) k (k dt d? 2 2 1 + = + + (3) The road spectrum of vehicle vibration are collected in time domain and usually transferred to frequency domain before application and analysis. To facilitate the further application of current model, the harmonic excitations are assumed in the following discussion, which leads to, ?? ik ? ? 1 ? ? k k ' ik' k' k + + + + = + = (5) Where the The experimental results of rubber components under static load and harmonic excitation are obtained to identify parameters of each element in the advanced rubber model. All the data measurement are finished using equipment dynamic mechanical analysis (DMA) which measures the influence of frequencies, time and temperature on materials mechanical properties. Harmonic excitation operated using this instrument can a) Friction Force (Hysteresis Damping) In the case of static loads, loading with lower velocity and excitation frequencies are utilized to simulate the quasi-static loading condition due to impractical static loading. Considering the accuracy and reliability, the excitation frequencies of quasi-static condition is chosen from F=0.01 to F=0.1HZ. The very low frequencies are required to keep the contribution of viscous damping in a negligible position. The stress is ramped from 0.001MPa to 0.05MPa to study the amplitude dependency when frequencies are fixed. The corresponding reports about stiffness, damping, tan delta and amplitude will output automatically as long as the pre-set working condition is reached and being stable. Fig. 4 shows the variation of amplitude of the rubber component is subjected to quasi-static low frequency excitation. The observed stiffness dependency upon excitation amplitudes in Fig. 3is the well known Payne effect commonly appearing in fillers strengthened rubber compounds. The reduction of stiffness of rubber components with the increasing excitation amplitudes is attributed to the weak van der Waals bonds existing in the agglomerates formed from carbon black particle (20). The recoverable rupture of physical bonds accumulates when imposing higher amplitudes is imposed, resulting in a significant relaxation of stiffness in the macro-level. For the very small amplitude, friction plays dominant contribution to the raise of stiffness. As the amplitude tends to become larger, part of molecular chains break in the micro-level which leads to friction release in macro-level. Considering the quasi-static hysteresis loop and the amplitude dependency of static stiffness, the tri-linear model has been suggested and its coefficients of springs can be estimated using multivariable constrained method. The equation formulated for parameterizing the hysteresis loop is based on the ratio of energy dissipation over total storage energy, Where S' is the area in side of hysteresis loop, and S is the storage energy during loading. As to the reported tan delta of every loading condition, corresponding strain can be collected and the area S' and S can be easily estimated using the Equation (1). There are three parameters to identify, which means at least two groups of data under the same or close excitation frequencies required to implement the optimization. Stiffness before and after the hysteresis harness during quasi-static measurement are set as the initial value for K 1 and K 2 . The optimization results are illustrated in Fig. 3b. Note that K 1 , K 2 and K 3 slightly decrease with increasing amplitudes and K 3 is much lower than K 1 , which also agrees well with the previous description of molecular chains. For further calculation, the K 1 , K 2 and K 3 are taken as 3.07MPa, 6.62MPa and 0.02 MPa separately. The tan? is taken as 0.05. It is also reasonable to observe that the fracture of molecular chains causes even more reduction of stiffness than the chains' relaxation. # b) Harmonic Excitation (Viscous Damping) For the dynamic stiffness besides of the hysteresis damping, the module is a combination of advanced Maxwell system and spring, which means three parameters k 1 , k 2 and c 1 are required to determine the frequency response of dynamic behaviors, see Equation (13). The dimension of the rubber sample used for the DMA test is less than 5mm in thickness and about 10mm in diameter, furthermore, stress higher than 0.1MPa could lead to shift of sample and fail to produce valid data. Thus, the corresponding output amplitude is pretty low in consideration of the sample size and loading restriction. In the modeling of railway and vehicle track interaction, vehicle suspensions are designed to make sure of the effective isolation of vibration when the frequency is below 10HZ. Track dynamics gradually play the dominated role than vehicle dynamics when the frequencies reached 20HZ. Under this condition, low frequencies force-deflection relations should be discussed here for stability and passenger comfort(3). In the current parameterization, frequencies of harmonic excitation are ramped from 1HZ to 100HZ to balance the time consuming data collection and representativeness over the large range of frequencies. In order to investigate the amplitude dependency of viscous damping, the stresses of each frequency is run from 0.01MPa to 0.1MPa. For low frequencies such as 0.001HZ and high frequency such as100HZ, the measurements are impeditive and some data are missing due to the testing condition of stresses and frequencies are not applicable to this sample. The overall variation of dynamic stiffness and damping over amplitudes and frequencies are plotted in Fig. 4. With the increase of frequencies, amplitude of the same stress will decrease while dynamic stiffness keeps decreasing. As a comparison, the damping presented a peak when the excitation frequency F=30HZ in spite of the variation of amplitude. For more in-depth discussion, outputs are sketched for more direct demonstration. It can be concluded from Fig. 5 that both the storage modulus and dynamic stiffness increase almost linearly with frequency and strain, which clearly illustrated the amplitude and frequencies dependency of dynamic behaviors. The almost linear increase of storage modulus also is proportional in line with energy input As discussed in Fig. 4, the energy dissipation presented with tan delta does not monotonously varying with the frequencies or amplitude. Similarly, note that the incongruous increase of loss modulus, damping and tan delta at frequencies F=30HZ in Fig. 6, that abrupt increase becomes even more obvious at a slightly higher amplitude. In sum, loss modulus increases with amplitude and damping variation displays the same tendency. In the range of lower frequencies, the increase with strain is slight and the increase of frequencies only slightly affects the position of those lines in Fig. 6. It is noticeable that the tan delta data at lower excitation frequencies can be almost fitted using straight horizontal lines. Similar phenomenon is also observed in the quasi-static deformation when harmonic excitation is as low as 0.01HZ. The linear fitting curve for tan delta upon quasi-static loading gives slope-2.12363×10 -4 , which illustrated that ratio of hysteresis damping to the energy input is independent of the excitation frequencies at quasi-static or lower frequency condition. It must be mentioned that quasi-static stiffness and tan delta should be excluded to identify the parameters in Equation ( 13) describing the dynamic behaviors. As to damping, it increases with frequencies gradually at the very beginning and abruptly reaches maximum at frequency around 30HZ. After this peak, damping decreases with continually rising frequencies. This dependency become seven more obvious with the increase of amplitude. The dependency of damping on frequencies can be split into three stages according to the mechanism in the macro-level. All molecular chains of rubber response to harmonic excitation of pretty low frequency changes almost simultaneously with the load and barely lag of phase. This stage is commonly named as high elastic rubber and display limited energy dissipation, which manifests as lower damping and tan delta. While, molecular chains can totally lag behind the deformation of rubber components when harmonic excitation frequencies are within a high range, such as more than 50HZ in the current case. Also, the rubber components upon high frequencies excitation are close to glass state and display much less energy dissipation. A more complicated dynamic behavior is the response of rubber upon the medium range frequencies excitation, for example 30HZ in Fig. 6. In that case, part of molecular chains can follow the macro-deflection of rubber, at the same time, resulting in a larger lag of angle delta and higher energy dissipation. In order to reproduce the variation of stiffness and damping response with the frequencies, there are three parameters k 1 , k 2 and c 1 needed to be identified with the experimental measurement subjected to medium and high-frequency harmonic excitation. From previous analysis, an optimization function should be formulated to satisfy the frequency-dependent stiffness and tan delta together. The multivariable equation is presented in Equation ( 7), ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? + + tan? tan? k' ' k' ? ? ? ? - ' k' k' min 0 0 0 0 2 2 (7) Dynamic stiffness is defined as 2 / 1 2 2 2 dyn ] tan? [1 k' ' k' k' k + = + = . Tan delta is less than 0.2 in almost all testing conditions, therefore, dyn k is approximate to k' with error less than 2%. Furthermore, k' is close to 2 k as long as the excitation frequencies are pretty low, finally, dynamic stiffness dyn k is selected to place the initial value of 2 k . Considering the friction effect of rubber upon the quasi-static deformation and its contribution to the hysteresis stiffness, the data used for the dynamic behavior's parameter identification should reduce the quasi-static stiffness and tan delta. The decrease of hysteresis stiffness becomes weak if the amplitude of dynamic excitation is set as pretty large since the reduction from the rapture of physical of bond between carbon black particles is only confined within a certain range of strain. The static stiffness variation upon amplitude is usually described exponentially or hyperbolically decreasing model. It can be concluded that the static stiffness with the auxiliary of friction effect can be treated as constant as the amplitude at higher frequencies is in hundreds. The initial value of k 1 usually starts from 0.7 dyn k to accelerate the optimization calculation. Damping coefficient is relative low and 1000N/m should be a good initial point. 7 shows the comparison of experiment curve against predicted curve plotted with the parameters obtained from optimization. The stress for this group of data is set as 0.05MPa. It is commonly to have the much better fitting for dynamic modulus than the tan delta fitting. Approximate correction can minimize the relative error of tan delta fitting, while that barely changes the parameter identification results. The fitting results in Fig. 7 give k 1 = 0.49MPa; k 2 = 1.009MPa and c=0.028MPa.t. The obtained stiffness k 2 is pretty close to the dynamic stiffness when excitation frequency f=1HZ, that result is verified against the proper assumption of initial values for optimization procedure. It can be seen easily that dynamic modulus increases with frequencies and there is a peak of damping with the increasing of frequency. As discussed from previous amplitude dependency, the dynamic stiffness and damping is generally amplitude dependent. Actually, DMA tests over a stress range 0.01MPa to 0.1MPa are carried out to found the best data for parameter identification, which is realized with the error analysis. To verify the optimal stress and frequency excitation condition to collect data for accurate model parameter identification, the three-dimensional diagram of dynamic stiffness and tan delta error varying with increasing amplitude and frequency are plotted in Fig. 8.As to the medium frequencies excitation, dynamic modulus and tan delta almost reach the maximum value, then, relative error of dynamics stiffness is minimized due to the absolute increment of stiffness. A better selection of excitation amplitude is located in the central part of the amplitude spectrum. The error analysis results indicated that optimization with experiment data upon lower stresses is poorer compared with that from medium or higher stress excitation. This may be attributed to the error of hysteresis stiffness estimated for quasi-static condition, which highly affects dynamic output upon lower frequencies excitation. While the tan delta of quasi-static loading is barely varied with amplitudes, thus, the Tan delta representing dynamic loss shows no obvious variation of relative error changing with amplitudes. As illustrated, the determinacy of frequency becomes more obvious in determining the error, which suggests the medium frequency as better experiment data range. Since the error distribution characteristics of dynamic stiffness and tan delta are consistent in the optimal selection of excitation stresses and frequencies for parameter identification, harmonic excitation are carried out with stresses set as 0.05MPa. Higher stress and higher frequency excitation condition is kind of out of the service scope of DMA test equipment, Again, excitation of too lower or too higher frequency leads to higher relative error, that makes the medium frequencies as much more cautious selection to collect the reliable experiment data. # IV. The fea Analysis about Temperature Distribution The FEA simulation with multi-physical COMSOL will be used as an instructive prediction about the heat generation magnitude and temperature distribution. The simulation includes two steps, the first is the force-frequency analysis to estimate the heat generation rate and the next step is the temperature distribution at steady state. The mechanical properties of viscoelastic materials are expressed in Prony series since the commonly recommended generalized Maxwell model introducing many parameters. The general Prony series expression representing shear stress relaxation modulus is Figure 9 : The effect of total strain and strain rate R (mm/s) on the stiffness of rubber It is important to claim that the heat generation rate comes from two aspects, one is the viscoelastic damping, and the other one is the hysteresis damping. From the specification of viscoelastic material, the heat generation results of force-frequency analysis merely presenting the viscoelastic damping. Thus, the hysteresis damping describing the friction effect can be set as heat source and the amplitude is determined from structure analysis. To get more precise prediction about temperature distribution, Hot Disk thermal analysis instrument based on the transient plane source (TPS) method is applied to measure the thermal diffusivity and thermal conductivity (Table 1).The probe of the instrument can provide heat source and record the temperature variation with time. The density and heat capacity are estimated using the mixing rule with the weight percentage of polymer and carbon black, which is estimated with thermal gravity analysis (TGA). The loading frequency is set as9.0134 HZ according to the time history of the loading from the field test. Fig. 10 shows the structure and thermal coupling analysis results produced by COMSOL software. The energy loss coming from the viscoelastic damping is presented as the power dissipation density (W/m 3 ) which plays as the major heat source for temperature rising. The secondary heat source is hysteresis damping which is defined as constant heat flow since the displacement amplitude is determined under specified load. The results display that the two faces of the rubber bushing accumulated maximum power dissipation density, in turn, the temperature of that two areas are higher. Gent (21) found the center of the bonded surface displays maximum compressive stress which is almost twice the average and again even much higher than the edges. The current FEA results display the same tendency of stress and temperature distribution which are decreasing from the central zone to the edge of the rubber block. # V. Conclusion A standard linear model parallel to a friction element is employed in this paper to represent the dynamic and quasi-static mechanical properties of rubber with high accuracy and credibility. Spring elements with different stiffness at corresponding stages of deformation are presented to frame the hysteresis damping subjected to the static loading. The tri-linear model is more specific in emphasizing the stiffness hardness and relaxation for each piece of the loop. DMA of different frequency and stress are measured and the model is parameterized using multivariable constrained optimization where data from medium stresses and frequencies are suggested to obtain the optimal dynamic parameters to controls the error. Observed from experimental measurements, the amplitude dependency of hysteresis stiffness illustrates exponential decrease while the dynamic stiffness of rubber increases with the amplitude. Parameter identification procedures are computationally inexpensive and easily applicable, along with the error erroneous analysis. The most adaptable experimental conditions are suggested to demonstrate the model that most precisely matches the real rubber mechanical performance. The powder dissipation density is obtained and its effect on temperature distribution has been analyzed. In the future, the thermal effect on life expectation of rubber bushing will be investigated. 1![Figure 1 : Berg's three branches model](image-2.png "Figure 1 :") 2![Figure 2 : Typical hysteresis curve of rubber under quasi-static harmonic excitation](image-3.png "Figure 2 :") ![from 0 to? 1 ; Stage 2for strain ranging from ? 1 to? 2 :](image-4.png "") ![and stage 3 for strain ranging from ? 2 to? 3:](image-5.png "") ![b) Dynamic Loading (Viscoelastic Damping) Besides static friction, the remaining branches form the standard linear solid model. The first and second branch of current module are placed together to formulate the kinetic equation characterizing the dynamic viscoelastic behaviors. The stress and strain of the two branches are represented in Equation (2),](image-6.png "") ![the Equation (4) into Equation (3) to get the complex modulus depending on frequency,](image-7.png "") 1![/k is relaxation time, representing the ratio of dashpot and spring stiffness in Maxwell model. In the case of viscous damping, the dynamic modulus of standard linear model is given as: be selected over a wide frequency range (0.001 -200HZ), but the frequency lower than 0.01HZ cannot be loaded in this experimental measurement. Parameterization is complicated and time consuming, furthermore, it is usually separated into two steps, static loads and high frequency harmonic excitation.](image-8.png "r ? = 1 C") 34![Figure 3 : Amplitude dependency of quasi-static (dynamic) excitation (a) Modulus; (b) modulus of K 1 , K 2 and K 3](image-9.png "Figure 3 : 4 Viscoelastic") 4![Figure 4 : Response of harmonic excitation versus amplitudes and frequencies: (a)Stiffness (N/m);(b)Tan delta](image-10.png "Figure 4 :") 5![Figure 5 : Illustration of the storage modulus (a, b) and stiffness(c, d) variation on frequencies and deformation amplitudes](image-11.png "Figure 5 :") 66![Figure 6 : Illustrationofamplitude and frequencies dependent: (a) loss modulus (MPa), (b)damping (N/m), and (c) tan delta](image-12.png "Figure 6 : 6 Viscoelastic") 7![Figure 7 : Experiment measurements against predicted behaviors: (a) Dynamic modulus (N/m 2 ) (b) Tan delta Fig.7shows the comparison of experiment curve against predicted curve plotted with the parameters obtained from optimization. The stress for this group of data is set as 0.05MPa. It is commonly to have the much better fitting for dynamic modulus than the tan delta fitting. Approximate correction can minimize the relative error of tan delta fitting, while that barely changes the parameter identification results. The](image-13.png "Figure 7 :") 8![Figure 8 : Erroneous analysis of dynamic stiffness and Tan delta upon with the variation of amplitudes and frequencies](image-14.png "Figure 8 :") 10![Figure 10 : The FEA result from COMSOL simulation (a) mesh, (b) displacement amplitude in x component(m), (c) Total power dissipation density (W/m3), (d) the temperature of rubber bushing at time=20s](image-15.png "Figure 10 :") ![](image-16.png "") 1DensityHeat CapacityThermal diffusivityThermal conductivityHeat transfer coefficient952.541611.440.2020.34310(Kg/m3)(J/Kg.K)(mm 2 /s)(W/m.K)(W=(m 2 /K)) © 2014 Global Journals Inc. (US) Viscoelastic Parameter Identification based Structure-Thermal Analysis of Rubber Bushing Year 2014 Viscoelastic Parameter Identification based Structure-Thermal Analysis of Rubber Bushing Year 2014 Viscoelastic Parameter Identification based Structure-Thermal Analysis of Rubber Bushing Where G i ? is the relaxation time for per series components, ) (? G is the long term shear modulus when t is approximated to infinite and . In the current standard linear model, N=1, and * Development of a rubber component model suitable for being implemented in railway dynamic simulation programs AAlonso NGil-Negrete JNieto JGiménez Journal of Sound and Vibration 2013 * Viscoelastic models for rubber mounts: influence on the dynamic behaviour of an elastomeric isolated system GarciaTarrago MGil-Negrete NVinolas J International Journal of Vehicle Design 49 4 2009 * Modelling of railway track and vehicle/track interaction at high frequencies. Vehicle system dynamics KKnothe SGrassie 1993 22 * A simplified methodology to predict the dynamic stiffness of carbon-black filled rubber isolators using a finite element code NGil-Negrete JVinolas LKari Journal of Sound and Vibration 296 4 2006 * Dynamic simulation studies and experiments on rubber structures used in rail vehicles RKLuo XWu WJMortel Proceedings of the Institution of Mechanical Engineers 227 1 2013 Part F: Journal of Rail and Rapid Transit * Elastomer bushing response: experiments and finite element modeling. Acta mechanica JKadlowec AWineman GHulbert 2003 163 * Constitutive modeling of the large strain time-dependent behavior of elastomers JBergström MBoyce Journal of the Mechanics and Physics of Solids 46 5 1998 * Determination of non-linear dynamic properties of carbon-filled rubbers GDean JDuncan AJohnson Polymer testing 4 2 1984 * Effect of carbon black on dynamic properties of rubber vulcanizates. Rubber Chemistry and Technology AMedalia 1978 51 * Dynamic mechanical behavior of filled rubber at small strains CRoland Journal of rheology 34 25 1990 * Modelling of rubber mounts and applications for time response analysis of dynamic systems including elastomerics X-YPan G-ZChai International Journal of Vehicle Design 49 4 2009 * A new model in railvehicles dynamics considering nonlinear suspension components behavior HSayyaadi NShokouhi International Journal of Mechanical Sciences 51 3 2009 * Constitutive approach to rateindependent properties of filled elastomers MKaliske HRothert International Journal of Solids and Structures 35 17 1998 * A non-linear rubber spring model for rail vehicle dynamics analysis. Vehicle system dynamics MBerg 1998 30 * Shakedown with nonlinear strain-hardening including structural computation using finite element method EStein GZhang JAKönig International Journal of Plasticity 8 1 1992 * Modeling and computation of shakedown problems for nonlinear hardening materials EStein GZhang YHuang Computer methods in applied mechanics and engineering 103 1 1993 * A model for rubber springs in the dynamic analysis of rail vehicles MBerg Proceedings of the Institution of Mechanical Engineers 211 2 1997 Part F: Journal of Rail and Rapid Transit * Softening of rubber by deformation. Rubber Chemistry and Technology LMullins 1969 42 * The effects of surface chemical interactions on the properties of filler-reinforced rubbers. Rubber Chemistry and Technology EDannenberg 1975 48 * Low strain dynamic properties of filled rubbers. Rubber chemistry and technology APayne RWhittaker 1971 44 * Interfacial stresses for bonded rubber blocks in compression and shear AGent RHenry MRoxbury Journal of Applied Mechanics 41 4 1974