Effects of Large Deformation and Velocity Impacts on the Mechanical Behavior of Filled Rubber: Microstructure-Based Constitutive Modeling and Mechanical Testing

Filled rubber has been extensively used in the repairing, retrofitting, and protecting of civil infrastructures due to its superior physical and mechanical properties. However, effects of large deformation and velocity impacts on the mechanical behavior of filled rubber are not well recognized, one of the major challenges in the past investigations is that the material exhibits significant nonlinearity and sensitivity to velocity. This paper presents a hyper-viscoelastic constitutive modeling and experimental study to capture both the hyperelastic and viscoelastic behaviors of filled rubber under large shear deformation and velocity impacts. Motivated by the micro-mechanism of filled rubber, the constitutive modeling consists of an equilibrium element in parallel with an improved Maxwell element to incorporate both nonlinear hyperelasticity and rate-dependent performance governed by the readjustment and rearrangement of molecular chains in the material. A new strain energy function is developed and the physical description of parameters in the strain energy function is highlighted. The Clausius-Duhem inequality is employed to consider the thermodynamic consistency of the model. Then, stress relaxation property and stress-strain response of filled rubber upon cyclic shear loading with different strain rates (ranging from 0.08 to 12.0 s−1) are experimentally studied, and some key observations are summarized. Subsequently, a “Gau-Poly” function is proposed based on the experimental data to describe the viscoelastic property of filled rubber versus strain and strain rate. Finally, stress-strain relationship and hysteretic area obtained from the experimental results were compared with the numerical results of the model, good agreement was achieved and the capacity of the model to accurately reproduce the mechanical behavior of filled rubber under a wide range of deformation and velocity impacts was verified.


Introduction
Filled rubber is an important class of materials due to its favorable flexibility, energy dissipation, self-centering, and other properties [1,2]. It has been extensively employed in isolation bearings, expansion joints, shock absorbers, and other civil engineering applications as a reliable method to mitigate the effects of vibration and dynamic impacts on structures, such as highway bridges and buildings [3][4][5]. The majority of these materials are filled with active fillers, such as carbon black, to remarkably enhance the properties, including stiffness, strength, damping, and abrasion resistance [6,7]. An illustration of the microstructure of a carbon black filled rubber is presented in Figure 1, it shows that the material comprises a rubber matrix, fillers, and rubber-filler interface.
interface. The carbon black particles are randomly distributed in the rubber matrix and are clustered into irregularly shaped aggregates with varying sizes and spacing. Moreover, these components further consist of three-dimensional networks, which include a large number of randomly oriented molecular chains with a broad range of lengths that are often cross-linked at junctions or become entangled among themselves [8,9]. Filled rubber-based structural devices are easily exposed to blast, earthquake, and other dynamic impacts ranging from low to high strain rates [10]. Previous experimental studies have revealed that the filled rubber has a complex nonlinear behavior under large deformation and high velocity impacts [11]. The stress-strain response of filled rubber is highly dependent on deformation rate and exhibits significant hysteresis upon cyclic shear loading, resulting in a series of undesired design problems in engineering practices [12]. Although the bilinear hysteretic model is a widely employed modeling approach, it excludes many of the advanced mechanical behaviors of filled rubber, such as strong nonlinear and hardening properties, thus it would yield inaccurate design results [13]. Therefore, the development of a high fidelity constitutive model of filled rubber under a wide range of strains and strains rates is imperative and contributes to the numerical simulation during the design process.
Currently, several sophisticated models and numerical computation are available in the literature to address different aspects of mechanical properties of filled rubber-e.g., Lion [14], Jankowski [15], Li et al. [16], Gjorgjiev et al. [17], Sahu et al. [18], Khajehsaeid et al. [19], and more recently Mokhireva et al. [20]. However, most available researches focused on the uniaxial tension or compression tests of the material, whereas limited experimental works regarding simple shear deformation have been reported, which essentially play an important role in analyzing the mechanical behavior of rubber-based structural devices under seismic excitation in civil engineering. It was found that the reasons are twofold: firstly, research interest in the shear deformation has merely grown in recent years. Secondly, additional measurement difficulty is encountered when filled rubber is under more challenging conditions like large deformation and high velocity impacts [21]. More importantly, those models have a drawback in that their estimations are unsatisfactory under large strains and high velocity circumstances. The viscosity effect of the material is roughly expressed via a Prony series approach [22], yet the nature and physical basis of the behavior of filled rubber is not well recognized.
On the other hand, there is another approach to develop the constitutive model based on the characteristics of microstructure of filled rubber and micro-mechanics owing to the fact that its mechanical behavior is mainly governed by the different cross-linking networks in the material [23]. Tomita et al. investigated the monotonic and cyclic behavior of carbon black filled rubber through a molecular-chain network model, but the time-dependent nature of the material is yet to be addressed [24]. Bergström et al. studied the microstructure of carbon black filled chloroprene rubber and concluded that the mechanical response of filled rubber incorporates equilibrium and rate-dependent components [25]. Pouriayevali et al. investigated the relationship between the readjustment of the molecular chain in rubber and its stress relaxation property in the macro level. Filled rubber-based structural devices are easily exposed to blast, earthquake, and other dynamic impacts ranging from low to high strain rates [10]. Previous experimental studies have revealed that the filled rubber has a complex nonlinear behavior under large deformation and high velocity impacts [11]. The stress-strain response of filled rubber is highly dependent on deformation rate and exhibits significant hysteresis upon cyclic shear loading, resulting in a series of undesired design problems in engineering practices [12]. Although the bilinear hysteretic model is a widely employed modeling approach, it excludes many of the advanced mechanical behaviors of filled rubber, such as strong nonlinear and hardening properties, thus it would yield inaccurate design results [13]. Therefore, the development of a high fidelity constitutive model of filled rubber under a wide range of strains and strains rates is imperative and contributes to the numerical simulation during the design process.
Currently, several sophisticated models and numerical computation are available in the literature to address different aspects of mechanical properties of filled rubber-e.g., Lion [14], Jankowski [15], Li et al. [16], Gjorgjiev et al. [17], Sahu et al. [18], Khajehsaeid et al. [19], and more recently Mokhireva et al. [20]. However, most available researches focused on the uniaxial tension or compression tests of the material, whereas limited experimental works regarding simple shear deformation have been reported, which essentially play an important role in analyzing the mechanical behavior of rubber-based structural devices under seismic excitation in civil engineering. It was found that the reasons are twofold: firstly, research interest in the shear deformation has merely grown in recent years. Secondly, additional measurement difficulty is encountered when filled rubber is under more challenging conditions like large deformation and high velocity impacts [21]. More importantly, those models have a drawback in that their estimations are unsatisfactory under large strains and high velocity circumstances. The viscosity effect of the material is roughly expressed via a Prony series approach [22], yet the nature and physical basis of the behavior of filled rubber is not well recognized.
On the other hand, there is another approach to develop the constitutive model based on the characteristics of microstructure of filled rubber and micro-mechanics owing to the fact that its mechanical behavior is mainly governed by the different cross-linking networks in the material [23]. Tomita et al. investigated the monotonic and cyclic behavior of carbon black filled rubber through a molecular-chain network model, but the time-dependent nature of the material is yet to be addressed [24]. Bergström et al. studied the microstructure of carbon black filled chloroprene rubber and concluded that the mechanical response of filled rubber incorporates equilibrium and rate-dependent components [25]. Pouriayevali et al. investigated the relationship between the readjustment of the molecular chain in Polymers 2020, 12, 2322 3 of 18 rubber and its stress relaxation property in the macro level. Based on this, they developed a relaxation time model to describe the strain rate sensitivity of rubber material [26].
Following the argument of Bergström et al. and inspired by the micro-mechanical behavior of filled rubber, the main objective of this study is to propose a constitutive model to capture both hyperelastic and viscoelastic properties of filled rubber under extensive ranges of strains and strain rates, especially attempting to consider the influence of fillers on altering the behavior of the material and to develop the viscosity law in a thermodynamically consistent way. The first part of the model comprised of a hyperelastic equation based on a newly developed strain energy function to represent the rate-independent response. The other part of the model incorporated an improved Maxwell element to consider the nonlinearity and rate-dependent response of filled rubber, multiplicative kinematics decomposition of deformation gradient tensor into elastic and inelastic parts was conducted to correlate the overstress tensor with strain rate. Besides, the Clausius-Duhem inequality was employed to ensure the thermodynamic consistency [27], the Helmholz free energy was separately stored in the above-mentioned two parts of the model, and an associated nonlinear viscosity coefficient was derived. Subsequently, a strain-controlled mechanical testing program covering large deformation (200% shear strain) and various velocity impacts (strain rate ranging from 0.08 to 12.0 s −1 ) was performed at room temperature to investigate the stress-strain relationship of the material. Then, hyperelastic parameters were identified straightforward based on the testing results and a "Gau-Poly" function was proposed to describe the viscoelastic property of filled rubber versus strain and strain rate. Finally, numerical results from the model were compared with experimental data.

Constitutive Modeling
Motivated by the micro-mechanism of filled rubber, a hyper-viscoelastic model is proposed to incorporate its hyperelastic and viscoelastic characteristics. The model is decomposable in two parallel parts, as depicted in Figure 2-part A is modeled by an equilibrium spring A to represent the rate-independent equilibrium response, while part B is consists of a generalized Maxwell element to describe the rate-dependent instantaneous response. It is an improvement on the traditional Maxwell element to incorporate nonlinearity in both spring and dashpot elements, contributing to model the rapid readjustments of short chains as well as the sluggish rearrangement and excessive entanglement of long chains in the material. Based on this, they developed a relaxation time model to describe the strain rate sensitivity of rubber material [26]. Following the argument of Bergström et al. and inspired by the micro-mechanical behavior of filled rubber, the main objective of this study is to propose a constitutive model to capture both hyperelastic and viscoelastic properties of filled rubber under extensive ranges of strains and strain rates, especially attempting to consider the influence of fillers on altering the behavior of the material and to develop the viscosity law in a thermodynamically consistent way. The first part of the model comprised of a hyperelastic equation based on a newly developed strain energy function to represent the rate-independent response. The other part of the model incorporated an improved Maxwell element to consider the nonlinearity and rate-dependent response of filled rubber, multiplicative kinematics decomposition of deformation gradient tensor into elastic and inelastic parts was conducted to correlate the overstress tensor with strain rate. Besides, the Clausius-Duhem inequality was employed to ensure the thermodynamic consistency [27], the Helmholz free energy was separately stored in the above-mentioned two parts of the model, and an associated nonlinear viscosity coefficient was derived. Subsequently, a strain-controlled mechanical testing program covering large deformation (200% shear strain) and various velocity impacts (strain rate ranging from 0.08 to 12.0 s −1 ) was performed at room temperature to investigate the stress-strain relationship of the material. Then, hyperelastic parameters were identified straightforward based on the testing results and a "Gau-Poly" function was proposed to describe the viscoelastic property of filled rubber versus strain and strain rate. Finally, numerical results from the model were compared with experimental data.

Constitutive Modeling
Motivated by the micro-mechanism of filled rubber, a hyper-viscoelastic model is proposed to incorporate its hyperelastic and viscoelastic characteristics. The model is decomposable in two parallel parts, as depicted in Figure 2-part A is modeled by an equilibrium spring A to represent the rate-independent equilibrium response, while part B is consists of a generalized Maxwell element to describe the rate-dependent instantaneous response. It is an improvement on the traditional Maxwell element to incorporate nonlinearity in both spring and dashpot elements, contributing to model the rapid readjustments of short chains as well as the sluggish rearrangement and excessive entanglement of long chains in the material. When loading rate is very slow, no stress is transferred through the dashpot, and the equilibrium spring A is the only source of stress tensor. Conversely, when an infinitely high rate is applied, the dashpot C does not have enough time to complete the relaxation processes-all the deformation is undergone by the intermediate spring B. These two boundary states, respectively, correspond to equilibrium and instantaneous responses of filled rubber, which can be described by hyperelastic equations based on the strain energy function. Between the limiting states, the behavior When loading rate is very slow, no stress is transferred through the dashpot, and the equilibrium spring A is the only source of stress tensor. Conversely, when an infinitely high rate is applied, the dashpot C does not have enough time to complete the relaxation processes-all the deformation is undergone by the intermediate spring B. These two boundary states, respectively, correspond to equilibrium and instantaneous responses of filled rubber, which can be described by hyperelastic equations based on the strain energy function. Between the limiting states, the behavior of the generalized Maxwell element is rate-dependent and the stress-strain behavior of the material is governed by a hyper-viscoelastic constitutive law.

Hyperelasticity
Considering a particle located at position X in a material, when simple shear deformation is applied on the body, as illustrated in Figure 3, new position x of the particle in the deformed configuration can be described by the reference configuration as: where γ is the amount of shear strain along the X 1 direction. of the generalized Maxwell element is rate-dependent and the stress-strain behavior of the material is governed by a hyper-viscoelastic constitutive law.

Hyperelasticity
Considering a particle located at position X in a material, when simple shear deformation is applied on the body, as illustrated in Figure 3, new position x of the particle in the deformed configuration can be described by the reference configuration as: where γ is the amount of shear strain along the X1 direction. Using Equation (1), total deformation gradient tensor F in the finite strain kinematic framework is defined as: Further, the left Cauchy-Green deformation tensor B is obtained by: where F T is the transposed matrix of F. If Equation (3) holds, the first, second, and third strain invariants of B are: where tr(·) and det(·) are trace and determinant operators, respectively.
Deformed configuration Legend: Using Equation (1), total deformation gradient tensor F in the finite strain kinematic framework is defined as: Further, the left Cauchy-Green deformation tensor B is obtained by: where F T is the transposed matrix of F. If Equation (3) holds, the first, second, and third strain invariants of B are: where tr(·) and det(·) are trace and determinant operators, respectively.
The bulk modulus of filled rubber is much greater than its shear modulus, suggesting that a hypothesis of Poisson's ratio ν= 0.5 is reasonable. Thus, applying isotropy and incompressibility yields: where S and T are weighted Cauchy and Cauchy stress tensors, respectively. Subscript E represents the extra part of S. p is the hydrostatic pressure of S; 1 is a second-order identity tensor. Following the framework of hyper-viscoelastic model, S E at any time is the sum of equilibrium stress tensor S eq E and over stress tensor S ov E : The stress tensor is determined by the differentiation of strain energy function (SEF) with respect to B:  [28,29]. Although the abovementioned models provide acceptable agreement with the experimental results in extension and compression modes, several drawbacks of these models are revealed, especially the fact that they fail to accurately reproduce the highly nonlinear stress-strain response of filled rubber under large shear deformation and high velocity impacts [30]. Therefore, a new SEF should be developed based on the mechanical characteristic of filled rubber under considered conditions in the study.
Generally, the SEF is expressed as the sum of two functions associated with I 1 and I 2 [31,32], where In order to ensure the proposed SEF is capable of representing the basic situations of filled rubber, it must satisfy the following requirements. Firstly, the SEF should vanish for the undeformed configuration of filled rubber.
Secondly, the SEF and obtained stress tensors approached infinity when filled rubber is subjected to infinite deformation. lim Moreover, a good SEF should have a simple mathematical structure with minimal number of parameters to describe the essential features of the material and to ensure the easy application. Based on these considerations, a five-constant polynomial SEF is proposed in Equations (15) and (16) to consider the "S" shaped strong nonlinear behavior of filled rubber associated with the complex entropic rearrangement of molecular chains in the material.
Polymers 2020, 12, 2322 6 of 18 where C eq i (i = 1~5) are responsible for the material parameters of equilibrium stress, while C ov i (i = 1~5) are those of overstress.
It is seen that the proposed SEF is also a particular case of the well-known Rivlin's expression: When small deformation like initial linear range is applied on filled rubber, the SEF is simplified as: Noting that the Neo-Hookean model is calculated by the material chain density n, Boltzman's coefficient k, and absolute temperature T [33], which is shown in Equation (19) The Neo-Hookean model is based on the statistical thermodynamics of molecular chains in the material, proceeding to compare Equation (18) and Equation (19) leads to: therefore, the physical meaning of parameters in the SEF is highlighted, C 1 is related to parameters of n, k, and T.
On the other hand, the behavior of filled rubber under large deformation is mainly dominated by the parameters C 4 and C 5 in the proposed SEF. From the perspective of the microstructure of the material, large deformation involves a greater proportion of molecular chains and is determined by their limiting extensibility.
In addition, it should be noted that ignoring I 2 in the proposed SEF may induce some but not significant errors in describing the advanced aspect like poynting-type effect in rubber material [34], it suffices to illustrate the core concept in the study.
Substituting Equation (15) into Equation (8), the hyperelastic relation for the shear component of equilibrium stress tensor is given by: where C eq i (i = 1~5) are rate-independent material parameters. Substitution of Equation (16) into Equation (9) yields the shear component of overstress tensor: where C ov i (i = 1~5) are rate-dependent material parameters.

Hyper-Viscoelastic
As discussed previously, between the two limiting states (i.e., very slow and fast loading rates), the Maxwell element can be decomposed within the framework of a multiplicative kinematics decomposition of F [35], as shown in Equation (24) and Figure 4.
where F B and F C represent the deformation gradient tensors of spring B and dashpot C, respectively.
where C eq i (i=1~5) are rate-independent material parameters. Substitution of Equation (16) into Equation (9) yields the shear component of overstress tensor: where C ov i (i=1~5) are rate-dependent material parameters.

Hyper-Viscoelastic
As discussed previously, between the two limiting states (i.e., very slow and fast loading rates), the Maxwell element can be decomposed within the framework of a multiplicative kinematics decomposition of F [35], as shown in Equation (24) and Figure 4. (24) where FB and FC represent the deformation gradient tensors of spring B and dashpot C, respectively.
After taking the material time derivative of BB and FB, we have where L is the total velocity gradient tensor and LC is the corresponding velocity gradient tensor of dashpot C, which are derived as

Reference configuration
Current configuration Following Equation (24), the left Cauchy-Green deformation tensors of intermediate spring B and dashpot C are deduced, respectively.
After taking the material time derivative of B B and F B , we have .
where L is the total velocity gradient tensor and L C is the corresponding velocity gradient tensor of dashpot C, which are derived as where Substituting Equations (28)-(32) into Equation (27), the following equation is derived. For an isothermal viscoelastic process, the second law of thermodynamics in the form of the Clausius-Duhem inequality requires that [36] where ψ is the Helmholtz free energy. According to the hyper-viscoelastic model, it is noted that the energy can be only stored in spring A and B, thus, ψ has the following form: Taking the material time derivative of ψ, the rate of free energy is expressed as: considering . .
thus, Equation (36) can be rewritten as: On the other hand, according to Equation (24) and Equation (29), the velocity gradient tensor can be decomposed into: Considering Equation (7), S E is the sum of S eq E and S ov E , thus By the analogy with Equation (42), the second term on the right side of Equation (40) becomes: Substituting Equations (38)-(43) into Equation (34), the inequality function can be transformed into: With the definition of equilibrium stress tensor and overstress tensor expressed in Equations (8) and (9), Equation (44) is reduced to the following equation, implying that the power between overstress tensor and velocity gradient tensor of dashpot C should be non-negative. For an arbitrary process, the simplest sufficient condition satisfying Equation (45) is: where η is a positive viscosity coefficient. For the sake of computational simplicity, the validity of Equation (46) is equivalent to the following equation [35]: operating on the inelastic intermediate configuration, we have: where superscript D is the deviatoric component of overstress tensor S ov E . Substituting Equation (48) into Equation (33) yields: Overall, a thermodynamically consistent hyper-viscoelastic constitutive model is established from a combination of two classes of equations. The first corresponds to the limiting states of the material defined by Equations (21)- (23). The second is linked to the rate-dependent response of the material and associated with part B in the model.

Specimen Preparation and Test Setup
Commercially available nitrile rubber was chosen as the base rubber and carbon black N330 was employed as the reinforcing filler to fabricate the specimens. According to the ISO 1827:2007 standard [37], Figure 5 shows the quadruple-lap shear specimen which comprises of four filled rubber layers and rigid steel blocks, the rubber layers are glued to the steel blocks by adhesive material to prevent sliding on the contacting surfaces.  Figure 6 depicts the geometrical sizes of the filled rubber and steel block in the specimen. Each rubber layer has the following dimensions: 25 mm length, 20 mm width, and 4 mm thickness, the length and width are much greater than the thickness to restrict end rotation of rubber layers during the shear deformation and ensure homogeneous strain under high loading rate condition.
(a) Plan view 20 24.8 80 25 Figure 5. Illustration of the quadruple-lap shear specimen. Figure 6 depicts the geometrical sizes of the filled rubber and steel block in the specimen. Each rubber layer has the following dimensions: 25 mm length, 20 mm width, and 4 mm thickness, the length and width are much greater than the thickness to restrict end rotation of rubber layers during the shear deformation and ensure homogeneous strain under high loading rate condition.  Figure 6 depicts the geometrical sizes of the filled rubber and steel block in the specimen. Each rubber layer has the following dimensions: 25 mm length, 20 mm width, and 4 mm thickness, the length and width are much greater than the thickness to restrict end rotation of rubber layers during the shear deformation and ensure homogeneous strain under high loading rate condition.  As presented in Figure 7, a servo-hydraulic MTS 831 elastomer testing machine was employed to provide required loading schemes and measure the corresponding shear force. The quadruple-lap shear specimen was mounted on the apparatus by the upper and lower grips, specified displacement was applied at the top steel block to allow it move freely in the vertical direction with a wide range of strain rates from 0.08 s -1 to 12.0 s -1 . The filled rubber deforms under simple shear deformation while the steel blocks are assumed to be undeformed in the tests. In order to remove the Mullins effect and obtain a relatively stable state of stress-strain curves of the specimens, they were subjected to ten cycles with a maximum shear strain of 200% before the formal tests. All the tests were performed at room temperature. As presented in Figure 7, a servo-hydraulic MTS 831 elastomer testing machine was employed to provide required loading schemes and measure the corresponding shear force. The quadruple-lap shear specimen was mounted on the apparatus by the upper and lower grips, specified displacement was applied at the top steel block to allow it move freely in the vertical direction with a wide range of strain rates from 0.08 s −1 to 12.0 s −1 . The filled rubber deforms under simple shear deformation while the steel blocks are assumed to be undeformed in the tests. In order to remove the Mullins effect and obtain a relatively stable state of stress-strain curves of the specimens, they were subjected to ten cycles with a maximum shear strain of 200% before the formal tests. All the tests were performed at room temperature.

Stress Relaxation and Cyclic Shear Tests
To assess the stress relaxation behavior of filled rubber, generally speaking, a prescribed step strain is exerted on the specimen and holds for a period ranging from 20-60 min to complete various forms of molecular chain readjustments in the material and reach a relatively stable state. Then, the holding period results in an asymptotical reduction in the stress to approximately one third or less of the stress at the beginning of relaxation, implying that the stress achieves an equilibrium state [38]. Once the holding period was finished, the procedure was repeated by increasing the strain amplitude.

Stress Relaxation and Cyclic Shear Tests
To assess the stress relaxation behavior of filled rubber, generally speaking, a prescribed step strain is exerted on the specimen and holds for a period ranging from 20-60 min to complete various forms of molecular chain readjustments in the material and reach a relatively stable state. Then, the holding period results in an asymptotical reduction in the stress to approximately one third or less of the stress at the beginning of relaxation, implying that the stress achieves an equilibrium state [38]. Once the holding period was finished, the procedure was repeated by increasing the strain amplitude.
In this study, the specimens were successively loaded to strain levels of 25%, 75%, 125%, 175%, and 200%, and held for 1200 s during each step. Figure 8 shows the stress relaxation results of filled rubber, a time shift of 150 s was introduced in the figure to highlight the comparison. It was observed that, during each holding step, the stress undergoes a rapid decrease at the very beginning of the process, almost 50% of the stress relaxation occurs in the first 5 s. Then, the stress gradually converges to a limit value, accompanied by the occurrence of rearrangement of microstructure in the filled rubber. At the end of the holding step, the stress at the termination point was identified as the approximation of equilibrium stress.

Stress Relaxation and Cyclic Shear Tests
To assess the stress relaxation behavior of filled rubber, generally speaking, a prescribed step strain is exerted on the specimen and holds for a period ranging from 20-60 min to complete various forms of molecular chain readjustments in the material and reach a relatively stable state. Then, the holding period results in an asymptotical reduction in the stress to approximately one third or less of the stress at the beginning of relaxation, implying that the stress achieves an equilibrium state [38]. Once the holding period was finished, the procedure was repeated by increasing the strain amplitude.
In this study, the specimens were successively loaded to strain levels of 25%, 75%, 125%, 175%, and 200%, and held for 1200 s during each step. Figure 8 shows the stress relaxation results of filled rubber, a time shift of 150 s was introduced in the figure to highlight the comparison. It was observed that, during each holding step, the stress undergoes a rapid decrease at the very beginning of the process, almost 50% of the stress relaxation occurs in the first 5 s. Then, the stress gradually converges to a limit value, accompanied by the occurrence of rearrangement of microstructure in the filled rubber. At the end of the holding step, the stress at the termination point was identified as the approximation of equilibrium stress. As shown in Figure 9, cyclic shear tests consist of saw-tooth wave displacement excitation with fully reversed cycles of loading at a maximum strain of 200% to consider general cases of filled As shown in Figure 9, cyclic shear tests consist of saw-tooth wave displacement excitation with fully reversed cycles of loading at a maximum strain of 200% to consider general cases of filled rubber-based structural devices used in civil engineering projects. Six different strain rates, including 0.08, 0.4, 0.8, 4.0, 8.0, and 12.0 s −1 were employed to investigate the cyclic mechanical performance of filled rubber under various velocity impacts, covering most typical dynamic scenarios that filled rubber undergoes during earthquake excitations. Among the loading rates, the strain rate = 12.0 s −1 was used as the approximation and alternative solution of the theoretical instantaneous response to calibrate the material parameters associated with the intermediate spring B owing to the fact that an infinitely fast rate was inaccessible because of the limitation of testing machine.   Figure 10 presents the cyclic shear test results of filled rubber, it is clearly found that the stress-strain responses not only exhibit pronounced nonlinearity with respect to small, moderate, and large strains, but also are very sensitive to the loading rate-specifically, the strain rate dependence is more obvious in the loading process than the unloading process. For any given strain level in the positive loading procedure, the shear stress shows a significant increasing trend when −1 Figure 9. Applied strain histories in the cyclic shear tests. Figure 10 presents the cyclic shear test results of filled rubber, it is clearly found that the stress-strain responses not only exhibit pronounced nonlinearity with respect to small, moderate, and large strains, but also are very sensitive to the loading rate-specifically, the strain rate dependence is more obvious in the loading process than the unloading process. For any given strain level in the positive loading procedure, the shear stress shows a significant increasing trend when the loading rate increases from 0.08 to 12.0 s −1 . More specifically, the maximum shear stress in the case of strain rate = 0.08 s −1 is 1.94 MPa, while that of strain rate = 12.0 s −1 is 2.86 MPa, the corresponding increment of the stress is 47.4%. These unique phenomena can be attributed to the presence of high amount of fillers which alter the microstructure of the material and eventually result in a viscoelastic behavior in the filled rubber. Figure 9. Applied strain histories in the cyclic shear tests. Figure 10 presents the cyclic shear test results of filled rubber, it is clearly found that the stress-strain responses not only exhibit pronounced nonlinearity with respect to small, moderate, and large strains, but also are very sensitive to the loading rate-specifically, the strain rate dependence is more obvious in the loading process than the unloading process. For any given strain level in the positive loading procedure, the shear stress shows a significant increasing trend when the loading rate increases from 0.08 to 12.0 s −1 . More specifically, the maximum shear stress in the case of strain rate =0.08 s −1 is 1.94 MPa, while that of strain rate =12.0 s −1 is 2.86 MPa, the corresponding increment of the stress is 47.4%. These unique phenomena can be attributed to the presence of high amount of fillers which alter the microstructure of the material and eventually result in a viscoelastic behavior in the filled rubber. In addition, geometric dimensions of the specimens were examined after the tests. It was observed that the filled rubber specimens almost retain their original geometric dimensions (i.e., the residual strains are negligible) after the loading tests, indicating the observed phenomenon is amenable to description by the proposed hyper-viscoelastic model. In addition, geometric dimensions of the specimens were examined after the tests. It was observed that the filled rubber specimens almost retain their original geometric dimensions (i.e., the residual strains are negligible) after the loading tests, indicating the observed phenomenon is amenable to description by the proposed hyper-viscoelastic model.

Parameter Identification
The identification for the model parameters consists of two procedures, the first is related to the hyperelastic parameters and the second focuses on the viscoelastic parameters.

Hyperelastic Parameters
Parameters of equilibrium stress were directly identified by minimizing the error between the stress relaxation test results and simulated results associated with Equation (21). It should be noted that C eq 5 was initially taken to be equal to zero because the first four items on the right side of Equation (21) were sufficiently to fit the stress relaxation results of filled rubber. While the parameter calibration for overstress was based on three steps, firstly, the case of strain rate = 12.0 s −1 in the cyclic shear tests was selected considering the approximation of infinitely fast rate. Then, the overstress was extracted from the total stress by subtracting the equilibrium stress. Finally, parameters could be estimated in accordance with Equation (23). By applying nonlinear least-squares method, the identified parameters are listed in Table 1.  Figure 11 presents the experimental and fitted results of the viscosity coefficient. In particular, the experimental results are first focused-it is easily observed that the viscosity coefficient is positively correlated with the strain rate. For a specified strain level, the viscosity coefficient is increased with the increment of strain rate-this variation is more pronounced in high rate conditions (i.e., strain rate > 0.8 s −1 ). While, for a given strain rate, it is also seen that the viscosity coefficient achieves the maximum value at moderate strain levels, whereas the value is decreased at small and large strains.  Figure 11 presents the experimental and fitted results of the viscosity coefficient. In particular, the experimental results are first focused-it is easily observed that the viscosity coefficient is positively correlated with the strain rate. For a specified strain level, the viscosity coefficient is increased with the increment of strain rate-this variation is more pronounced in high rate conditions (i.e., strain rate>0.8 s −1 ). While, for a given strain rate, it is also seen that the viscosity coefficient achieves the maximum value at moderate strain levels, whereas the value is decreased at small and large strains.  Therefore, it is concluded that the viscosity coefficient η is closely related to the strain rate v and strain γ (i.e., η = η (γ, v)), this phenomenon is also experimentally demonstrated for styrene butadiene rubber by Fatt et al. [39]. Motivated by the experimental results, a fitting function is introduced to characterize how the viscosity coefficient was affected by the strain and strain rate. In this function, a Gaussian function of γ and v was implemented in conjunction with a second-order polynomial function of v to develop an accurate and robust model for the reproduction of multi-dependent behavior of η, which is called "Gau-Poly" function and expressed in the following form:

Viscoelastic Parameters
where parameters z 0 , A, x c , w 1 , y c , and w 2 in the Gaussian function govern the multi-dependence of η, while the rest of the parameters b, c, and d are served to describe the strong variation of η with strain rate. It is noted that x c , w 1 , and b are non-unit parameters. The parameter identification of viscosity coefficient was done by the following double sub-processes. Firstly, the rate-dependent characteristic of filled rubber was temporarily neglected under a special condition of strain rate = 0.08 s −1 and it was reasonable to deem that the η was only determined by the Gaussian function in this case [40]. Once the parameters belonging to the Gaussian function were determined, the second step was conducted to identify the rest parameters by approximating the rate-dependent variation in η. This arises from the fact that the viscosity coefficient is nonlinearly proportional to the strain rate. The identified parameters and their units are summarized in Table 2-reasonable agreement between the experimental and fitted results is shown in Figure 11, demonstrating that the function has a good generalization ability to predict the variation in viscosity coefficient with varying strains and strain rates.

Model Application, Validation, and Discussion
In this section, numerical results of the proposed model are compared with the experimental tests to validate how it works in predicting the mechanical behavior of the material.

Monotonic Shear Tests and Model Predictions
To validate the applicability of the proposed constitutive model in predicting the rate-dependent and nonlinear properties of filled rubber, monotonic shear tests with three levels of strain rates, including 0.08, 0.8, and 8.0 s −1 were simulated, roughly corresponding to slow, moderate, and high velocity scenarios. By incorporating the proposed model in a MATLAB software code, Figure 12 presents the comparison of experimental results and model predictions. It is seen that the numerical simulations can accurately reproduce the stress-strain response of filled rubber-a good agreement including a rapid stiffness growth at vary small strain amplitudes and rate-dependent performance for different strain rate cases is obtained. Among all the loading schemes, the discrepancy is within the range of 6% and mainly exists at very small strains, resulting from the fact that the fillers substantially increase the initial stiffness of filled rubber.  Figure 13 shows the comparison of experimental results and model predictions for the cyclic shear tests with two levels of strain rates, including 0.4 and 4.0 s −1 -coherent agreement in terms of the stress-strain responses of the material between the experimental and numerical results is observed in both cases, indicating the favorable capacity and reliability of the proposed model in predicting the rate-dependent and cyclic behavior of filled rubber. Notably, the proposed model presents an error in small deformations. For instance, the stress predicted by the model was approximately 7% greater than that of experimental results at appropriately 40% shear strain in the case of strain rate =4.0 s −1 . It should be noted that the error merely exists in local regions and is not  Figure 13 shows the comparison of experimental results and model predictions for the cyclic shear tests with two levels of strain rates, including 0.4 and 4.0 s −1 -coherent agreement in terms of the stress-strain responses of the material between the experimental and numerical results is observed in both cases, indicating the favorable capacity and reliability of the proposed model in predicting the rate-dependent and cyclic behavior of filled rubber. Notably, the proposed model presents an error in small deformations. For instance, the stress predicted by the model was approximately 7% greater than that of experimental results at appropriately 40% shear strain in the case of strain rate = 4.0 s −1 . It should be noted that the error merely exists in local regions and is not enlarged after the initial deformation. The reasons for this discrepancy are twofold: one is that the strain energy function of model yields a less satisfactory agreement with the strong nonlinear behavior of the material at small strains; the other is the appropriation of the viscosity coefficient.  Figure 13 shows the comparison of experimental results and model predictions for the cyclic shear tests with two levels of strain rates, including 0.4 and 4.0 s −1 -coherent agreement in terms of the stress-strain responses of the material between the experimental and numerical results is observed in both cases, indicating the favorable capacity and reliability of the proposed model in predicting the rate-dependent and cyclic behavior of filled rubber. Notably, the proposed model presents an error in small deformations. For instance, the stress predicted by the model was approximately 7% greater than that of experimental results at appropriately 40% shear strain in the case of strain rate =4.0 s −1 . It should be noted that the error merely exists in local regions and is not enlarged after the initial deformation. The reasons for this discrepancy are twofold: one is that the strain energy function of model yields a less satisfactory agreement with the strong nonlinear behavior of the material at small strains; the other is the appropriation of the viscosity coefficient.  Considering the hysteretic behavior (i.e., energy dissipation capacity) of the filled rubber plays a crucial role in civil engineering, the hysteretic area (unit of the area is ignored for simplicity) of stress-strain curves of the specimens obtained from Figure 13 was calculated and is given in Table 3 for a further comparison. It is found that the predicted hysteretic area is in good agreement with the experimental results, the relative error for both cases is within 5%. Overall, the agreement is sufficiently to substantiate that the proposed model is capable of yielding a reasonable description in engineering practices when a filled rubber-based structural device is incorporated.

Conclusions
This paper presented a micro-mechanism-based hyper-viscoelastic constitutive model to describe the behavior of filled rubber under large shear deformation and various velocity impacts. Parameters in the model were identified by conducting experimental tests and the accuracy of the model was verified by comparing the experimental and numerical results. Major conclusions can be summarized as below: 1.
The constitutive model comprises two parts: the first captures the equilibrium and instantaneous responses of filled rubber; the second incorporates the decomposition of the deformation gradient tensor to correlate the overstress tensor with strain rate. The proposed model covers a wide range of conditions, including small to large shear deformations as well as low to high velocity impacts that the filled rubber is expected to undergo in engineering practices.

2.
Considering the similar nature (i.e., nonlinear hyperelasticity) of equilibrium and instantaneous responses of filled rubber, a newly-developed polynomial strain energy function is applied for equilibrium and intermediate springs in the constitutive model. The proposed strain energy