Viscoelasticity

The response of a viscoelastic material includes both elastic (instantaneous) and viscous (time-dependent) behavior. The instantaneous elastic response of the material is followed by creep when subjected to a fixed applied stress, and stress-relaxation when subjected to a fixed applied strain.

This page discusses:

Reduced states of stress
Automatic time stepping procedure

See Also

In Other Guides

Time Domain Viscoelasticity

ProductsAbaqus/StandardAbaqus/Explicit

The basic hereditary integral formulation for linear isotropic viscoelasticity is

σ (t) = \int_{0}^{t} 2 G (τ - τ^{'}) \dot{e} d t^{'} + I \int_{0}^{t} K (τ - τ^{'}) \dot{ϕ} d t^{'} .

Here $e$ and $ϕ$ are the mechanical deviatoric and volumetric strains; K is the bulk modulus and G is the shear modulus, which are functions of the reduced time $τ$ ; and $\dot{}$ denotes differentiation with respect to $t^{'}$ .

The reduced time is related to the actual time through the integral differential equation

τ = \int_{0}^{t} \frac{d t^{'}}{A_{θ} (θ (t^{'}))}, \frac{d τ}{d t} = \frac{1}{A_{θ} (θ (t))},

where $θ$ is the temperature and $A_{θ}$ is the shift function. (Hence, if $A_{θ} = 1$ , $τ = t$ .) A commonly used shift function is the Williams-Landel-Ferry (WLF) equation, which has the following form:

- \log A_{θ} = h (θ) = \frac{C_{1}^{g} (θ - θ_{g})}{C_{2}^{g} + (θ - θ_{g})},

where $C_{1}^{g}$ and $C_{2}^{g}$ are constants and $θ_{g}$ is the “glass” transition temperature. This is the temperature at which, in principle, the behavior of the material changes from glassy to rubbery. If $θ \leq θ_{g} - C_{2}^{g}$ , deformation changes will be elastic. $C_{1}^{g}$ and $C_{2}^{g}$ were once thought to be “universal” constants whose values were obtained at $θ_{g}$ , but these constants have been shown to vary slightly from polymer to polymer.

Abaqus allows the WLF equation to be used with any convenient temperature, other than the glass transition temperature, as the reference temperature. The form of the equation remains the same, but the constants are different. Namely,

- \log A_{θ} = h (θ) = \frac{C_{1} (θ - θ_{0})}{C_{2} + (θ - θ_{0})},

where $θ_{0}$ is the reference temperature at which the relaxation data are given, and $C_{1}$ and $C_{2}$ are the calibration constants at the reference temperature. The “universal” constants $C_{1}^{g}$ and $C_{2}^{g}$ are related to $C_{1}$ and $C_{2}$ as follows:

\begin{matrix} C_{1} & ​ = \frac{C_{1}^{g}}{1 + (θ_{0} - θ_{g}) / C_{2}^{g}}, \\ C_{2} & ​ = C_{2}^{g} + θ_{0} - θ_{g} . \end{matrix}

Other forms of $h (θ)$ are also used, such as a power series in $θ - θ_{0}$ . Abaqus allows a general definition of the shift function with user subroutine UTRS.

The relaxation functions $K (t)$ and $G (t)$ can be defined individually in terms of a series of exponentials known as the Prony series:

K (τ) = K_{\infty} + \overset{n_{K}}{\sum_{i = 1}} K_{i} e^{- τ / τ_{i}^{K}} G (τ) = G_{\infty} + \overset{n_{G}}{\sum_{i = 1}} G_{i} e^{- τ / τ_{i}^{G}},

where $K_{\infty}$ and $G_{\infty}$ represent the long-term bulk and shear moduli. In general, the relaxation times $τ_{i}^{K}$ and $τ_{i}^{G}$ need not equal each other; however, Abaqus assumes that $τ_{i} = τ_{i}^{K} = τ_{i}^{G}$ . On the other hand, the number of terms in bulk and shear, $n_{K}$ and $n_{G}$ , need not equal each other. In fact, in many practical cases it can be assumed that $n_{K} = 0$ . Hence, we now concentrate on the deviatoric behavior. The equations for the volumetric terms can be derived in an analogous way.

The deviatoric integral equation is

\begin{matrix} S & ​ = \int_{0}^{t} 2 (G_{\infty} + \sum_{i = 1}^{n_{G}} G_{i} e^{(τ^{'} - τ) / τ_{i}}) \dot{e} d t^{'} \\ ​ = \int_{0}^{τ} 2 (G_{\infty} + \sum_{i = 1}^{n_{G}} G_{i} e^{(τ^{'} - τ) / τ_{i}}) \frac{d e}{d τ^{'}} d τ^{'} . \end{matrix}

We rewrite this equation in the form

\begin{matrix} (1) & S = 2 G_{0} (e - \sum_{i = 1}^{n} α_{i} e_{i}), \end{matrix}

where $G_{0} = G_{\infty} + \sum_{i = 1}^{n} G_{i}$ is the instantaneous shear modulus, $α_{i} = G_{i} / G_{0}$ is the relative modulus of term i, and

\begin{matrix} (2) & e_{i} = \int_{0}^{τ} (1 - e^{(τ^{'} - τ) / τ_{i}}) \frac{d e}{d τ^{'}} d τ^{'} \end{matrix}

is the viscous (creep) strain in each term of the series. For finite element analysis this equation must be integrated over a finite increment of time. To perform this integration, we will assume that during the increment $e$ varies linearly with $τ$ ; hence, $d e / d τ^{'} = Δ e / Δ τ$ . To use this relation, we break up Equation 2 into two parts:

\begin{matrix} e_{i}^{n + 1} & ​ = \int_{0}^{τ^{n}} (1 - e^{(τ^{'} - τ^{n + 1}) / τ_{i}}) \frac{d e}{d τ^{'}} d τ^{'} \\ + \int_{τ^{n}}^{τ^{n + 1}} (1 - e^{(τ^{'} - τ^{n + 1}) / τ_{i}}) \frac{d e}{d τ^{'}} d τ^{'} . \end{matrix}

Now observe that

1 - e^{(τ^{'} - τ^{n + 1}) / τ_{i}} = 1 - e^{- Δ τ / τ_{i}} + e^{- Δ τ / τ_{i}} (1 - e^{(τ^{'} - τ^{n}) / τ_{i}}) .

Use of this expression and the approximation for $d e / d τ^{'}$ during the increment yields

\begin{matrix} e_{i}^{n + 1} = (1 - e^{- Δ τ / τ_{i}}) & \int_{0}^{τ^{n}} \frac{d e}{d τ^{'}} d τ^{'} \\ + e^{- Δ τ / τ_{i}} & \int_{0}^{τ^{n}} (1 - e^{(τ^{'} - τ^{n}) / τ_{i}}) \frac{d e}{d τ^{'}} d τ^{'} \\ + \frac{Δ e}{Δ τ} & \int_{τ^{n}}^{τ^{n + 1}} (1 - e^{(τ^{'} - τ^{n + 1}) / τ_{i}}) d τ^{'} . \end{matrix}

The first and last integrals in this expression are readily evaluated, whereas from Equation 2 follows that the second integral represents the viscous strain in the $i^{th}$ term at the beginning of the increment. Hence, the change in the $i^{th}$ viscous strain is

\begin{matrix} (3) & \begin{matrix} Δ e_{i} & ​ = (1 - e^{- Δ τ / τ_{i}}) e^{n} + (e^{- Δ τ / τ_{i}} - 1) e_{i}^{n} + (Δ τ - τ_{i} (1 - e^{- Δ τ / τ_{i}})) \frac{Δ e}{Δ τ} \\ ​ = \frac{τ_{i}}{Δ τ} (\frac{Δ τ}{τ_{i}} + e^{- Δ τ / τ_{i}} - 1) Δ e + (1 - e^{- Δ τ / τ_{i}}) (e^{n} - e_{i}^{n}) . \end{matrix} \end{matrix}

If $Δ τ / τ_{i}$ approaches zero, this expression can be approximated by

\begin{matrix} (4) & Δ e_{i} = \frac{Δ τ}{τ_{i}} (\frac{1}{2} Δ e + e^{n} - e_{i}^{n}) . \end{matrix}

The last form is used in the computations if $Δ τ / τ_{i} < 10^{- 7}$ .

Hence, in an increment, Equation 3 or Equation 4 is used to calculate the new value of the viscous strains. Equation 1 is then used subsequently to obtain the new value of the stresses.

The tangent modulus is readily derived from these equations by differentiating the deviatoric stress increment, which is

Δ S = 2 G_{0} (Δ e - \overset{n_{G}}{\sum_{i = 1}} α_{i} (e_{i}^{n + 1} - e_{i}^{n}))

with respect to the deviatoric strain increment $Δ e$ . Since the equations are linear, the modulus depends only on the reduced time step:

G^{T} = {\begin{matrix} G_{0} [1 - \sum_{i = 1}^{n} α_{i} \frac{τ_{i}}{Δ τ} (\frac{Δ τ}{τ_{i}} + e^{- Δ τ / τ_{i}} - 1)] & if Δ τ / τ_{i} > 10^{- 7} \\ G_{0} [1 - \sum_{i = 1}^{n} \frac{1}{2} α_{i} \frac{Δ τ}{τ_{i}}] & if Δ τ / τ_{i} < 10^{- 7} \end{matrix}

Finally, one needs a relation between the reduced time increment, $Δ τ$ , and the actual time increment, $Δ t$ . To do this, we observe that $A_{θ}$ varies very nonlinearly with temperature; hence, any direct approximation of $A_{θ}$ is likely to lead to large errors. On the other hand, $h (θ)$ will generally be a smoothly varying function of temperature that is well approximated by a linear function of temperature over an increment. If we further assume that incrementally the temperature $θ$ is a linear function of time t, one finds the relation

h (θ) = - \ln A_{θ} (θ (t)) = a + b t

A_{θ}^{- 1} (θ (t)) = e^{a + b t}

with

\begin{matrix} a & ​ = \frac{1}{Δ t} [t^{n + 1} h (θ^{n}) - t^{n} h (θ^{n + 1})] \\ b & ​ = \frac{1}{Δ t} [h (θ^{n + 1}) - h (θ^{n})] . \end{matrix}

This yields the relation

\begin{matrix} Δ τ & ​ = \int_{t^{n}}^{t^{n + 1}} e^{a + b t} d t \\ ​ = \frac{1}{b} (e^{a + b t^{n + 1}} - e^{a + b t^{n}}) . \end{matrix}

This expression can also be written as

Δ τ = \frac{A_{θ}^{- 1} (θ^{n + 1}) - A_{θ}^{- 1} (θ^{n})}{h (θ^{n + 1}) - h (θ^{n})} Δ t .

Viscoelasticity

Reduced states of stress

Automatic time stepping procedure