Time-Stepping Algorithm

In hypoelasticity, we required the use of the so-called “Jaumann rate of Cauchy Stress”, $\mathring{\boldsymbol{\sigma}}$ , because $\dot{\boldsymbol{\sigma}}$ is not work conjugate with $\boldsymbol{\sigma}$ or $\mathbf{D}$ , while $\mathring{\boldsymbol{\sigma}}$ is. However, we still need $\dot{\boldsymbol{\sigma}}$ , because $\boldsymbol{\sigma}_n=\boldsymbol{\sigma}_{n+1}+\dot{\boldsymbol{\sigma}}\Delta t$ is the stress tensor that we need. The Cauchy stress (and its rate) have a well-understood physical meaning (whereas $\mathring{\boldsymbol{\sigma}}$ and $\mathring{\boldsymbol{\sigma}}\Delta t$ do not have physical meaning) – namely, $\boldsymbol{\sigma}$ is defined in the “spatial” coordinate system. Recall once more that this is necessary in order to have a consistent frame of reference for all of the elements in the finite element simulation. Thus, before we can consider the problem “solved” for any particular instance in time, we need to add an additional “step” to our solution procedure that obtains $\dot{\boldsymbol{\sigma}}$ from $\mathring{\boldsymbol{\sigma}}$ . We can do this by simply solving for $\dot{\boldsymbol{\sigma}}$ from our previous expression for the Jaumann rate of Cauchy Stress (eq. 3 in Section 7: Hypoelasticity). Thus, during each time step, the Jaumann rate is the stress that is used for the constitutive relationship, but the Cauchy stress is what is needed for the equation of motion (e.x. stress equilibrium).

It turns out that this same multi-step time-stepping “algorithm” will work for linear infinitesimal elasticity as well. We already mentioned that $\hat{\boldsymbol{\sigma}}$ is analogous to the linear infinitesimal stress tensor, and $\mathbf{E}$ is analogous to the linear infinitesimal strain tensor, $\boldsymbol{\epsilon}$ . In other words, the work-conjugate pair – defined in material coordinates – that is used in linear infinitesimal elasticity, can be thought of as a simplified version of $\hat{\boldsymbol{\sigma}} \leftrightarrow \mathbf{E}$ . $\dot{\hat{\boldsymbol{\sigma}}}$ and $\dot{\mathbf{E}}$ we know to be similarly work-conjugate and invariant to rigid body rotation.

Thus, one possible expression in linear infinitesimal elasticity that is analogous to our hypoelastic expression (eq. 7 in Section 7: Hypoelasticity) could be:

$\dot{\hat{\sigma}}_{ij}=\mathcal{C}_{ijkl} \dot{\epsilon}_{kl}$

where $\mathcal{C}_{ijkl}$ is a fourth-order tensor that relates linear infinitesimal stress, $\hat{\boldsymbol{\sigma}}$ , to linear infinitesimal strain, $\boldsymbol{\epsilon}$ , as we will see in the next chapter.

Obtaining the stress in spatial coordinates can be accomplished via eq. 2 in Section 4: Alternative Measures of Stress, namely, $\dot{\hat{\boldsymbol{\sigma}}}=\frac{d}{dt}\left[(det\mathbf{F})\mathbf{F}^{-1}\cdot\boldsymbol{\sigma}\cdot\mathbf{F}^{-T}\right]$ . After performing this time derivative, we could then rearrange to find the rate of Cauchy stress, $\dot{\boldsymbol{\sigma}}$ , in terms of $\dot{\hat{\boldsymbol{\sigma}}}$ . However, it may be desirable to take a different approach in order to avoid having to explicitly calculate the deformation gradient, $\mathbf{F}$ , or its time-derivative $\dot{\mathbf{F}}$ .

Additionally, let’s keep in mind that the linear infinitesimal stress and strain are not exactly equal to $\boldsymbol{\hat{\sigma}}$ and $\mathbf{E}$ , respectively. Thus, obtaining the Cauchy stress from the expression for $\boldsymbol{\hat{\sigma}}$ (eq. 2 in Section 4: Alternative Measures of Stress) may not be appropriate here because the linear infinitesimal stress is not exactly equal to $\boldsymbol{\hat{\sigma}}$ .

note: We cannot always simply replace the linear infinitesimal stress and strain with $\boldsymbol{\hat{\sigma}}$ and $\mathbf{E}$ . Although this is what we did, above, it was not derived because it is not mathematically correct. The above expression is shown primarily for illustrative purposes. It turns out that it is possible to relate $\boldsymbol{\hat{\sigma}}$ and $\mathbf{E}$ using $\boldsymbol{\mathcal{C}}$ , exactly, but this requires that we use the so-called Saint Venant-Kirchhoff strain energy density function.

Let’s instead consider our previous definitions of the Truesdell and the Jaumann rate, which, for infinitesimal elasticity, state that:

$\overset{\bigtriangledown}{\boldsymbol{\sigma}}\approx\mathring{\boldsymbol{\sigma}}\approx\mathbf{F}\cdot\dot{\hat{\boldsymbol{\sigma}}}\cdot\mathbf{F}^{T}$

Here, we should note that the term $\frac{1}{det\mathbf{F}}$ was eliminated since $det\mathbf{F}\approx1$ for infinitesimal deformations.

Since $\dot{\hat{\boldsymbol{\sigma}}}=\boldsymbol{\mathcal{C}}: \dot{\boldsymbol{\epsilon}}$ and $\mathring{\boldsymbol{\sigma}}\approx\mathbf{F}\cdot\dot{\hat{\boldsymbol{\sigma}}}\cdot\mathbf{F}^{T}$ , we can see that:

(1) $\begin{equation*} \mathring{\boldsymbol{\sigma}}=\boldsymbol{\mathcal{C}}:(\mathbf{F}\cdot\dot{\boldsymbol{\epsilon}}\cdot\mathbf{F}^{T}) \end{equation*}$

From Chapter 3, we know that eq. 1 conveniently reduces (with reasonable approximation) to:

(2) $\begin{equation*} \mathring{\boldsymbol{\sigma}}=\boldsymbol{\mathcal{C}}:\mathbf{D} \end{equation*}$

Or,

(3) $\begin{equation*} \mathring{\sigma}_{ij}=\mathcal{C}_{ijkl}D_{kl} \end{equation*}$

In eq. 2, we recall that $\mathbf{D}$ is the “Rate of Deformation Tensor.”

In linear infinitesimal elasticity, $\mathbf{D}$ is often simply referred to as the “strain rate.”

Some authors write the constitutive relationship for linear infinitesimal elasticity as $\boldsymbol{\sigma}=\boldsymbol{\mathcal{C}}:\mathbf{e}$ . The physical nature of $\mathbf{E}$ and $\mathbf{e}$ was illustrated in the example at the end of Section 4: Alternative Measures of Stress, where we saw that the two strain measures gave different results, even for the infinitesimal deformation case. Ultimately, however, one can argue that we are free to choose either stress-strain pair in linear infinitesimal elasticity. Applying the Jaumann operator to either constitutive expression results in the equation 2. Thus, our “time-stepping algorithm,” which we will see shortly, is independent of our “interpretation” of linear infinitesimal stress and strain.

The Kirchhoff stress, $\boldsymbol{\tau}$ , is defined as $\boldsymbol{\tau}=(det\mathbf{F})\boldsymbol{\sigma}$ . Sometimes $\boldsymbol{\tau}$ is used in place of $\boldsymbol{\sigma}$ , in eq. 2 for example, so that the small deformation approximation of $det\mathbf{F}\approx1$ does not need to be used.

Obtaining the stress rate in spatial coordinates simply requires that we solve for $\dot{\boldsymbol{\sigma}}$ from eq. 3 in Section 7: Hypoelasticity, namely,

(4) $\begin{equation*} \dot{\boldsymbol{\sigma}}=\mathring{\boldsymbol{\sigma}}+\mathbf{W} \cdot \boldsymbol{\sigma}+\boldsymbol{\sigma} \cdot \mathbf{W}^{T} \end{equation*}$

An alternative proof of eq. 4, for infinitesimal elasticity, is given in Appendix D.1

Thus, our constitutive expressions, in rate-form, for both linear infinitesimal elasticity and hyperelasticity (hypoelasticity), involve the Jaumann rate. For hypoelasticity, we had eq. 7 in Section 7: Hypoelasticity, while for linear infinitesimal elasticity, we have eq. 3. The following time-stepping algorithm, which allows us to obtain Cauchy (“spatial”) stress is identical for both hyperelasticity and linear infinitesimal elasticity.

note that the bases do not always remain orthogonal when there is shear (recall the second figure in Section 6: Phenomenological and Micromechanical Models, for example). The Jaumann rate in eq. 4 assumes orthogonal bases, which is a simplification.

General Framework For Problem Solving

1. Use an appropriate work-conjugate constitutive relationship:
$\mathring{\sigma}_{ij}=\Lambda_{ijkl}D_{kl}$ or $\mathring{\sigma}_{ij}=\mathcal{C}_{ijkl} D_{kl}$

2. Obtain a spatial representation of stress by first finding $\dot{\boldsymbol{\sigma}}$ from eq. 4, and then finding $\boldsymbol{\sigma}$ from $\boldsymbol{\sigma}_{n+1}=\boldsymbol{\sigma}_{n}+\dot{\boldsymbol{\sigma}}\Delta t$

3. Solve the equation of motion for the next increment in time:
$\sigma_{ij,j} + \rho b_{i}=\rho \frac{dv_{i}}{dt}$

The majority of this text is devoted to “1” and “2” (in a theoretical sense). Good course sequences on FEA spend a great deal of time on the implementation of these steps – in particular on step “3.” Recall that the expression given in step “3” was derived in Chapter 4.

The above “time-stepping algorithm” is just one possibility. The purpose here is merely to show one particular methodology that FEA software use for solving real-life problems. This methodology has only been shown to apply to isotropic elastic materials. Different time-stepping algorithms are surely used for other situations.

Plasticity

Obtaining Cauchy stress from the Jaumann rate is what is most often done by FEA software, even though the Truesdell rate form is better (exact) for elasticity. The Jaumann rate is known to have some problems with shear behavior due to the assumption of orthogonal bases [Valanis]. In addition, for highly compressible materials, the Truesdell rate can give significantly better results [Bazant]. Perhaps, then, the main reason that the Jaumann rate is preferred over the Truesdell rate (and other rates) has to do with plasticity.

In plasticity, many assumptions are made. The Truesdell rate turns out not to be good for plasticity, in general, because the dilatational (volume-changing) terms that relate $\dot{\boldsymbol{\sigma}}$ and the Truesdell rate are ignored in plasticity. So, suffice to say, plasticity is a complex topic that will not be covered in this introductory text, but the Jaumann rate, $\mathring{\boldsymbol{\sigma}}$ , is often used in FEA for solid elements, because it is computationally efficient and handles plasticity better than the Truesdell rate.

A general overview of plasticity equations used in FEA can be found in [Belytschko] and in [Simo]. Where finite strain elasticity and plasticity (e.x. rubber elastoplasticity) are considered, good theoretical references include [Lee], [Lubarda], and [Volokh].

V. A. Lubarda, Elastoplasticity Theory, Boca Raton, FL: CRC Press, 2002.
[Bibtex]

@BOOK{Lubarda,
    Address        = {Boca Raton, FL},
    Author         = {Vlado A. Lubarda},
    Edition        = {},
    Publisher      = {CRC Press},
    Title          = {Elastoplasticity Theory},
    Year           = {2002}
}

J. C. Simo and T. J. R. Hughes, Computational inelasticity, Springer, 1998.
[Bibtex]

@book{hughes,
  title={Computational inelasticity},
  author={Simo, J.C. and Hughes, T.J.R.},
  Publisher = {Springer},
  year={1998}
}

K. C. Valanis, “Back stress and Jaumann rates in finite plasticity,” International Journal of Plasticity, vol. 6, iss. 3, pp. 353-367, 1990.
[Bibtex]

@article{Valanis,
title = "Back stress and Jaumann rates in finite plasticity",
journal = "International Journal of Plasticity",
volume = "6",
number = "3",
pages = "353 - 367",
year = "1990",
note = "",
issn = "0749-6419",
doi = "DOI: 10.1016/0749-6419(90)90007-2",
author = "K.C. Valanis"
}

Z. P. Bazant, M. Gattu, and J. Vorel, “Work conjugacy error in commercial finite-element codes: its magnitude and how to compensate for it,” Proceedings of the Royal Society A: Mathematical, Physical and Engineering Science, 2012.
[Bibtex]

@article{Bazant,
author = {Bazant, Zdenek P. and Gattu, Mahendra and Vorel, Jan}, 
title = {Work conjugacy error in commercial finite-element codes: its magnitude and how to compensate for it},
year = {2012}, 
doi = {10.1098/rspa.2012.0167}, 
abstract ={}, 
URL = {}, 
eprint = {}, 
journal = {Proceedings of the Royal Society A: Mathematical, Physical and Engineering Science} 
}

E. H. Lee, Elastic-plastic deformation at finite strains, Stanford University, Divsion of Engineering Mechanics, 1968.
[Bibtex]

@book{Lee,
  title={Elastic-plastic deformation at finite strains},
  author={Lee, Erastus H},
  year={1968},
journal = "Transactions of the {A}{S}{M}{E}",
  publisher={Stanford University, Divsion of Engineering Mechanics}
}

T. Belytschko, B. Moran, and W. K. Liu, Nonlinear finite element analysis for continua and structures, Wiley, 1999, vol. 1.
[Bibtex]

@book{Belytschko,
  title={Nonlinear finite element analysis for continua and structures},
  author={Belytschko, Ted and Moran, Brian and Liu, Wing Kam},
  volume={1},
  year={1999},
  publisher={Wiley}
}

K. Volokh, “An approach to elastoplasticity at large deformations,” European Journal of Mechanics-A/Solids, 2012.
[Bibtex]

@article{Volokh,
  title={An approach to elastoplasticity at large deformations},
  author={Volokh, KY},
  journal={European Journal of Mechanics-A/Solids},
  year={2012},
  publisher={Elsevier}
}

Leave a Reply Cancel reply