A New Formula of Impact Stiffness in Linear Viscoelastic Model for Pounding Simulation

The phenomenon of earthquake-induced structural pounding was extensively studied by some researchers using different models for the impact force. The aim of this paper is to provide a new formula of impact stiffness in the linear viscoelastic contact model, based on the assumption that the maximum impact deformation from the distributed mass model should be equal to that from the equivalent lumped mass model. The correctness and accuracy of the proposed formula have been confirmed by comparing the pounding simulation using the present formula of impact stiffness with those using the existing formulae.


Introduction
During earthquakes, the adjacent structures with different dynamic characteristics can vibrate out-of-phase, which could result in pounding if the reserved separation space is insufficient to accommodate the relative displacement between them.This pounding can generate high magnitude and short duration stress pulse, which could lead to serious damage or even the collapse of buildings and bridges [1][2][3].Pounding damage was identified in several severe earthquakes, for example, the San Fernando earthquake in 1971 [4], the Loma Prieta earthquake in 1989 [5], the Northridge earthquake in 1994 [6], the Kobe earthquake in 1995 [7], and the Taiwan Chi-Chi earthquake in 1999 [8].Therefore, it is valuable to study the collision of structures under earthquakes.Several researchers studied the pounding of adjacent R/C structures [9,10], the pounding of adjacent buildings considering soil-structure interaction [11], and the pounding of bridge deck under bidirectional earthquake excitations [12].Licari et al. [13] proposed a special "multilink viscoelastic" finite element model for seismic pounding between R/C frame buildings, which could simulate the nonlinear time-dependent damping coefficients.
Structural seismic pounding is a high nonlinear phenomenon, which involves plastic deformation, local cracking or crushing, and fracturing as well as friction damage at the contact area.Owing to such complexity, performing the exact mathematical analysis is very difficult for this type of problem.Therefore, based on idealizations and assumptions, several simplified contact models have been developed to study the structural pounding during earthquakes, such as the linear spring model [14], the Hertz model [15], the linear viscoelastic model [16], the Hertz-damp model [17], and the nonlinear viscoelastic model [18].The disadvantage of the linear spring model and the Hertz model is that they are both of full elasticity, so they cannot permit energy dissipation during contact occurrence.It has been shown that the other three impact models, the linear viscoelastic model, the Hertz-damp model, and the nonlinear viscoelastic model, can provide sufficiently accurate results for the overall structural response if the parameters in these impact models could be effectively selected [17].
The linear viscoelastic model has been widely used in the simulation of structural pounding because it is not only simple but also easy to implement in commercial software programs, even if it does have some drawbacks: an initial jump of the pounding force presents at the beginning of the impact and the tensile force could occur before the separation of two pounding components [17].In order to avoid the tensile impact forces arising during detachment, a variation on the linear viscoelastic model was proposed by Komodromos et al. [19].However, it still suffers a major limitation: the proper value of the impact stiffness to be used in the simulation cannot be known in advance.An impact spring with understiffness could lead to the overlap of adjacent structures while an impact spring with overstiffness could result in unrealistically high impact force and unstable numerical simulation.
Existing studies have demonstrated that impact stiffness has a significant influence on the structural pounding response, especially on the structural acceleration response and the inertia force response [19][20][21].However, the selection of the impact stiffness in the pounding analysis is a complex issue that could need considering the unknown geometry of the impact surfaces, the uncertain material properties under impact load, the variable impact velocities, and so forth.Very limited research has been carried out for the selection of the impact stiffness.Typically, the impact stiffness in the linear viscoelastic model was empirically taken proportional to the axial stiffness of the colliding structures.This approach has been commonly used to select the impact stiffness in the analysis of structural seismic pounding [17].In general, the impact stiffness used in the bridge pounding analysis is taken as the axial stiffness of girder [22] while the impact stiffness used in the building pounding analysis is taken as the axial stiffness of building slab [14].Cole et al. [23] derived a formula to determine the impact stiffness of the linear viscoelastic model, based on the assumption of the same durations from the collision model of distributed mass and the collision model of equivalent lumped masses.However, the effectiveness of Cole's formula has not been verified sufficiently.
The purpose of this paper is to derive a new formula of impact stiffness for the linear viscoelastic model, based on the assumption of the same maximum impact deformations from the collision model of distributed masses and the collision model of equivalent lumped masses.The effectiveness of the present analytical formula is verified through comparison with the results of numerical analysis.

Formula for Impact Stiffness
2.1.The Maximum Impact Deformation from the Distributed Mass Model.In order to present a better physical understanding of the structural seismic pounding, the distributed mass model has been used to analyze the pounding phenomenon based on the one-dimensional wave propagation theory [23][24][25].The collision of two adjacent components during an earthquake can be simplified as a collision of two distributed masses, as shown in Figure 1.
The distributed mass  ( = 1, 2, resp.) has a section   , a length   , an elastic modulus   , and a mass density   .Based on the results of a single distributed mass colliding with a rigid wall, the motion for the collision of two distributed masses was derived by Cole et al. [24] as follows: in which where   (, ) ( = 1, 2) is the displacement at an arbitrary point  in the distributed mass  at time , which takes zero value at the initial time  = 0; V 0 is the initial velocity of the distributed mass  before impact occurrence; and V  is the intermediate velocity of the collision interface.The wave propagation velocity   of the distributed mass  is as follows: Equations ( 1) are only valid during the period of contact occurrence between two distributed masses.The term V   in (1) is used to describe the rigid body displacement of each mass during the collision and the series terms are used to describe the impact deformation of each mass.Namely, the impact displacements in two distributed masses are taken as  1 (, ) and  2 (, ), respectively.
The natural collision period of the distributed mass  can be written as follows: When two distributed masses collide, the collision duration takes the smaller one of natural collision periods of two distributed masses.The reason is that once  reaches min( 1 ,  2 ), the mass with the shorter period moves away while the other mass internally oscillates because of the incomplete reflection of the shock wave.
Without losing generality, we assume  1 ≤  2 when the pounding between two distributed masses is studied.Therefore, the collision duration is  1 .Furthermore, when a collision occurs, two compression waves from the collision interface will propagate within the distributed masses in opposite directions.As the compression waves pass through the masses, the deformation in each mass occurs.Because the natural collision period of distributed mass 1 is less than that of mass 2, the compression wave is firstly reflected when it reaches the free end of distributed mass 1.At this moment, mass 1 presents its maximum compression deformation.However, the compression deformation in mass 2 is  1 / 2 times the maximum compression deformation itself.The total impact deformation is the sum of the deformations in mass 1 and mass 2. When the compression wave reaches the free end of distributed mass 1, the maximum impact deformation is  dm , which can be obtained from the mathematical sum of the infinite series in (2) as follows: Substituting ( 4) into (5) yields the following:

The Maximum Impact Deformation from the Equivalent
Lumped Mass Model.Two distributed masses given in Figure 1 can be simplified as two equivalent lumped masses during the analysis of seismic pounding, as shown in Figure 2. Therefore, mass   in Figure 2 can be expressed as follows: The pounding phenomenon between two lumped masses as shown in Figure 2 has been analyzed, based on the linear viscoelastic model by Anagnostopoulos [26].The contact force for the linear viscoelastic model can be expressed as follows: where   ( = 1, 2) and u  are the displacements and velocity of the lumped mass , respectively, where the dot denotes the differentiation with respect to time .  is the initial separation distance between two lumped masses. is the impact stiffness of the linear viscoelastic model, and  is the corresponding impact damping coefficient which can be expressed by the following formula [26]: where  denotes an impact damping ratio correlated with the restitution coefficient , which accounts for the energy dissipation during the pounding.Without loss of generality, it is assumed that the beginning of the impact occurrence is at  = 0. Therefore, we have  1 (0) =  2 (0) = 0, u 1 (0) = V 10 , and u 2 (0) = V 20 .When the contact touches off, the motion equation of two lumped masses are as follows: where ü  ( = 1, 2) is the acceleration of the lumped mass .Equation ( 10) can be solved by using the mode superposition method.With the initial above-mentioned conditions, the solutions of the displacements are as follows: in which where  is the radial frequency and   is the damped radial frequency.
The corresponding velocities are as follows: When the maximum deformation occurs during the impact, the relative velocity between two masses should be zero; namely, where   is the time corresponding to the occurrence of maximum deformation during the impact.Substituting ( 13) into ( 14), the expression for   can be obtained as follows: It should be mentioned that   is always smaller than /  which is the time when the two lumped masses begin to separate [26].

The Formula of Impact
Stiffness.Real structural seismic pounding can only have a maximum deformation during a pounding period.However, both the distributed mass model shown in Figure 1 and the lumped mass model shown in Figure 2 can be used to analyze the seismic pounding between the adjacent structures.The maximum impact deformation from the distributed mass model is equal to that from the equivalent lumped mass model; namely, Substituting ( 6) and ( 17) into (18) leads to the following: ⋅  ((ln )/) arcsin(/√ 2 +(ln ) 2 ) . ( Solving the above equation for  yields the following: The axial stiffness of the distributed mass  given in Figure 1 can be expressed as follows: Substituting ( 3), (7), and ( 21) into (20) yields the following: It should be mentioned that formula (22) is valid only when  1 ≤  2 .When  1 >  2 , a similar formula can be obtained by taking  1 ,  2 , and  2 instead of  2 ,  1 , and  1 in (22), respectively.

Numerical Verification
During the derivation of the formula for impact stiffness, some assumptions have been set.Therefore, it is necessary to verify the correctness and the accuracy of the proposed formula (22).The values of impact stiffness from the present formula have been compared with those determined by other methods.Furthermore, the precision of the present formula is verified through numerical simulations for pounding response.
The comparison study among three kinds of methods for impact stiffness has been conducted; that is, taking the impact stiffness follows, respectively, the axial stiffness formula (see (21)), the present formula (see ( 22)), and Cole's formula which is expressed as follows [23]: Moreover, the effect of different parametric values such as the mass ratio of two impact bodies and the restitution coefficient on impact stiffness has been investigated.We take the impact between two simplified concrete bodies as an example.The following parameters in the analysis have been used: and  1 = 40 m, and the mass ratio of two bodies ( 1 ≤  2 ) is taken as follows: In the case of constant mass ratio  = 2, the restitution coefficient  ranged from 0.3 to 1 with an increment of 0.01.Moreover, in the case of fixed restitution coefficient  = 0.65, the mass ratio  has been set from 1 to 5 with an increment of 0.05 (the change of the mass ratio results from the length variation of body 2 in Figure 1).The values of the impact stiffness from three kinds of methods are normalized as follows: where  is the stiffness ratio relative to the axial stiffness of impact body 1.
The comparisons of results from three kinds of methods for impact stiffness are presented in Figures 3 and 4, respectively.It can be seen from both figures that the impact stiffness calculated by the present formula (see (22)) is obviously different from those determined by the other two methods.Moreover, it can be seen from Figure 3 that the impact stiffness given by the proposed formula is much smaller than those from Cole's formula (see (23)) and the axial stiffness formula (see (21)) at the level of restitution coefficient  varying from 0.3 to 1. Similarly, this trend has been observed in Figure 4 where the values given by the proposed formula are much smaller than those given by Cole's formula and the axial stiffness formula.Guo et al. [27,28] and Wang et al. [29] commented that, for the seismic pounding of a bridge, the impact stiffness obtained from experiment is significantly smaller than that from theoretical analysis.Khatiwada et al. [30] studied five different impact force models and carried out shake table tests for pounding between two steel portal frames.It was found that the Hertz-damp model always overestimated the responses while the other four models also frequently overestimated the amplifications.Thus, the predictions from the four models were not significantly different.Since the linear viscoelastic model requires substantially less computation, this model is more suitable for numerical modelling of pounding responses, compared with the other models.It is seen from Figures 3 and 4 that the impact stiffness from the present formula is just lower than that from the axial stiffness formula and greatly lower than that from Cole's formula.Furthermore, it can be seen from Figures 3  and 4 that both restitution coefficient and mass ratio have a significant influence on the impact stiffness described by the proposed formula.The impact stiffness increases with the increase of restitution coefficient or mass rate.
In order to verify the accuracy of the proposed formula for impact stiffness, a comparison of the maximum impact forces from pounding simulations using the linear viscoelastic model has been made.In the numerical analysis, the impact stiffness determined by three kinds of methods, respectively, was used and the same model parameters as given above were taken.In the case of mass ratio  = 2 and the initial relative velocity V 10 − V 20 = 2 m/s, the coefficient of restitution r ranged from 0.4 to 1 with an increment of 0.1.In the case of mass ratio  = 2 and the restitution coefficient  = 0.65, the initial relative velocity (V 10 − V 20 ) was taken from 1 m/s to 5 m/s with an increment of 1 m/s.In the case of the restitution coefficient  = 0.65 and the initial relative velocity V 10 − V 20 = 2 m/s, the mass ratio  varies from 1 to 5 with an increment of 1.The fourth-order Runge-Kutta method with the constant time step Δ = 0.0001 s was applied in the numerical computations.The referenced maximum impact force  ,max used for comparison comes from the following formula [18,19]: The relative error of the maximum impact force from the numerical simulations was assessed by the following: where  is the relative error. max is the maximum impact force from the numerical simulations.The drift envelopes of the relative errors of maximum impact forces from numerical simulations, based on the linear viscoelastic model with respect to impact stiffness given by three kinds of methods, are presented in Figures 5, 6, and 7.It can be seen from Figures 5-7 that the linear viscoelastic model with the impact stiffness coming from the proposed formula (see (22)) gives the least relative errors on the average meaning, compared to the results based on the impact stiffness given by Cole's formula (see (23)) and that given by the axial stiffness formula (see (21)).Therefore, it may be inferred that the linear viscoelastic model using the impact stiffness given by the proposed formula could provide a more precise simulation to the structural pounding than those given by the other two formulae.It should be mentioned that, in order to derive the present formula of impact stiffness, the maximum impact deformation in the distributed mass model is assumed to be the same as that from the equivalent lumped mass model.This equivalence is essentially approximate.However, the numerical comparison has verified the acceptable accuracy of the present formula, at least for the actual engineering application.

Conclusions
In this paper, a new formula to calculate the impact stiffness of linear viscoelastic model for pounding simulation has  been derived by equating the maximum impact deformation from the distributed mass model to that from the equivalent lumped mass model.In order to verify the accuracy of the proposed formula, the comparison of pounding simulations from three different formulae has been conducted.The analysis indicates that the impact stiffness of linear viscoelastic model, provided by the proposed formula, is much smaller than the axial stiffness of the colliding body, which is closer to the experimental observation.The accuracy of the proposed formula has been verified through comparing the numerical simulation results with the analytical solutions.Therefore, more reliable results for pounding simulation in structural engineering could be provided by using the linear viscoelastic model with the present formula of impact stiffness.

Figure 3 :Figure 4 :
Figure 3: The impact stiffness from three kinds of methods with respect to different restitution coefficients.

Figure 5 :Figure 6 :
Figure5: The relative errors of maximum impact forces from numerical simulations based on three kinds of methods to calculate impact stiffness with respect to different restitution coefficients.

Figure 7 :
Figure7: The relative errors of maximum impact forces from numerical simulations based on three kinds of methods to calculate impact stiffness with respect to different mass ratios of two bodies.