HEX20 — 20-node quadratic hexahedral solid ========================================== *MAPDL-compatible alias:* ``SOLID186``. Degenerate forms ``WEDGE15`` and ``PYR13`` (MAPDL ``SOLID186W`` / ``SOLID186P``) auto-dispatch from this entry. **Kinematics.** 20-node serendipity hexahedron; three translational DOFs per node (:math:`u_x, u_y, u_z`); 60 DOFs per element. Also dispatches degenerate wedge (15-node) and pyramid (13-node) shapes automatically when a CDB-stored element has collapsed corners. **Stiffness.** :math:`K_e = \int_{\Omega_e} B^\top C B\,\mathrm{d}V` integrated on the reference cube with 3×3×3 Gauss-Legendre (full) or 2×2×2 (reduced). **Mass.** Consistent :math:`M_e = \int \rho\, N^\top N\,\mathrm{d}V` with 3×3×3 Gauss. **MAPDL compatibility.** Reproduces SOLID186 at KEYOPT(2) = 1 (full 3×3×3 integration) and KEYOPT(2) = 0 (uniform-reduced 2×2×2 without explicit hourglass control). .. contents:: On this page :local: :depth: 2 ---- Shape functions (serendipity 20-node hex) ----------------------------------------- The 20-node brick is the quadratic serendipity extension of the trilinear SOLID185 [ZTZ2013]_ [Bathe2014]_. Eight corner nodes sit at :math:`(\xi_i, \eta_i, \zeta_i) \in \{\pm 1\}^3`; twelve mid-edge nodes sit where exactly one natural coordinate is zero. There are **no** mid-face or body-centre nodes — that's the hallmark of the serendipity family [Cook2001]_. * **Corner node** :math:`i` (all three :math:`|\xi_i^c| = 1`): .. math:: N_i(\xi, \eta, \zeta) = \tfrac{1}{8}\,(1 + \xi_i \xi)\,(1 + \eta_i \eta)\, (1 + \zeta_i \zeta)\,\bigl(\xi_i \xi + \eta_i \eta + \zeta_i \zeta - 2\bigr). * **Mid-edge node** on a :math:`\xi`-directed edge (:math:`\xi_i^m = 0`): .. math:: N_i(\xi, \eta, \zeta) = \tfrac{1}{4}\,(1 - \xi^2)\,(1 + \eta_i \eta)\, (1 + \zeta_i \zeta), with the symmetric forms for :math:`\eta`- and :math:`\zeta`-directed mid-edges. Jacobian, :math:`B` matrix, and the isotropic constitutive :math:`C` are the same six-component Voigt objects as for SOLID185 — only the shape-function and derivative entries change. 3×3×3 Gauss quadrature ---------------------- Full integration uses 27 points at natural coordinates :math:`\xi_g \in \{-\sqrt{3/5},\, 0,\, +\sqrt{3/5}\}` (and similarly for :math:`\eta`, :math:`\zeta`) with weights :math:`w_g` from the tensor product of the one-dimensional 3-point rule (:math:`5/9,\ 8/9,\ 5/9`). The rule is exact for polynomials up to degree 5 in each direction, which integrates :math:`B^\top C B` on an undistorted hex without approximation. Reduced (``"reduced"``) integration drops to the 2×2×2 rule (same points as SOLID185). On a single element this produces six extra zero-energy modes (hourglass modes); in a patch of elements the modes usually don't connect and the global :math:`K` stays stable. MAPDL's KEYOPT(2) = 0 option applies explicit hourglass-control stabilisation that femorph-solver does *not* yet implement, so use ``"full"`` in any production run. Numpy walk-through ------------------ The following listing reuses the SOLID185 pattern with the 20-node shape functions. On a unit cube with generic steel it produces the expected 6 rigid-body + 54 elastic eigenvalues: .. code-block:: python import numpy as np # Natural coordinates of the 20 nodes, in MAPDL / VTK order: # corners (bottom 4, top 4), bottom-face midsides, top-face # midsides, vertical midsides. xi_n = np.array([ (-1,-1,-1), (+1,-1,-1), (+1,+1,-1), (-1,+1,-1), (-1,-1,+1), (+1,-1,+1), (+1,+1,+1), (-1,+1,+1), ( 0,-1,-1), (+1, 0,-1), ( 0,+1,-1), (-1, 0,-1), ( 0,-1,+1), (+1, 0,+1), ( 0,+1,+1), (-1, 0,+1), (-1,-1, 0), (+1,-1, 0), (+1,+1, 0), (-1,+1, 0), ], dtype=float) # 3x3x3 Gauss-Legendre. p1d = np.array([-np.sqrt(3/5), 0.0, +np.sqrt(3/5)]) w1d = np.array([5/9, 8/9, 5/9]) gp = np.array([(a, b, c) for a in p1d for b in p1d for c in p1d]) w = np.array([a * b * c for a in w1d for b in w1d for c in w1d]) def shape_and_derivs(xi, eta, zeta): N = np.empty(20) dN = np.empty((20, 3)) for i, (xi_i, eta_i, zeta_i) in enumerate(xi_n): one_xi = 1.0 + xi_i * xi one_eta = 1.0 + eta_i * eta one_zeta = 1.0 + zeta_i * zeta if abs(xi_i) + abs(eta_i) + abs(zeta_i) == 3: # corner combo = xi_i * xi + eta_i * eta + zeta_i * zeta - 2.0 N[i] = 0.125 * one_xi * one_eta * one_zeta * combo dN[i, 0] = 0.125 * one_eta * one_zeta * (xi_i * combo + one_xi * xi_i) dN[i, 1] = 0.125 * one_xi * one_zeta * (eta_i * combo + one_eta * eta_i) dN[i, 2] = 0.125 * one_xi * one_eta * (zeta_i * combo + one_zeta * zeta_i) elif xi_i == 0: # mid-edge on ξ q = 1.0 - xi * xi N[i] = 0.25 * q * one_eta * one_zeta dN[i, 0] = 0.25 * (-2.0 * xi) * one_eta * one_zeta dN[i, 1] = 0.25 * q * eta_i * one_zeta dN[i, 2] = 0.25 * q * one_eta * zeta_i elif eta_i == 0: # mid-edge on η q = 1.0 - eta * eta N[i] = 0.25 * q * one_xi * one_zeta dN[i, 0] = 0.25 * q * xi_i * one_zeta dN[i, 1] = 0.25 * (-2.0 * eta) * one_xi * one_zeta dN[i, 2] = 0.25 * q * one_xi * zeta_i else: # mid-edge on ζ q = 1.0 - zeta * zeta N[i] = 0.25 * q * one_xi * one_eta dN[i, 0] = 0.25 * q * xi_i * one_eta dN[i, 1] = 0.25 * q * one_xi * eta_i dN[i, 2] = 0.25 * (-2.0 * zeta) * one_xi * one_eta return N, dN # Assemble K_e the same way as the SOLID185 listing: loop over Gauss # points, compute J = dN^T X, detJ, dN/dx = J^{-1} dN/dξ, build the # 6x60 B, accumulate K += B^T C B detJ w. The element passes the patch test to machine precision: a linear displacement field ramps exactly through every node of any distorted hex. Degenerate shapes — wedge and pyramid ------------------------------------- MAPDL stores pyramid and wedge elements as 20-slot hex connectivity with repeated node IDs at collapsed corners: * **Wedge (15 nodes)** — corner collapse ``K = L`` and ``O = P`` on adjacent vertical edges. * **Pyramid (13 nodes)** — corner collapse ``M = N = O = P`` on all four top-face nodes (the apex). * **Degenerate tet** would collapse further, but in practice MAPDL writes those as SOLID187. :mod:`mapdl_archive` parses the degenerate connectivity and emits ``VTK_QUADRATIC_WEDGE`` (15 pts) or ``VTK_QUADRATIC_PYRAMID`` (13 pts) cells. femorph-solver's :class:`SOLID186W ` and :class:`SOLID186P ` wrappers restore the full 20-slot layout, delegate :math:`K_e` / :math:`M_e` to the plain SOLID186 kernel, and fold the 60 × 60 result back down to 45 × 45 (wedge) or 39 × 39 (pyramid) with a constant indicator matrix :math:`T`: .. math:: K_\text{deg} = T^\top K_\text{hex}\, T. This is exactly what the global assembler would compute if the element record had 20 node IDs with repeats — bit-for-bit identical to MAPDL's degenerate-hex formulation. Validated on a 5 166-element rotor-blade CDB (3 069 regular hex + 497 pyramid + 19 wedge + 1 581 tet): 12-mode modal parity to 0.08 % relative vs MAPDL (the residual is the plain-Gauss vs B-bar offset, not the degenerate dispatch). Formulation flag ---------------- .. list-table:: :widths: 25 25 50 :header-rows: 1 * - Name - MAPDL KEYOPT(2) - Description * - ``"full"`` - 1 - 3×3×3 Gauss (27 points). Fully ranked single-element :math:`K_e`. * - ``"reduced"`` - 0 - 2×2×2 Gauss (8 points). 6 hourglass modes per element. No stabilisation yet — caller's responsibility at mesh level. .. code-block:: python m.materials[1]["_SOLID186_FORMULATION"] = "reduced" Validation ---------- **Rigid-body modes.** 54 non-zero elastic modes + 6 rigid-body zeros on an undistorted hex (``"full"``); 48 + 12 (6 rigid + 6 hourglass) for ``"reduced"``. **Degenerate-shape correctness.** SOLID186W on a 15-node wedge and SOLID186P on a 13-node pyramid each have exactly 6 rigid-body zero eigenvalues on the reference element. **MAPDL parity.** Plain modal on a quarter-arc blade sector (fixed root, 12 lowest modes) agrees with MAPDL to 0.08 % relative. API reference ------------- .. autoclass:: femorph_solver.elements.solid186.SOLID186 :members: .. autoclass:: femorph_solver.elements.solid186_degenerate.SOLID186W :members: .. autoclass:: femorph_solver.elements.solid186_degenerate.SOLID186P :members: References ---------- .. [ZTZ2013] Zienkiewicz, O. C., Taylor, R. L., and Zhu, J. Z. (2013). *The Finite Element Method: Its Basis and Fundamentals*, 7th ed. Butterworth-Heinemann. https://doi.org/10.1016/C2009-0-24909-9 .. [Bathe2014] Bathe, K.-J. (2014). *Finite Element Procedures*, 2nd ed. Prentice-Hall / Klaus-Jürgen Bathe. Ch. 5 (isoparametric elements). https://www.klausjurgenbathe.com/fepbook/ .. [Cook2001] Cook, R. D., Malkus, D. S., Plesha, M. E., and Witt, R. J. (2001). *Concepts and Applications of Finite Element Analysis*, 4th ed. Wiley. ISBN 978-0-471-35605-9. https://www.wiley.com/en-us/Concepts+and+Applications+of+Finite+Element+Analysis%2C+4th+Edition-p-9780471356059