HEX8 — 8-node trilinear hexahedron
==================================
The workhorse 3D solid kernel in femorph-solver. Eight corner
nodes map the reference cube
:math:`\hat\Omega = [-1, +1]^{3}` to a physical-space hex via
the isoparametric scheme of :doc:`../theory/isoparametric`.
24 DOFs per element (three translations per node).
* **Spec:** ``ELEMENTS.HEX8`` (default integration), or
``ELEMENTS.HEX8(integration="enhanced_strain")`` /
``ELEMENTS.HEX8(integration="plain_gauss")``
.. raw:: html
Verified vs MAPDL VM-1, Lamé closed form, NAFEMS LE5 — all benchmarks within 5 %
Cross-vendor mapping
--------------------
========================== ============================== ============================================
Solver Element name Notes
========================== ============================== ============================================
femorph-solver ``ELEMENTS.HEX8`` this page; ``integration="enhanced_strain"`` for thin-bending solids
ANSYS Mechanical APDL ``SOLID185`` KEYOPT(2)=2 → enhanced strain; KEYOPT(2)=3 → simplified-enhanced strain
NX / MSC Nastran ``CHEXA(8)`` PSOLID property, 8-node hex
Abaqus ``C3D8`` / ``C3D8I`` ``C3D8`` plain-Gauss; ``C3D8I`` adds incompatible modes (Simo-Rifai EAS analogue)
LS-DYNA ``ELEMENT_SOLID`` ELFORM=1 (under-integrated) / ELFORM=2 (full)
========================== ============================== ============================================
Restrictions
------------
Use a different element when:
* **Bending-dominated thin geometry** with :math:`L/h > 50` —
HEX8 even with EAS still admits some shear-locking residue.
Prefer ``QUAD4_SHELL`` or the planned ``SHELL181`` for
plate-bending dominated structures.
* **Highly-distorted hexes** (Jacobian determinant ratio > 10
across the element). HEX20 (20-node serendipity) tolerates
distortion better; for unstructured meshes prefer ``TET10``.
* **Singular fields** (cracks, sharp re-entrant corners).
Singular elements / quarter-point hexes are not yet
supported; use a refinement gradient instead.
Shape functions
---------------
Trilinear Lagrange on the reference cube, with corner-node
sign vector :math:`(\xi_i, \eta_i, \zeta_i) \in \{-1, +1\}^{3}`:
.. math::
N_i(\xi, \eta, \zeta) =
\tfrac{1}{8}
(1 + \xi_i\, \xi)(1 + \eta_i\, \eta)(1 + \zeta_i\, \zeta).
These satisfy the Kronecker-delta property
:math:`N_i(\boldsymbol{\xi}_j) = \delta_{ij}` at the eight
corner nodes and the partition-of-unity
:math:`\sum_i N_i \equiv 1` on :math:`\hat\Omega`. See
:ref:`sphx_glr_gallery_elements_solid185_example_hex8_shape_function_surfaces.py`
for slice-by-slice plots of the eight basis functions.
Integration
-----------
Stiffness uses the **2 × 2 × 2 Gauss-Legendre** rule (8 points,
degree 3 exact) — see :doc:`../theory/quadrature`. Three
integration-mode variants are exposed:
* ``ELEMENTS.HEX8`` (default) — full 2 × 2 × 2 Gauss with
**B-bar volumetric stabilisation** (Hughes 1980). The
volumetric component of the strain field is averaged across
the element to suppress *volumetric locking* when
:math:`\nu \to 0.5`. Recommended for nearly-incompressible
elasticity.
* ``ELEMENTS.HEX8(integration="enhanced_strain")`` — Simo–Rifai
9-parameter **enhanced assumed strain** (Simo & Rifai 1990).
Adds nine static-condensed internal modes that fix bending-
mode shear-locking in thin-bending problems. The natural
choice for plate / shell-bending solid kernels and for every
bending-dominated benchmark in
:doc:`/verification/index`.
* ``ELEMENTS.HEX8(integration="plain_gauss")`` — full 2 × 2 × 2
with **no** B-bar / EAS stabilisation. Useful only for
cross-solver comparisons against engines that ship the
vanilla 2 × 2 × 2 hex (e.g., scikit-fem); not recommended
for production runs because it locks both volumetrically
(large :math:`\nu`) and in shear (thin bending).
Mass
----
Consistent mass at 2 × 2 × 2 Gauss is the default. Row-sum
(HRZ) lumping (:func:`Model.solve_modal(lumped=True)
`) returns a diagonal mass
matrix; off-diagonal contributions are summed onto the diagonal
entry of each row, scaled to preserve total element mass. HRZ
is preferred over naive lumping for explicit dynamics — see
Cook Table 11.6-1.
Stress recovery
---------------
:func:`femorph_solver.result._stress_recovery.compute_nodal_stress`
evaluates :math:`\mathbf{B}(\boldsymbol{\xi}_\mathrm{node})\,
\mathbf{u}_e` at each of the eight element-local node positions,
applies the 6 × 6 isotropic elasticity matrix
:math:`\mathbf{C}` to get Voigt stress, then averages across
every element that touches a node — the standard nodal-averaging
stress-recovery scheme (Cook §6.14).
Performance characteristics
---------------------------
The HEX8 kernel is the most heavily-tuned in femorph-solver:
* Per-element ``ke``, ``me``, ``strain_batch`` route through a
C++ batch entry point so the assembly loop runs in native
code with vectorised inner products.
* Stiffness assembly on a 100 k-cell mesh is ~0.3 s on a
single core; mass assembly is comparable.
Verification cross-references
-----------------------------
* :ref:`sphx_glr_gallery_elements_solid185_example_solid185.py` — single-hex uniaxial.
* :ref:`sphx_glr_gallery_verification_example_verify_cantilever_eb.py` — Euler-Bernoulli cantilever (EAS).
* :ref:`sphx_glr_gallery_verification_example_verify_lame_cylinder.py` — Lamé thick cylinder.
* :ref:`sphx_glr_gallery_verification_example_verify_clamped_plate_static.py` — NAFEMS LE5 clamped plate.
* :ref:`sphx_glr_gallery_verification_example_verify_ss_plate_modes.py` — SS plate Kirchhoff fundamental.
Implementation: :mod:`femorph_solver.elements.hex8`.
References
----------
* Zienkiewicz, O. C. and Taylor, R. L. (2013) *The Finite
Element Method: Its Basis and Fundamentals*, 7th ed.,
Butterworth-Heinemann, §6 (8-node hex Lagrange family),
§6.6.2 (volumetric locking + B-bar).
* Cook, R. D., Malkus, D. S., Plesha, M. E., Witt, R. J.
(2002) *Concepts and Applications of Finite Element
Analysis*, 4th ed., Wiley, §6.2 (hex shape functions),
§6.14 (nodal stress recovery), Table 11.6-1 (HRZ lumping).
* Hughes, T. J. R. (1980) "Generalization of selective
integration procedures to anisotropic and nonlinear media,"
*International Journal for Numerical Methods in
Engineering* 15 (9), 1413–1418.
* Simo, J. C. and Rifai, M. S. (1990) "A class of mixed
assumed strain methods and the method of incompatible
modes," *International Journal for Numerical Methods in
Engineering* 29 (8), 1595–1638.
* Bathe, K.-J. (2014) *Finite Element Procedures*, 2nd ed.,
§5.3.1 (8-node hex), §5.3.4 (incompatible modes / EAS).