Finite-strain viscoelasticity

It is assumed that the instantaneous response of the material follows from the hyperelastic constitutive equations:

for a compressible material and

for an incompressible material. In the above, ${\mathit{\tau }}_0^D$ and ${\mathit{\tau }}_0^H$ are, respectively, the deviatoric and the hydrostatic parts of the instantaneous Kirchhoff stress ${\mathit{\tau }}_0$. $\overline{\mathbf{F}}$ is the “distortion gradient” related to the deformation gradient $\mathbf{F}$ by

where

is the volume change.

Using integration by parts and a variable transformation, the basic hereditary integral formulation for linear isotropic viscoelasticity can be written in the form

or entirely in terms of stress quantities,

where $\tau$ is the reduced time, $\dot{G}\left({\tau }^{\prime }\right)=dG\left({\tau }^{\prime }\right)/d{\tau }^{\prime }$, and $\dot{K}\left({\tau }^{\prime }\right)=dK\left({\tau }^{\prime }\right)/d{\tau }^{\prime }$. $G_0$ and $K_0$ are the instantaneous small-strain shear and bulk moduli, and $G\left(t\right)$ and $K\left(t\right)$ are the time-dependent small-strain shear and bulk relaxation moduli. Recall that the reduced time represents a shift in time with temperature and is related to the actual time through the differential equation

where $\theta$ is the temperature and $A_{\theta }$ is the shift function.

Using the volumetric/deviatoric-split hereditary integral in the reference configuration for large strain (hyperelastic) materials, and then using a standard push-forward operator (see Simo, 1987), one obtains the following set of equations in the current configuration:

where $\mathrm{dev}(\bullet )=(\bullet )-\frac{1}{3}(\bullet ):\mathbf{I}\mathbf{I}$ and ${\overline{\mathbf{F}}}_t\left(t-t^{\prime }\right)$ is the distortional deformation gradient of the state at $t-t^{\prime }$ relative to the state at t . A transformation is performed on the stress relating the state at time $t-t^{\prime }$ to the state at time t.

As in small-strain viscoelasticity, we represent the relaxation moduli in terms of the Prony series

where $g_i$ and $k_i$ are the relative moduli of terms i. Note that $g_{\mathrm{\infty }}+{\sum }_{i=1}^{N_G}{g}_i=k_{\mathrm{\infty }}+{\sum }_{i=1}^{N_K}{k}_i=1$. Abaqus assumes that the relaxation times ${\tau }_i={\tau }_i^K={\tau }_i^G$ are the same so that from here on, we will sum on N terms for both bulk and shear behavior. In reality, the number of nonzero terms in bulk and shear, $N_K$ and $N_G$, need not be equal, unless the instantaneous behavior is based on the hyperfoam model. In the latter case, the two deformation modes are closely related and are then assumed to relax equally and simultaneously.

Substituting Equation 2 and Equation 3 in Equation 1, we obtain

Next, we introduce the internal stresses, associated with each term of the series

These stresses are stored at each material point and are integrated forward in time. We will assume that the solution is known at time t, and we need to construct the solution at time $t+\mathrm{\Delta }t$.

The internal hydrostatic stresses ${\mathit{\tau }}_i^H$ at time $t+\mathrm{\Delta }t$ follow from

With $\overline{\tau }={\tau }^{\prime }-\mathrm{\Delta }\tau$ and $\overline{t}=t^{\prime }-\mathrm{\Delta }t$, it follows that

which yields with Equation 6

To integrate the first integral in Equation 7, we assume that ${\mathit{\tau }}_0^H\left(t-\overline{t}\right)$ varies linearly with the reduced time $\overline{\tau }$ over the increment

Substituting into Equation 7 yields

The integrals are readily evaluated, providing the solution at the end of the increment

or, in a slightly different form

with

Observe that for $\mathrm{\Delta }t=\mathrm{\Delta }\tau =0$, ${\gamma }_i=1$ and ${\alpha }_i={\beta }_i=0$. For $\mathrm{\Delta }t=\mathrm{\Delta }\tau =\mathrm{\infty }$, ${\alpha }_i=1$ and ${\gamma }_i={\beta }_i=0$.

The internal deviatoric stresses ${\mathit{\tau }}_i^D$ at time $t+\mathrm{\Delta }t$ follow from

Observe that

and the inverse of this is

which—when substituted into Equation 10 with $\mathrm{\Delta }\overline{\mathbf{F}}\equiv {\overline{\mathbf{F}}}_t\left(t+\mathrm{\Delta }t\right)$ and $\mathrm{\Delta }{\overline{\mathbf{F}}}^{-1}\equiv {\overline{\mathbf{F}}}_{t+\mathrm{\Delta }t}\left(t\right)$—gives

With $\overline{\tau }={\tau }^{\prime }-\mathrm{\Delta }\tau$ and $\overline{t}=t^{\prime }-\mathrm{\Delta }t$, it follows:

Now introduce the variable

Note that

and

Then we can also introduce

Substitution of Equation 5, Equation 12, and Equation 15 into Equation 11 yields

To integrate the first integral in Equation 16, we assume that ${\widehat{\mathit{\tau }}}_0^D\left(t-\overline{t}\right)$ varies linearly with the reduced time $\overline{\tau }$ over the increment:

which with Equation 14 becomes

Equation 16 and Equation 17 for the deviatoric stress have exactly the same form as Equation 7 and Equation 8 for the hydrostatic stress. Hence, after integration we obtain

with

Equation 13, Equation 15, and Equation 18, thus, provide a straightforward integration scheme.

The total stress at the end of the increment becomes

which with Equation 9 and Equation 18 can also be written as

To solve the system of nonlinear equations generated by the constitutive equations, we need to generate the corotational constitutive rate equations. From Equation 20 it follows

where $\check{\mathit{\tau }}\left(t\right)$ is the corotational (Jaumann) stress rate. The rate form of the above equation is used to compute the Jacobian.

All equations have been worked out in terms of the Kirchhoff stress. However, the implementation in Abaqus uses the Cauchy stress. To transform to Cauchy stress, we use the relations

With $\mathrm{\Delta }J\equiv J\left(t+\mathrm{\Delta }t\right)/J\left(t\right)$, this allows us to write Equation 9, Equation 13, Equation 15, Equation 18, and Equation 19 in the following form:

The virtual work and rate of virtual work equations are written with respect to the current volume. Therefore, the corotational stress rates are rates of Kirchhoff stress mapped into the current configuration and transformed in the same way as the stresses themselves.

This set of equations—combined with the expressions for ${\alpha }_i$, ${\beta }_i$, and ${\gamma }_i$—describe the full implementation of the hyper-viscoelasticity model in a displacement formulation.

The rate equations can be written in a form similar to Hyperelastic material behavior. Introduce

and

where ${\mathbf{C}}_0^S$ and $K_0$ are the instantaneous moduli, corresponding to ${\mathbf{C}}^S$ and $K$ of Hyperelastic material behavior. Thus, all rate equations can be obtained by substitution of ${\mathbf{C}}_0^S$ by ${\mathbf{C}}_v^S$ and $K_0$ by $K_v$.