QUAD4_SHELL — 4-node Mindlin–Reissner shell =========================================== *MAPDL-compatible alias:* ``SHELL181``. **Kinematics.** Flat four-node first-order shear-deformation shell (Mindlin–Reissner theory). Six DOFs per node (:math:`u_x, u_y, u_z, \theta_x, \theta_y, \theta_z`); 24 DOFs per element in the local frame. **Stiffness.** :math:`K_e = K_m + K_b + K_s + K_\text{drill}` — membrane + bending + transverse shear + small drilling stabilisation. **Mass.** Consistent :math:`M_e` includes both translational and (through-thickness) rotational inertia. **MAPDL compatibility.** Matches SHELL181 at default settings (flat shell, selective-reduced transverse shear, drilling-DOF stabilisation). .. contents:: On this page :local: :depth: 2 ---- Mindlin–Reissner kinematics --------------------------- First-order shear-deformation theory [Hughes1987]_ [Bathe2014]_ describes the displacement of a point at through-thickness coordinate :math:`z` in terms of the mid-surface displacement :math:`(u, v, w)` and two rotation DOFs :math:`(\theta_x, \theta_y)`: .. math:: u_z(x, y, z) = u(x, y) + z\,\theta_y(x, y), \\ v_z(x, y, z) = v(x, y) - z\,\theta_x(x, y), \\ w_z(x, y, z) = w(x, y). The transverse normal remains straight but *not* necessarily perpendicular to the deformed mid-surface — the angle between the normal and the mid-surface is the transverse shear strain, captured by .. math:: \gamma_{xz} = \partial_x w + \theta_y, \qquad \gamma_{yz} = \partial_y w - \theta_x. The in-plane strains split into membrane (:math:`z^0`) and bending (:math:`z^1`) contributions after through-thickness integration: .. math:: \varepsilon_m = \begin{bmatrix} \partial_x u \\ \partial_y v \\ \partial_x v + \partial_y u \end{bmatrix}, \qquad \varepsilon_b = \begin{bmatrix} \partial_x \theta_y \\ -\partial_y \theta_x \\ \partial_x \theta_x - \partial_y \theta_y \end{bmatrix}. Local frame and bilinear interpolation -------------------------------------- The element is assembled in a local 2-D frame with :math:`\hat{x}` along edge 1→2, :math:`\hat{z}` the best-fit normal, and :math:`\hat{y} = \hat{z} \times \hat{x}`. Nodes project to this plane and the element uses standard bilinear shape functions on :math:`(\xi, \eta) \in [-1, 1]^2`: .. math:: N_i(\xi, \eta) = \tfrac{1}{4}\,(1 + \xi_i \xi)\,(1 + \eta_i \eta), with corner nodes at :math:`\{\pm 1\}^2` in CCW order — same as PLANE182 / VTK_QUAD. Selective reduced integration (membrane, bending, shear) -------------------------------------------------------- Locking management splits the stiffness into three integrals with different Gauss rules [Hughes1987]_ [ZTZ2013]_: * **Membrane**: 2 × 2 Gauss (full). * **Bending**: 2 × 2 Gauss (full). * **Transverse shear**: 1 × 1 Gauss (reduced). The reduced transverse-shear integration is the classical cure for **shear locking** on thin shells: full integration over-stiffens the element dramatically as the thickness-to-span ratio :math:`t/L \to 0`. MITC4 [BatheDvorkin1985]_ is the alternative cure (via tying of the shear strain at specific points); femorph-solver uses selective reduced integration, which matches MAPDL's default SHELL181 behaviour. Through-thickness integration uses the Reissner–Mindlin shear correction factor :math:`\kappa = 5/6` on the transverse-shear stiffness. Plane-stress constitutive ------------------------- SHELL181 runs in the plane-stress regime (:math:`\sigma_{zz} = 0` in the local frame). Membrane and bending share the 3 × 3 plane-stress elastic matrix: .. math:: C_\text{ps} = \frac{E}{1 - \nu^2} \begin{bmatrix} 1 & \nu & 0 \\ \nu & 1 & 0 \\ 0 & 0 & (1 - \nu) / 2 \end{bmatrix}. Through-thickness constants: .. math:: D_m = C_\text{ps}\, t, \qquad D_b = \frac{t^3}{12}\, C_\text{ps}, \qquad D_s = \kappa\, G\, t\, I_2, with :math:`G = E / (2 (1 + \nu))`. Drilling-DOF stabilisation -------------------------- A flat plate has no physical stiffness for the local :math:`\theta_z` rotation — the drilling DOF. Without stabilisation, the 24 × 24 local :math:`K` has four extra zero-energy modes (one per corner) and global assembly can become singular. femorph-solver adds a small penalty term .. math:: K_\text{drill} = \alpha\, G\, t\, A_e \, \text{diag}_i\!\bigl(\delta_{i, \theta_z^\text{local}}\bigr), with :math:`\alpha \approx 10^{-3}` — enough to stabilise, small enough not to bleed into physical bending / membrane stiffness. A block-diagonal 3 × 3 direction-cosine matrix then rotates every node's 6-DOF block from the local frame back to global. Numpy walk-through (skeleton) ----------------------------- .. code-block:: python import numpy as np # Bilinear shape functions on the reference quad. xi_n = np.array([(-1,-1), (+1,-1), (+1,+1), (-1,+1)], dtype=float) def shape_and_derivs(xi, eta): N = 0.25 * (1 + xi_n[:, 0] * xi) * (1 + xi_n[:, 1] * eta) dN = np.empty((4, 2)) dN[:, 0] = 0.25 * xi_n[:, 0] * (1 + xi_n[:, 1] * eta) dN[:, 1] = 0.25 * xi_n[:, 1] * (1 + xi_n[:, 0] * xi) return N, dN # Full 2x2 for membrane + bending, single centre point for shear. g = 1.0 / np.sqrt(3.0) gp_full = np.array([(-g,-g), (+g,-g), (+g,+g), (-g,+g)]) w_full = np.ones(4) gp_shear = np.array([(0.0, 0.0)]) w_shear = np.array([4.0]) # full square area in natural coords # Assemble K_m and K_b at gp_full with (D_m, D_b); # assemble K_s at gp_shear with D_s; sum the three; # add K_drill; rotate the 24x24 matrix to global. Real constants and KEYOPT ------------------------- .. list-table:: :widths: 25 75 :header-rows: 1 * - Slot - Meaning * - ``real[0]`` - Through-thickness dimension :math:`t` (mandatory). The broader MAPDL SHELL181 KEYOPT table (composite lay-ups, shell offsets, incompatible-modes extension) is not yet exposed; the present implementation targets KEYOPT(3) = 0 (flat, uniform thickness). Validation ---------- **Rigid-body modes.** 6 rigid-body zeros + 18 elastic eigenvalues on an undistorted square shell. **Bending patch.** A slender cantilever strip (length :math:`L`, width :math:`b`, thickness :math:`t`) with tip moment :math:`M` should deflect to :math:`u_\text{tip} = M L^2 / (2 E I)` with :math:`I = b t^3 / 12`. One-element-wide mesh reaches this limit to within a few percent via selective-reduced shear; plain full integration would be off by orders of magnitude on the thin-shell limit. **MacNeal–Harder benchmarks.** The twisted-beam and pinched-cylinder problems from [MacNealHarder1985]_ are the community-standard shell stress tests; they're on the roadmap for the verification gallery. API reference ------------- .. autoclass:: femorph_solver.elements.shell181.SHELL181 :members: References ---------- .. [Hughes1987] Hughes, T. J. R. (1987 / Dover 2000). *The Finite Element Method: Linear Static and Dynamic Finite Element Analysis*. Ch. 5 (Mindlin plates / shells), §5.2 (selective integration). https://store.doverpublications.com/0486411818.html .. [Bathe2014] Bathe, K.-J. (2014). *Finite Element Procedures*, 2nd ed. Ch. 5 & 6 on plates and shells. https://www.klausjurgenbathe.com/fepbook/ .. [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 .. [BatheDvorkin1985] Bathe, K.-J. and Dvorkin, E. N. (1985). "A four-node plate bending element based on Mindlin/Reissner plate theory and a mixed interpolation." *International Journal for Numerical Methods in Engineering* **21** (2), 367–383. https://doi.org/10.1002/nme.1620210213 .. [MacNealHarder1985] MacNeal, R. H. and Harder, R. L. (1985). "A proposed standard set of problems to test finite element accuracy." *Finite Elements in Analysis and Design* **1** (1), 3–20. https://doi.org/10.1016/0168-874X(85)90003-4