# Influence of large strain preloads on the viscoelastic response of rubber-like materials under small oscillations

Anne-Sophie Lectez, Erwan Verron

► To cite this version:

Anne-Sophie Lectez, Erwan Verron. Influence of large strain preloads on the viscoelastic response of rubber-like materials under small oscillations. International Journal of Non-Linear Mechanics, 2016, 81, 10.1016/j.ijnonlinmec.2015.12.003. hal-03286600

HAL Id: hal-03286600

https://hal.science/hal-03286600

Submitted on 28 Jul 2021

HAL is a multi-disciplinary open access
archive for the deposit and dissemination of sci-
entific research documents, whether they are pub-
lished or not. The documents may come from
teaching and research institutions in France or
abroad, or from public or private research centers.

L'archive ouverte pluridisciplinaire HAL, est
destinée au dépôt et à la diffusion de documents
scientifiques de niveau recherche, publiés ou non,
émanant des établissements d'enseignement et de
recherche français ou étrangers, des laboratoires
publics ou privés.

<figure><img src="image_1.png" /></figure>

Distributed under a Creative Commons Attribution 4.0 International License
# Influence of large strain preloads on the viscoelastic response of rubber-like materials under small oscillations

A.-S. Lectez*, E. Verron

Institut de Recherche en Génie Civil et Mécanique (GeM), UMR CNRS 6183, Ecole Centrale de Nantes, BP 92101, 44321 Nantes cedex 3,
France

The viscoelastic response predicted by linearized internal variables models in the case of small oscilla-
tions superimposed on a large static preload is investigated, comparing simple forms of the Zener and
the Poynting–Thomson models. It is shown that both of them predict a preload dependency of the
equilibrium linearized stress, but only the latter take into account such a dependency on the out-of-
equilibrium part. Yet, the Zener model is much more frequently linearized as the Poynting–Thomson
model. The formulation of each model for finite deformations is quickly reminded, before linearizing
them around a large static preload. Finally, a comparison of the influence of preload on each model is
proposed for uniaxial extension, before discussing which kind of model has to be chosen regarding
theoretical and practical aspects.

**1. Introduction**

Elastomers are often used in industrial applications for which
vibrations must be filtered, as in automotive engine mounts. These
parts are subjected to a large strain preload due to engine weight,
on which some oscillations caused by engine vibrations and/or
road conditions are superimposed. Simulating accurately the
mechanical response of engine mounts is an industrial ongoing
challenge which necessitates to predict the change in dynamic
properties of elastomers due to the preload. Authors interested in
this problem mostly choose to linearize viscoelastic models.

This linearization is justified by several arguments. Firstly, the
amplitude of oscillations is small for a large part of industrial
applications. In this case, a complex nonlinear viscoelastic model
including several superimposed large strain steps is unnecessary.
Secondly, linearized models are sufficient to study the influence of
a preload on dynamic properties, as it is shown in this paper.
Finally, characterizing a material subjected to small oscillations
enables the calculation of well-known intrinsic mechanical
quantities, like storage and loss moduli. These quantities are only
rigorously defined for small strain, as highlighted by Govindjee
and Simo [1]. Moreover, linearization leads to a relatively simple
closed-form solution, which ensure a faster identification
procedure.

Viscoelastic models are classically derived following two
approaches, referred to as integral and internal variables approa-
ches. Integral models are based on the assumption that the stress
is an explicit function of strain history. This approach generalizes
the Boltzmann principle (see Ward [2] for example). Green and
Rivlin [3] and Coleman and Noll [4] derived such early models,
while Christensen [5] and K-BKZ models [6,7] were developed
later. Internal variables models consider that the stress depends on
state variables, including internal ones. These internal variables
depend on strain history through evolution laws. Rheological
models are parts of this approach. The simplest constitutive
equations of this type are Maxwell and Kelvin-Voigt models
which simply associate a spring and a dashpot, respectively in
series and in parallel. To correctly predict the behavior of a vis-
coelastic solid, a supplementary branch is needed, which leads to
Zener (parallel branches) and Poynting-Thomson (series bran-
ches) models (see for example Wineman and Rajagopal [8]). One
of the first large strain models issued from this approach is due to
Green and Tobolsky [9]. Further developments have been pro-
posed by Lubliner [10] and Sidoroff [11], from which comes the
multiplicative decomposition of the deformation gradient [12] as
applied to viscoelasticity.

Several authors studied the linearization of viscoelastic models,
most of them choosing the integral approach, and particularly
models based on Coleman and Noll [4] proposal [13–18] and on
Pipkin and Rogers [19] proposal [20]. These models, expressed in
their most general form, contain a preload dependency of both
equilibrium and out-of-equilibrium stresses. However, application

<footnote>* Corresponding author. Current affiliation: IRT Jules Verne, Chemin du Chaffault, 44340 Bouguenais, France. Tel.: +33 02 40 37 16 56.

E-mail address: anne-sophie.lectez@ec-nantes.fr (A.-S. Lectez).</footnote>
<figure><img src="image_2.png" /><figcaption>Fig. 1. Notations used in the Zener model and the Poynting–Thomson model for large strain formulation.</figcaption></figure>

of these models requires to consider restrictions of material functions as order of expansion or decoupling between time and strain influences on relaxation functions. Sometimes, these hypotheses lead to a loss of dependency of out-of-equilibrium stress on preload (see Lion et al. [17] for example). This can be corrected by adding a dependency through correction factors [21] or preload dependent moduli [18]. Thus, thanks to their general and flexible forms, it has been proven by several authors that linearized integral models provide an efficient way to explicitly introduce dependency of equilibrium and out-of-equilibrium stress on preload: a multitude of functions of various order and coupling is available to be compared with experimental data. However, this strategy has to be employed carefully, since too flexible models with an abusive large number of material parameters would directly lead to non-robust models and irrelevant mechanical predictions out of the range of deformations considered during the identification of material functions. Moreover, integral models are generally not based on the second principle of thermodynamics (see Petiteau et al. [22]). In contrary to that, internal variables models are derived from the second principle. Several authors have exploited this advantage in the context of linearization around a large preload. For example, Reese and Govindjee [23], Haupt and Sedlan [24], Amin et al. [25], Lejeunes et al. [26] studied linearized Zener-like models. In the present paper, it is chosen to focus on linearization of internal variable models. Particularly, it will be shown that linearized Zener-like models, in their simpler form, do not take into account a preload dependency on the out-of-equilibrium stress. As well as integral models, it is possible to modify the model introducing such a dependency, but to the authors opinion, this solution is more artificial than using a model which already include the preload dependency, even in its simplest form. The Poynting–Thomson model presents this property and yet is very few considered in linearization processes. Thus, both Zener and Poynting–Thomson models are considered: while their small strain formulation are equivalent, it has been shown by Huber and Tsakmakis [27] that they differ for large strain. In this paper, the comparison between both models is carried further, studying the difference due to large strain after linearization. First, large strain formulations are reminded, before introducing the linearized ones. Then, an analysis of how both models account for modifications of dynamic properties due to a preload is performed. Finally their differences, advantages and drawbacks are discussed in the last section.

## 2. Linearization of large strain formulations

Based on Huber and Tsakmakis [27], large strain formulations of Zener and Poynting–Thomson models are briefly recalled, before being linearized.

### 2.1. Large strain

In addition to the framework of incompressible hyperelasticity, both models are based on three major hypotheses:

* The multiplicative decomposition of the deformation gradient. This assumption is related to the choice of internal variables which describe the problem and define the intermediate configuration.

* The additive decomposition of the strain energy density into two strain energy densities. On the rheological schemes, one strain energy density is associated with each spring: it defines the nonlinear stress–strain relationship.

* The choice of an evolution law which respects the second principle of thermodynamics. Here, a simple case is considered with constant viscosity for both models.

From these three hypotheses follow both large strain formulations, whose major components are introduced in Fig. 1.

#### 2.1.1. Constitutive equations: Zener model

The stress–strain relationship of the large strain formulation of the Zener model based on the multiplicative split of the deformation gradient is described by the following four equations:

$$ 
\left\{ 
\begin{aligned}
\boldsymbol{\sigma} &= \boldsymbol{\sigma}_{\infty} + \boldsymbol{\sigma}_{OE} \\
\boldsymbol{\sigma}_{\infty} &= -p_{\infty}\boldsymbol{I} + 2\boldsymbol{B}_{e}\frac{\partial W_{a}}{\partial\boldsymbol{B}_{e}} \\
\boldsymbol{\sigma}_{OE} &= -(p - p_{\infty})\boldsymbol{I} + 2\boldsymbol{B}_{e}\frac{\partial W_{b}}{\partial\boldsymbol{B}_{e}} \\
\boldsymbol{B}_{e} &= -\frac{2}{\eta}\boldsymbol{B}_{e}(\boldsymbol{\sigma}_{OE})^{D} + \boldsymbol{B}_{e}\boldsymbol{L}^{T} + \boldsymbol{L}\boldsymbol{B}_{e}
\end{aligned}
\right. 
\qquad (1) 
$$

where:

* $\boldsymbol{\sigma}$, $\boldsymbol{\sigma}_{\infty}$ and $\boldsymbol{\sigma}_{OE}$ are respectively the total Cauchy stress tensor, the equilibrium stress tensor and the out-of-equilibrium stress tensor,

* $\boldsymbol{B}$ and $\boldsymbol{L}$ are respectively the total left Cauchy–Green strain tensor and spatial velocity gradient, and $\boldsymbol{B}_{e}$ is the left Cauchy–Green strain tensor associated with spring b,

* $W_a$ and $W_b$ are the strain energy densities which define the nonlinear stiffness of respectively spring a and spring b (hyper-elasticity), $\eta$ is the viscosity constant,

* $p$ and $p_{\infty}$ are hydrostatic pressures related to the incompressibility constraint.

#### 2.1.2. Constitutive equations: Poynting–Thomson model

The stress–strain relationship of the large strain formulation of the Poynting–Thomson model based on the multiplicative split of the deformation gradient is described by the following three
equations:

$$
\left\{
\begin{array}{l}
\boldsymbol{\sigma} = -p\boldsymbol{I} + 2\boldsymbol{B}_e^* \frac{\partial W_1}{\partial \boldsymbol{B}_e^*} \\[1em]
\hat{\boldsymbol{S}}_S = -p_s\boldsymbol{I} + 2\boldsymbol{B}_i^* \frac{\partial W_2}{\partial \boldsymbol{B}_i^*} \\[1em]
\hat{\boldsymbol{B}}_e^* = -\frac{2}{\eta^*} (\boldsymbol{B}_e^* \boldsymbol{\sigma}^D - \hat{\boldsymbol{F}}_e^* \hat{\boldsymbol{S}}_S^D \hat{\boldsymbol{F}}_e^{*T}) + \boldsymbol{B}_e^* \boldsymbol{L}^T + \boldsymbol{L}\boldsymbol{B}_e^*
\end{array}
\right.
\qquad (2)
$$

where:

* $\hat{S}_S$ is the static stress tensor defined in the intermediate configuration, which is the stress in the top branch,

* $\mathbf{B}_e^*$ and $\mathbf{B}_i^*$ are the left Cauchy-Green strain tensor associated respectively with spring 1 and spring 2 (also with the dashpot),

* $W_1$ and $W_2$ are the strain energy densities which define the nonlinear stiffness of respectively spring 1 and spring 2 (hyperelasticity), $\eta^*$ is the viscosity constant,

* $p$ and $p_s$ are hydrostatic pressures related to the incompressibility constraint.

## 2.2. Linearization

The large strain formulations of both models are then linearized. The linearization of the Zener model is classic. For example, Lion [28], Haupt et al. [29], Lejeunes et al. [26] proposed linearization of this model, more or less enriched. The linearization of the Poynting-Thomson model is more original.

### 2.2.1. Method

A static deformation is considered, characterized by its deformation gradient $\mathbf{F}_0$. A displacement increment $\boldsymbol{\Delta}\mathbf{u}$ is superimposed to this initial deformation, and the total deformation gradient is noted $\mathbf{F}$. The directional derivative of $\mathbf{F}$ is written [30, p. 375]:

$$ \boldsymbol{\Delta}\mathbf{F} = D_{\boldsymbol{\Delta}\mathbf{u}}\mathbf{F} = \mathbf{Grad}(\boldsymbol{\Delta}\mathbf{u}), \qquad (3) $$

where $\mathbf{Grad}$ is the gradient operator defined with respect to the undeformed configuration. With the differential composition, it leads to:

$$ \mathbf{Grad}(\boldsymbol{\Delta}\mathbf{u}) = \mathbf{grad}(\boldsymbol{\Delta}\mathbf{u})\mathbf{F}_0, \qquad (4) $$

where $\mathbf{grad}$ is the gradient operator defined with respect to the preloaded configuration. In what follows, $\mathbf{grad}(\boldsymbol{\Delta}\mathbf{u})$ will be noted $\mathbf{h}$. Finally, the total deformation gradient is:

$$ \mathbf{F} = \mathbf{F}_1\mathbf{F}_0, \qquad (5) $$

$$ \mathbf{F} \approx (\mathbf{I} + \mathbf{h})\mathbf{F}_0, \qquad (6) $$

$$ \mathbf{F} \approx \mathbf{F}_0 + \boldsymbol{\Delta}\mathbf{F}, \qquad (7) $$

where $\mathbf{F}_1$ is the deformation gradient due to small displacements $\boldsymbol{\Delta}\mathbf{u}$. This decomposition leads to the linearization of all strain quantities involved in Eqs. (1) and (2). Details on linearization about large strain can be found in Lianis [14], Casey [31] and Holzapfel [30]. For the sake of brevity, only linearization of the Cauchy stress tensor is considered here. The total Cauchy stress tensor is written as:

$$ \boldsymbol{\sigma} = -p\boldsymbol{I} + \mathbf{F} \frac{\partial W}{\partial \boldsymbol{E}} \mathbf{F}^T, \qquad (8) $$

with $\boldsymbol{E} = \frac{\boldsymbol{F}^T\boldsymbol{F} - \boldsymbol{I}}{2}$ being the Green-Lagrange strain tensor. Linearizing this expression gives:

$$ \boldsymbol{\Delta}\boldsymbol{\sigma} = -\boldsymbol{\Delta}p\boldsymbol{I} + 2\,2h\boldsymbol{B}_0 \left. \frac{\partial W}{\partial \boldsymbol{B}} \right|_{\boldsymbol{B}_0}^S + \boldsymbol{F}_0 \boldsymbol{\Delta} \left. \left( \frac{\partial W}{\partial \boldsymbol{E}} \right) \right|_{\boldsymbol{E}}^T \boldsymbol{F}_0^T, \qquad (9) $$

where $\bullet^S = \bullet^{+} \bullet^{-}$ stands for the symmetric operator. Linearization of

the strain energy derivative gives:

$$ \boldsymbol{\Delta} \left( \frac{\partial W}{\partial \boldsymbol{E}} \right) \bigg|_{\boldsymbol{E}_0} = \left. \frac{\partial^2 W}{\partial \boldsymbol{E}^2} \right|_{\boldsymbol{E}_0} : \boldsymbol{\Delta}\boldsymbol{E} = \boldsymbol{C}^4 : \mathbf{F}_0^T \boldsymbol{\epsilon} \boldsymbol{F}_0, \qquad (10) $$

$\boldsymbol{C}^4$ being the material tangent tensor and $\boldsymbol{\epsilon} = \frac{1}{2}(\boldsymbol{h} + \boldsymbol{h}^T)$ the infinitesimal strain tensor. A push-forward operation of the previous equation gives:

$$ \mathbf{F}_0 \boldsymbol{\Delta} \left( \frac{\partial W}{\partial \boldsymbol{E}} \right) \bigg|_{\boldsymbol{E}_0} \mathbf{F}_0^T = \boldsymbol{c}^4 : \boldsymbol{\epsilon}, \qquad (11) $$

where $\boldsymbol{c}^4$ is the spatial tangent tensor. It is defined as the push-forward operation of $\boldsymbol{c}^4$, and using the relationship $\boldsymbol{\Delta}\boldsymbol{\sigma} = \frac{\partial\boldsymbol{\sigma}}{\partial\boldsymbol{B}}|_{\boldsymbol{B}_0} : \boldsymbol{\Delta}\boldsymbol{B}$, it can be shown that $\boldsymbol{c}^4$ also satisfies the following equation (see for example Holzapfel [30, p. 264]):

$$ \boldsymbol{c}^4 = 4\boldsymbol{B}_0 \left. \frac{\partial^2 W}{\partial \boldsymbol{B}^2} \right|_{\boldsymbol{B}_0} \boldsymbol{B}_0 = \boldsymbol{c}^4(\boldsymbol{B}_0). \qquad (12) $$

Finally, the expression of the total Cauchy stress tensor becomes:

$$ \boldsymbol{\sigma} = \boldsymbol{\sigma}_0 + \boldsymbol{\Delta}\boldsymbol{\sigma} = \boldsymbol{\sigma}_0 - \boldsymbol{\Delta}p\boldsymbol{I} + 2 \left( 2h\boldsymbol{B}_0 \left. \frac{\partial W}{\partial \boldsymbol{B}} \right|_{\boldsymbol{B}_0} \right)^S + \boldsymbol{c}^4(\boldsymbol{B}_0) : \boldsymbol{\epsilon}. \qquad (13) $$

In the following, this result is applied to both Zener and Poynting-Thomson models; the subscript $\bullet_0$ always denotes the preloading state.

### 2.2.2. Zener model

The response of the material during a static preload with the Zener model is obtained by applying the condition of zero strain rate in Eqs. (1):

$$
\left\{
\begin{array}{l}
\boldsymbol{\sigma}_0 = \boldsymbol{\sigma}_{\infty 0} + \boldsymbol{\sigma}_{OE0} \\[1em]
\boldsymbol{\sigma}_{\infty 0} = -p_0 \boldsymbol{I} + 2 \boldsymbol{B}_0 \frac{\partial W_a}{\partial \boldsymbol{B}} \bigg|_{\boldsymbol{B}_0} \\[1em]
\boldsymbol{\sigma}_{OE0} = \mathbf{0} \\[1em]
\boldsymbol{B}_{E0} = \mathbf{I}
\end{array}
\right.
\qquad (14)
$$

Then, from Eqs. (7) and (13) follows the total response of the material to a static preload and a small increment of displacement:

$$
\left\{
\begin{array}{l}
\boldsymbol{\sigma} = \boldsymbol{\sigma}_0 + \boldsymbol{\Delta}\boldsymbol{\sigma}_{\infty} + \boldsymbol{\Delta}\boldsymbol{\sigma}_{OE} \\[1em]
\boldsymbol{\Delta}\boldsymbol{\sigma}_{\infty} = -\boldsymbol{\Delta}p_{\infty}\boldsymbol{I} + \boldsymbol{c}_a^4(\boldsymbol{B}_0) : \boldsymbol{\epsilon} + 2 \, 2h\boldsymbol{B}_0 \left. \frac{\partial W_a}{\partial \boldsymbol{B}} \right|_{\boldsymbol{B}_0}^S \\[1em]
\boldsymbol{\Delta}\boldsymbol{\sigma}_{OE} = -(\boldsymbol{\Delta}p - \boldsymbol{\Delta}p_{\infty})\boldsymbol{I} + \boldsymbol{c}_b^4(\boldsymbol{I}) : \boldsymbol{\epsilon}_e \\[1em]
\dot{\boldsymbol{\epsilon}}_e = -\frac{1}{\eta} \boldsymbol{\Delta}\boldsymbol{\sigma}_{OE}^D + \dot{\boldsymbol{\epsilon}}
\end{array}
\right.
\qquad (15)
$$

where:

* $\boldsymbol{\epsilon}_e = \frac{1}{2}(\boldsymbol{h}_e + \boldsymbol{h}_e^T)$, $\boldsymbol{h}_e$ being the small displacement gradient introduced by the linearization of the multiplicative decomposition of the deformation gradient,

* $\boldsymbol{c}_a$ and $\boldsymbol{c}_b$ are the spatial tangent tensors associated respectively with strain energy densities $W_a$ and $W_b$; they represent the nonlinear stiffness of each spring of the model around a given position, respectively $\boldsymbol{B}_0$ (deformed position) and $\mathbf{I}$ (undeformed position).

### 2.2.3. Poynting-Thomson model

The response of the material during a static preload with the Poynting-Thomson model is obtained by applying the condition of
<figure><img src="image_3.png" /><figcaption>Fig. 2. Comparison of the state of each elements in rheological schemes of Zener and Poynting-Thomson models after a static deformation.</figcaption></figure>

zero strain rate in Eq. (2):

$$ \left\{ \begin{array}{l} \sigma_0 = -p_0^* \mathbf{I} + 2 \mathbf{B}_{e0}^* \frac{\partial W_1}{\partial \mathbf{B}_e^*} \bigg|_{\mathbf{B}_{e0}} \\ \hat{\mathbf{S}}_{S0} = -p_{S0} \mathbf{I} + 2 \mathbf{B}_{i0}^* \frac{\partial W_2}{\partial \mathbf{B}_i^*} \bigg|_{\mathbf{B}_{i0}} \\ \mathbf{B}_{e0}^* \sigma_0^D = \mathbf{F}_{e0}^* \hat{\mathbf{S}}_{S0}^D \mathbf{F}_{e0}^{*T} \end{array} \right. \qquad (16) $$

In contrary to the Zener model, none of the internal variables are zero in the case of a static preload. The static condition is expressed differently, through the second equation of the system: in the static case, the stresses in both springs are equal. As it will be seen in the following, having non-null internal variables describing the preload is a good start for obtaining a natural influence of the preload on all stress components. However, less simplification occurs, since the only terms which disappear are those which contain $\mathbf{L}$ or $\hat{\mathbf{B}}_e^*$.

Then, from Eqs. (7) and (13) follows the total response of the material to a static preload and a small increment of displacement:

$$ \left\{ \begin{array}{l} \sigma = \sigma_0 + \Delta\sigma \\ \Delta\sigma = -\Delta p \mathbf{I} + c_1^4 (\mathbf{B}_{e0}^*) : \epsilon_e^* + 2 \mathbf{h}_e^* \mathbf{B}_{e0}^* \frac{\partial W_1}{\partial \mathbf{B}_e^*} \bigg|_{\mathbf{B}_{e0}} \bigg|^S \\ \Delta \hat{\mathbf{S}}_S = -\Delta p \mathbf{I} + c_2^4 (\mathbf{B}_{i0}^*) : \epsilon_i^* + 2 \mathbf{h}_i^* \mathbf{B}_{i0}^* \frac{\partial W_2}{\partial \mathbf{B}_i^*} \bigg|_{\mathbf{B}_{i0}} \bigg|^S \\ (\hat{\mathbf{h}}_e^* \mathbf{B}_{e0}^*)^S = -\frac{2}{\eta^*} \left[ (\hat{\mathbf{h}}_e^* \mathbf{B}_{e0}^*)^S \sigma_0^D + \frac{1}{2} \mathbf{B}_{e0}^* \Delta\sigma^D - (\hat{\mathbf{h}}_e^* \mathbf{F}_{e0}^* \mathbf{S}_{S0}^D \mathbf{F}_{e0}^{*T})^S - \mathbf{F}_{e0}^* \Delta \hat{\mathbf{S}}_S^D \mathbf{F}_{e0}^{*T} \right] + (\hat{\mathbf{h}} \mathbf{B}_{e0}^*)^S \end{array} \right. \qquad (17) $$

In this system:

* $\epsilon_e^* = \frac{1}{2}(\hat{\mathbf{h}}_e^* + \hat{\mathbf{h}}_e^{*T})$, $\epsilon_i^* = \frac{1}{2}(\hat{\mathbf{h}}_i^* + \hat{\mathbf{h}}_i^{*T})$, $\hat{\mathbf{h}}_e^*$ and $\hat{\mathbf{h}}_i^*$ being the small displacement gradients introduced by the linearization of the multiplicative decomposition of the deformation gradient,

* $c_1$ and $c_2$ are the spatial tangent tensors associated respectively with strain energy densities $W_1$ and $W_2$: they represent the nonlinear stiffness of each spring of the model around a given position, respectively $\mathbf{B}_{e0}^*$ (deformed position) and $\mathbf{B}_{i0}^*$, which both represent a deformed position.

The linearization of the Zener model around a static preload has been recalled, and one of the Poynting-Thomson models has been derived. Both systems are now compared in detail, focusing on the influence of the preload.

### 3. Analysis of the influence of a large static preload

#### 3.1. Comparison of linearized formulations

Observing the system of Eqs. (15) for the Zener model, the linearized equilibrium stress tensor depends on the preload. Obviously, $c_a^4$ and $\frac{\partial W_a}{\partial \mathbf{B}}|_{\mathbf{B}_0}$ are expressed in the preloaded state. In contrary to that, none of the terms of the linearized out-of-equilibrium stress tensor depend on the preload: $c_b^4$ is expressed in the undeformed state (the left Cauchy-Green strain tensor is identity).

One spring is not activated during a static preload for the Zener model, because one part of the multiplicative decomposition of the deformation gradient is equal to the identity tensor. It is due to the definition of the intermediate configuration of this model, and leads to the independency of the linearized overstress on the preload. In the case of Poynting-Thomson model, prescribing an infinitely slow preload does not lead to one of the part of the multiplicative decomposition being equal to the identity tensor: both springs are deformed. The static condition imposes another type of constraint, which is the equality of stresses in both springs after the preload (see Eq. (16)). As it can be seen in Eqs. (17), all terms of the linearized total and static stresses are expressed with respect to a deformed state. This difference of behavior is summarized in Fig. 2. For the Zener model, the only spring which response depends on strain rate is associated in series with a dashpot and is thus undeformed in the static case. Thus, the overstress during oscillations around a static preload is linearized around an undeformed state for the Zener model. Concerning the Poynting-Thomson model, the response of both springs depends on strain rate and both springs are deformed, even in the static case. Thus, each linearized stress of this model depends on the strain rate as well as on the preload.

#### 3.2. Equivalence for a simple case

At this point, it has been clearly shown that:

* the equilibrium linearized stress depends on preload for the Zener model,

* the out-of-equilibrium linearized stress does not depend on preload for the Zener model (in a natural way),

* each stress component of the Poynting-Thomson model depends on preload.

Thus, we cannot establish any rigorous comparison on the effect of the preload on both models as long as the equilibrium and out-of-
<figure><img src="image_4.png" /><figcaption>Fig. 3. Evolution of equilibrium stiffness with respect to extension of a static preload, predicted by Zener and Poynting–Thomson models with C₁=0.3 MPa, C₂=0.75 MPa and Cₐ=0.7 MPa.</figcaption></figure>

equilibrium parts of the linearized stress have not been expressed for the Poynting–Thomson model, and the influence of the preload on these precise parts not studied. As the equilibrium/out-of-equilibrium decomposition is inherent to Zener model, it means that a strict equivalence has to be done between both models. For small strains without any preload, the equivalence is straightforward (see for example Huber and Tsakmakis [27]). In the case of small strains superimposed on a preload, the equivalence becomes more complex to establish.

We choose here to study a simple theoretical case of a large static preload followed by small oscillations, both in uniaxial tension–compression loading conditions. The large preload is characterized by its extension λ₀. The values considered here for the preload vary between a compressive strain of 50% (i.e. λ₀ = 0.5) to a tensile strain of 300% (i.e. λ₀ = 4). In the case of the Poynting–Thomson model, the extension is decomposed following the multiplicative decomposition of the deformation gradient:

$$\lambda_0 = \lambda_{e0}^* \lambda_{i0}^* \qquad (18)$$

A uniaxial small strain ε is prescribed during oscillations, and is also decomposed:

$$\epsilon = \epsilon_e^* + \epsilon_i^* \qquad (19)$$

Moreover, the simplest strain energies are chosen. For the Poynting–Thomson model, neo-Hookean strain energies are selected, characterized by parameters C₁, associated with spring 1, and C₂, associated with spring 2. For the Zener model, a neo-Hookean strain energy is also selected for spring a, characterized by parameter Cₐ. Choosing this type of law for spring b would lead to zero out-of-equilibrium response, that is why a Mooney–Rivlin type is preferred:

$$W_b = C_{b1}(I_1 - 3) + C_{b2}(I_2 - 3), \qquad (20)$$

I₁ and I₂ being the two first invariants of Cauchy–Green strain tensors, Cₐ₁ and Cₐ₂ are parameters. First, the equilibrium stresses due to the static preload and to superimposed cyclic oscillations of the Poynting–Thomson model are compared with responses of the Zener model. Then, the out-of-equilibrium linearized stress is calculated, and it is shown that it depends on the preload, while Zener one does not.

### 3.2.1. Equilibrium stress

(a) Preload: For an infinitely slow motion, Eqs. (16) hold for the Poynting–Thomson model. Projecting in the tension direction and replacing the hydrostatic pressure with a free lateral surface

condition, the system becomes:

$$\left\{ \begin{array}{l} \sigma_0 = 2C_1 \lambda_{e0}^{*2} - \frac{1}{\lambda_{e0}^{*}} \\ \frac{4}{3} C_1 \lambda_{e0}^{*2} \lambda_{e0}^{*2} - \frac{1}{\lambda_{e0}^{*}} = \frac{4}{3} C_2 \lambda_{e0}^{*2} \lambda_{i0}^{*2} - \frac{1}{\lambda_{i0}^{*}} \end{array} \right. \qquad (21)$$

where σ₀ is the only non-null component of the Cauchy stress tensor of the large static preload. Taking into account that λᵢ₀* = λ₀/λᵢ₀*, the evolution law becomes:

$$\frac{C_1}{C_2} \lambda_{e0}^{*2} - \frac{1}{\lambda_{e0}^{*}} = \left( \frac{\lambda_{e0}^{*2}}{\lambda_{e0}^{*2}} - \frac{\lambda_{e0}^{*}}{\lambda_{0}} \right) \qquad (22)$$

To solve Eqs. (21), λᵢ₀* has to be expressed as a function of λ₀. This is equivalent to finding a real positive root to the polynomial:

$$P(\lambda_{e0}^*) = G\lambda_{e0}^{*4} + \frac{\lambda_{e0}^{*3}}{\lambda_0} - G\lambda_{e0}^* - \lambda_0^2 \qquad (23)$$

where G = C₁/C₂. This problem can be solved using Ferrari and Cardan methods: λᵢ₀* can be found for each value of G and λ₀.

(b) Oscillations: In order to obtain the equilibrium stress from the Poynting–Thomson model, static oscillations superimposed on a static preload are considered. The equilibrium linearized stress is found by linearizing system of Eqs. (21):

$$\left\{ \begin{array}{l} \Delta\sigma_{\infty} = 2C_1\epsilon_e^* \left( 2\lambda_{e0}^{*2} + \frac{1}{\lambda_{e0}^{*}} \right) \\ \epsilon_e^* = \frac{C_2 \left( 2\lambda_{i0}^{*2} + \frac{1}{\lambda_{i0}^{*}} \right)}{C_1 \left( 2\lambda_{e0}^{*2} + \frac{1}{\lambda_{e0}^{*}} \right) + C_2 \left( 2\lambda_{i0}^{*2} + \frac{1}{\lambda_{i0}^{*}} \right)} \epsilon \end{array} \right. \qquad (24)$$

First, εᵢ* is calculated through the linearized evolution law. It depends on λᵢ₀* and λᵢ₀*, and thus on λ₀, which are known from the previous paragraph. Knowing εᵢ* and λᵢ₀*, Δσᵢ* is determined. The response predicted by the Zener model for uniaxial static oscillations follows from Eq. (15)₂ and is given by:

$$\Delta\sigma_{\infty} = 2C_a\epsilon\left(2\lambda_0^2 + \frac{1}{\lambda_0}\right) \qquad (25)$$

Fig. 3 shows the evolution of equilibrium stiffness (taken as the ratio of Δσᵢ* on ε) predicted by both models with respect to the preloading extension. As expected, both Zener and Poynting–Thomson models include an influence of the preload on the equilibrium part of the linearized stress. Next step consists in studying the out-of-equilibrium part.

### 3.2.2. Out-of-equilibrium stress

The out-of-equilibrium stress is simply the difference between the total and equilibrium stresses. Preload being considered as static, the linearized out-of-equilibrium stress predicted by the Poynting–Thomson model is thus:

$$\Delta\sigma_{oe} = \Delta\sigma - \Delta\sigma_{\infty} = 2C_1 \left( 2\lambda_{e0}^{*2} + \frac{1}{\lambda_{e0}^{*}} \right) \left( \epsilon_e^* - \epsilon_{e\infty}^* \right), \qquad (26)$$

where εᵢ* is calculated with Eq. (24)₂. It is chosen to consider the simple case of infinitely fast oscillations around a static preload. Then, εᵢ* = ε for the Poynting–Thomson model. Replacing εᵢ* by ε and εᵢ* by Eqs. (24)₂, Eq. (26) becomes:

$$\Delta\sigma_{oe} = M_0\epsilon, \qquad (27)$$
<figure><img src="image_5.png" /><figcaption>Fig. 4. Evolution of out-of-equilibrium stiffness during infinitely fast oscillations with respect to extension of a static preload, predicted by Zener and Poynting-Thomson models with C₁=0.3 MPa, C₂=0.75 MPa and C₉₂=2 MPa.</figcaption></figure>

where:

$$M_0 = \frac{2C_1^2 2\lambda_{e0}^{*2} + \frac{1}{\lambda_{e0}^*}}{C_1 2\lambda_{e0}^{*2} + \frac{1}{\lambda_{e0}^*} + C_2 \left(2\lambda_{i0}^{*2} + \frac{1}{\lambda_{i0}^*}\right)} \quad (28)$$

depends on the static preload extension.

For the Zener model, Eq. (15)₃ applied to the case of infinitely fast oscillations (e = e₀) becomes:

$$\Delta\sigma_{oe} = \frac{C_{b2}}{2}\epsilon. \quad (29)$$

Fig. 4 shows the evolution of the out-of-equilibrium stiffness with respect to the preload extension for both models, i.e. M₀ for the Poynting-Thomson model and C₉₂/₂ for the Zener model. Clearly, the Zener model is unable to predict any dependency of the out-of-equilibrium stiffness on the static preload, while the Poynting-Thomson model does, even with less parameters.

# 4. Discussion

The objective here is to answer the question: Which is the better viscoelastic approach in the context of a linearization around a static preload? For this, several elements have to be considered. The first one is of course the influence of the preload.

## 4.1. Influence of the preload

In the previous section, it has been shown that both models induce an influence of the preload on the linearized equilibrium stress. However, only the Poynting-Thomson model induces an influence of the preload on the linearized overstress. For Zener-like models, some authors choose to modify the evolution law and introduce a strain-dependent viscosity. For example, Holzapfel [32], Reese and Govindjee [23], Miehe and Keck [33], Haupt and Sedlan [24], Amin et al. [25] adopted this strategy. This approach can be considered as artificial, since it supposes to modify the model. The question then stands to know which dependency of the viscosity on the preload has to be chosen, which can be determined with experimental tests. With Poynting-Thomson-like models, this question does not arise since the influence of the preload on the linearized overstress already exists with the simplest version of the model. Of course, it will also be needed to check if this influence is in agreement with experimental data, but it is clearly an advantage to naturally include this influence of preload.

## 4.2. Identification procedure

To be easier, the identification procedure has to be based on experimental loading conditions identical to those under which is studied the material behavior, that is to say small oscillations around a static preload. Thus, the prediction of both models has first to be accurate in the case of static loading conditions.

However, the best way to identify parameters of models is to decompose the identification with several loading conditions, each of them corresponding to a particular set of parameters. For example, for the Zener model, in order to isolate the parameters of the top spring a, static loading tests have to be performed. So, the parameters of the top spring can be identified with the preload data, and the parameters of bottom branch (spring b and dashpot viscosity) are identified with data collected during oscillating loading conditions. The procedure is more complicated with the Poynting-Thomson model. Indeed, in order to isolate one of the springs and be able to identify only its parameters with particular loading conditions, it would be necessary to consider infinitely fast loading conditions. Performing high speed tests is possible but not easy to perform, as it can be seen in Aloui et al. [34] and Petiteau et al. [35]. Moreover, it supposes to add some tests to the already existing experimental campaign. Thus, one solution is to identify all parameters simultaneously on all tests data, i.e. both preload and oscillations. This procedure is simple since other tests are not needed, but less robust than identifying parameters from adapted tests.

## 4.3. Finite element simulation

Industrial applications which involve elastomers subjected to large static preloads followed by small oscillations often concern geometrically complex parts. Thus, the prediction of the behavior often requires finite element simulations, possibly with commercial codes. It can be quite easy to use Zener model in finite element simulations, since it can be either yet implemented, or it is easy to find an equivalent, for example with Prony series models. In contrary to that, a user-defined implementation would be required to use Poynting-Thomson-like models.

## 4.4. Advantages and drawbacks

Finally, choosing between both models is not obvious. The Poynting-Thomson model has indeed the large advantage to include naturally an influence of the preload on both linearized equilibrium and out-of-equilibrium parts of the stress. This is due to its intermediate configuration which is not an equilibrium one. However, what makes a big advantage induces also a major drawback, since the identification procedure based on preload and oscillation tests is not well adapted. Finally, in an industrial application, it will be more difficult to use, since its use in finite element commercial codes requires to implement it.

In contrary, the Zener model already exists or has direct equivalents in finite element commercial codes, and the identification procedure is well established, since the preload data correspond only to the parameters of one spring. Another advantage of the Zener model is the easiness to be generalized: it is well known that adding branches to rheological models, which means increasing the number of internal variables, allows to improve their accuracy on a wider relaxation spectrum (see for example Kaliske and Rothert [36] or Haupt et al. [29]). For the Zener model, the generalization results in N Maxwell elements are in parallel. Then, N independent evolutions laws have to be solved to determine the N multiplicative decompositions of the deformation gradient. The Poynting-Thomson model is more complex to generalize, since it results in N Kelvin elements in series, and thus N
coupled evolution laws. However, the Zener model does not have an influence of the preload on the out-of-equilibrium part of the linearized stress, unless an artificial influence is introduced through a particular form of a strain-dependent viscosity that has to be postulated.

## 5. Conclusion

Some challenging industrial applications require to master modeling the behavior of elastomers under large static preloads on which are superimposed small oscillations. In this paper, attention is focused on linearization of internal variables models. Particularly, both well-known Zener and Poynting-Thomson are considered. While the linearization of Zener model is classic, the linearization of the Poynting-Thomson model has, to the authors knowledge, never been presented. It has been shown that the Poynting-Thomson model naturally produces an influence of preload on both parts of the linearized stress, while until now this influence has been forced through a modified viscosity in Zener-like models. However, this characteristic is not the only one to take into account in the comparison of both models. Thus, it has been highlighted that the identification procedure as well as the implementation are easier with the Zener model. Finally, the choice depends on context, time and resources: the Poynting-Thomson model is more interesting but requires time for its implementation and probably more parameters than the Zener model to accurately predict the behavior during preload, while the Zener model is simpler to use in finite element simulations; nevertheless if the influence of the preload on the overstress has to be taken into account, the introduction of a strain dependent viscosity remains an open question. Of course, it is essential to compare the influence of the preload predicted by both models with experimental data to pronounce in favor of one or the other model. This subject is currently still investigated, but to the authors opinion, the theoretical reasoning and possibilities exposed in this paper justify not to limit studies on linearized viscoelasticity around large preload to Zener-like models.

## Acknowledgments

The authors are grateful to Renault S.A.S. for having financed this project.

## References

[1] S. Govindjee, J.C. Simo, Mullins effect and the strain amplitude dependence of the storage modulus, Int. J. Solids Struct. 29 (14) (1992) 1737-1751.

[2] I.M. Ward, Mechanical Properties of Solid Polymers, 2nd edn., John Wiley and Sons Ltd., Chichester/New York, 1983.

[3] A.E. Green, R.S. Rivlin, The mechanics of non-linear materials with memory, Arch. Ration. Mech. Anal. 1 (1) (1957) 1-21.

[4] B.D. Coleman, W. Noll, Foundations of linear viscoelasticity, Rev. Mod. Phys. 33 (2) (1961) 239-249.

[5] R.M. Christensen, A nonlinear theory of viscoelasticity for application to elastomers, J. Appl. Mech. 47 (1980) 762-768.

[6] B. Bernstein, E.A. Kearsley, L.J. Zapas, A study of stress relaxation with finite strain, Trans. Soc. Rheol. 7 (1) (1963) 391-410.

[7] A. Kaye, Non-Newtonian Flow in Incompressible Fluids, College of Aeronautics, Cranford, UK, Note No. 134, 1962.

[8] A.S. Wineman, K.R. Rajagopal, Mechanical Response of Polymers, An Introduction, Cambridge University Press, Cambridge, 2000.

[9] M.S. Green, A.V. Tobolsky, A new approach for the theory of relaxing polymeric media, J. Chem. Phys. 14 (1946) 87-112.

[10] J. Lubliner, A model of rubber viscoelasticity, Mech. Res. Commun. 12 (1985) 93-99.

[11] F. Sidoroff, Un modèle viscoélastique non linéaire avec configuration intermédiaire, J. Méc. 13 (1974) 679-713.

[12] E.H. Lee, Elastic-plastic deformation at finite strains, J. Appl. Mech. 36 (1) (1969) 1-6.

[13] A.C. Pipkin, Small finite deformations of viscoelastic solids, Rev. Mod. Phys. 36 (1964) 1034-1041.

[14] G. Lianis, Application of irreversible thermodynamics in finite viscoelastic deformations, J. Appl. Mech. 32 (1965) 623-629.

[15] W. Goldberg, G. Lianis, Behavior of viscoelastic media under small sinusoidal oscillations superposed on finite strain, J. Appl. Mech. 35 (3) (1968) 433-440.

[16] K.N. Morman, J.C. Nagtegaal, Finite element analysis of sinusoidal small-amplitude vibrations in deformed viscoelastic solids. Part I: theoretical development, Int. J. Numer. Methods Eng. 19 (7) (1983) 1079-1103.

[17] A. Lion, J. Retka, M. Rendek, On the calculation of predeformation-dependent dynamic modulus tensors in finite nonlinear viscoelasticity, Mech. Res. Commun. 36 (6) (2009) 653-658.

[18] D. Wollscheid, A. Lion, Predeformation-and-frequency-dependent material behaviour of filler-reinforced rubber: experiments, constitutive modelling and parameter identification, Int. J. Solids Struct. 50 (9) (2013) 1217-1225.

[19] A. Pipkin, T. Rogers, A non-linear integral representation for viscoelastic behaviour, J. Mech. Phys. Solids 16 (1) (1968) 59-72.

[20] A. Wineman, Nonlinear viscoelastic solids—a review, Math. Mech. Solids 14 (3) (2009) 300-366.

[21] B.-K. Kim, S.-K. Youn, A viscoelastic constitutive model of rubber under small oscillatory load superimposed on large static deformation, Arch. Appl. Mech. 71 (11) (2001) 748-763.

[22] J.-C. Petiteau, E. Verron, R. Othman, H. Le Sourne, J.-F. Sigrist, G. Barras, Large strain rate-dependent response of elastomers at different strain rates: convolution integral vs. internal variable formulations, Mech. Time-Depend. Mater. 17 (3) (2013) 349-367.

[23] S. Reese, S. Govindjee, A theory of finite viscoelasticity and numerical aspects, Int. J. Solids Struct. 35 (26) (1998) 3455-3482.

[24] P. Haupt, K. Sedlan, Viscoplasticity of elastomeric materials: experimental facts and constitutive modelling, Arch. Appl. Mech. 71 (2-3) (2001) 89-109.

[25] A. Amin, A. Lion, S. Sekita, Y. Okui, Nonlinear dependence of viscosity in modeling the rate-dependent response of natural and high damping rubbers in compression and shear: experimental identification and numerical verification, Int. J. Plast. 22 (9) (2006) 1610-1657.

[26] S. Lejeunes, T. Nguyen Van, A. Boukamel, D. Eyheramendy, Linearized behavior of a dissipative rubber with large static preloads, Comput. Struct. 96 (2012) 46-53.

[27] N. Huber, C. Tsakmakis, Finite deformation viscoelasticity laws, Mech. Mater. 32 (1) (2000) 1-18.

[28] A. Lion, Thixotropic behaviour of rubber under dynamic loading histories: experiments and theory, J. Mech. Phys. Solids 46 (5) (1998) 895-930.

[29] P. Haupt, A. Lion, E. Backhaus, On the dynamic behaviour of polymers under finite strains: constitutive modelling and identification of parameters, Int. J. Solids Struct. 37 (26) (2000) 3633-3646.

[30] G.A. Holzapfel, Nonlinear Solid Mechanics. A Continuum Approach for Engineering, J. Wiley and Sons, Chichester, 2000.

[31] J. Casey, On infinitesimal deformation measures, J. Elast. 28 (3) (1992) 257-269.

[32] G.A. Holzapfel, On large strain viscoelasticity: continuum formulation and finite element applications to elastomeric structures, Int. J. Numer. Methods Eng. 39 (22) (1996) 3903-3926.

[33] C. Miehe, J. Keck, Superimposed finite elastic-viscoelastic-plastoelastic stress response with damage in filled rubbery polymers. Experiments, modelling and algorithmic implementation, J. Mech. Phys. Solids 48 (2) (2000) 323-365.

[34] S. Aloui, R. Othman, E. Verron, P. Guégan, Semi-analytic inverse method for rubber testing at high strain rates, Mech. Res. Commun. 47 (2013) 97-101.

[35] J.-C. Petiteau, R. Othman, P. Guégan, H.L. Sourne, E. Verron, Dynamic uniaxial extension of elastomers at constant true strain rate, Polym. Test. 32 (2) (2013) 394-401.

[36] M. Kaliske, H. Rothert, Formulation and implementation of three-dimensional viscoelasticity at small and finite strains, Comput. Mech. 19 (3) (1997) 228-239.