Biomedical Engineering Reference
In-Depth Information
¼
R
e
K
e
R
e
x
x
0
f
ext
e
M
e
€
x
þ
C
e
_
x
þ
(6)
Although the method cannot model material nonlinearities, it offers a very
efficient way to achieve a geometric nonlinear formulation. This is due to the fact
that only one linear system has to be solved for each time step if the semidiscrete
equation of motion is discretized using an implicit time integration technique.
However, it is important to point out that, in contrast to the linear FEM, the rotation
matrices have to be computed and assembled into the stiffness matrix every time step.
In summary, the corotational FEM requires significantly less computational
effort than a fully nonlinear formulation. Also, it is very stable as the extraction
of the rotational component changes the condition number of the element matrices
only marginally.
2.2 Corotated Quadratic Tetrahedra
A FE formulation based on quadratic shape functions can offer several advantages
over linear elements. First of all the convergence properties are much better if the
solution is sufficiently smooth. Furthermore, the accuracy of linear elements can be
severely degraded by volume locking during the simulation of soft tissue. This is in
particular true if linear tetrahedra are used. Quadratic tetrahedra significantly
reduce this artificial stiffness. Finally, isoparametric quadratic tetrahedra are better
suited to capture curved geometries than linear elements.
One integration point is needed for the linear tetrahedron in order to evaluate the
stiffness integral. Consequently, the element's rotation has to be determined only
once per element. This approach can be naturally extended to quadratic elements. In
accordance with the formulation by Mezger et al. [
10
], we perform a polar decom-
position at each of the four integration points of the quadratic tetrahedron.
In order to speed up the simulation, the element rotation matrices
R
e
are not
explicitly built. Instead, at each integration point we calculate the rotation matrix
R
of the deformation gradient
r'
. Furthermore, the tangential stiffness matrix
K
iIjJ
is
re-arranged into ten 3
3 matrices
K
IJ
that each represent the derivative of the
I
-th
nodal force in the direction of the
J
-th nodal spatial coordinates. For each integra-
tion point, the elemental co-rotated nodal forces then read
f
int
I
RK
IJ
R
T
x
J
x
J
¼
(7)
and the co-rotated elemental stiffness matrix can be assembled by summing up
K
CR
IJ
RK
IJ
R
T
¼
(8)
over each integration point. The polar decomposition can be obtained in a fast and
robust way by applying the iterative scheme proposed by Higham et al. [
11
].
