Continuum deformation

Note that the conservative elastic material forces \(f_{\textrm{e}}\) in general require some rather involved computation. Since this computation is not straightforward, we briefly explain the derivation of the nonlinear hyper-elastic material forces, which follows standard nonlinear finite element procedures ¹²³. \(\renewcommand{\p}{\partial}\) \(\renewcommand{\R}{\mathbb{R}}\) \(\renewcommand{\vec}[1]{\mathbb{#1}}\) \(\renewcommand{\FB}{\mathbf{F}}\) \(\renewcommand{\fB}{\mathbf{f}}\) \(\renewcommand{\x}{\mathbf{x}}\)

Intermezzo: deformation gradient

A fundamental measure of deformation in continuum mechanics is the deformation gradient, denoted by \(\mathbf{F}\). The deformation gradient characterizes the local deformation for a neighborhood of the continuum body \(\Omega\). Since a subvolume of the continuum body cannot be reduced to a point, it follows that \(\det{\mathbf{F}} = J > 0\) and \(\mathbf{F}^{-1}\) exists. The term \(J\) is called the relative volume change and it is equal to 1 for isochoric deformations, such as rigid body deformations. Given these properties, the deformation gradient can then be factorized into \(\mathbf{F} = \mathbf{Q} \mathbf{V}\), where \(\mathbf{V} \succ 0\) is the right-handed stretch tensor and \(\mathbf{Q} \in \mathrm{SO}(3)\) is a rotation matrix belonging to the special orthogonal group ¹². For convenience, we summarize the derived quantities of \(\mathbf{F}\) in the table below that will be used throughout this section.

Intermezzo: Derivation of hyperelastic forces

Let \(\Omega_i\) denote the subspace spanned by the \(i\)-th element of the finite element mesh, and let \(\x_i\) denote its nodal displacement vector. The elasticity of the constitutive soft material can be described by a strain-energy density function \(\Psi: \FB \to \R_{\ge 0}\). A comprehensive discussion on common constitutive models for \(\Psi\) will be provided later in the subsequent paragraph. The elastic potential energy of the continuum body is given by \(\mathcal{U}_e = \int_\Omega \Psi(\cdot) \; dV\), and the conservative hyper-elastic force contribution can be computed as \(\fB_e := \nabla_{\x}\,\mathcal{U}_e\). This contribution can be approximated using piecewise finite element interpolation and integrated using the Gauss quadrature rule¹ as follows:

where the Gauss weights are denoted by \(w_j > 0\), and the number of finite elements and Gauss samples are represented by \(N_e\) and \(N_w\), respectively.

The term \({\partial \Psi}/{\partial \FB}\) is also referred to as the first Piolla-Kirchhoff (PK1) stress tensor, which can be represented in closed-form for many constitutive models. The term \({\partial \FB}/{\partial \x_e}\) denotes the deformation Jacobian, which can also be given in closed-form but depends on the choice of element type. In addition to the first Piolla stress tensor, we also introduce the Cauchy stress tensor (i.e., true stress) \(\vec{\sigma} := J^{-1} \frac{\partial \Psi}{\partial \FB} \FB^\top\), which is a symmetric second-order tensor whose components represent the true stress. It should be noted that these tensor calculations are highly nonlinear, making their computation the most time-consuming aspect of the finite element assembly. To enhance computational efficiency, the toolkit employs \texttt{.mex} executable code, generated during installation (\texttt{Matlab Coder} is required).