Source code for femorph_solver.validation._problems.cantilever_udl

"""Clamped cantilever under uniform distributed transverse load.

Reference geometry
------------------

Slender cantilever of length :math:`L`, rectangular cross-section
:math:`w \\times h`, isotropic elastic material with Young's
modulus :math:`E`.  Clamped at :math:`x = 0`, loaded by a uniform
transverse distributed load :math:`q` (force per unit length)
running along the upper surface, acting in the :math:`-z`
direction (gravity-like).

Closed-form quantities
----------------------

For a slender Euler-Bernoulli cantilever::

    EI w''''(x) = q,   w(0) = w'(0) = 0,
    EI w'''(L) = 0 (tip shear),  EI w''(L) = 0 (tip moment).

Integrating four times with the given BCs:

* **Tip deflection** (Timoshenko §5.4)

  .. math::

      \\delta_{\\text{tip}} = \\frac{q L^{4}}{8 E I}, \\qquad
      I = \\frac{w h^{3}}{12}.

* **Root reaction force**: by global equilibrium,
  :math:`R = q L` (carried by the clamped face).

References
----------
* Timoshenko, S. P.  *Strength of Materials, Part I*, 3rd ed.,
  Van Nostrand, 1955, §5.4 — standard cantilever deflection
  under uniformly distributed load.
* Gere & Goodno (2018), *Mechanics of Materials* 9th ed., §9.3
  Table 9-1, case 1 — tabulated tip deflection for a prismatic
  cantilever with UDL.

Cross-references (public verification manuals — fair-use
citation of problem IDs only; no vendor content vendored):

* **Abaqus AVM 1.5.x** cantilever-with-UDL family.

Implementation notes
--------------------

The UDL is applied as nodal forces distributed across the top
face of the hex mesh: the force per node is ``q * L / n_top``
where ``n_top`` is the number of top-face nodes.  This is the
standard consistent-load lumping for a uniformly-loaded
structured grid — the per-node force is constant along the
span, mirroring the physical load.  End-effects from the
discrete edge nodes are absorbed by the shape-function
integration and vanish on mesh refinement.
"""

from __future__ import annotations

from dataclasses import dataclass
from typing import Any

import numpy as np
import pyvista as pv

import femorph_solver
from femorph_solver import ELEMENTS
from femorph_solver.validation._benchmark import (
    BenchmarkProblem,
    PublishedValue,
)


[docs] @dataclass class CantileverUDL(BenchmarkProblem): """Cantilever tip deflection under uniform distributed load.""" name: str = "cantilever_udl" description: str = ( "Slender clamped cantilever under uniformly-distributed " "transverse load — Euler-Bernoulli closed-form " "(Timoshenko 1955 §5.4)." ) analysis: str = "static" #: Cantilever length [m]. L: float = 1.0 #: Cross-section width [m]. width: float = 0.05 #: Cross-section height [m]. height: float = 0.05 #: Young's modulus [Pa]. E: float = 2.0e11 #: Poisson's ratio. nu: float = 0.3 #: Density [kg/m^3] — irrelevant for the static solution but #: required by the material stamp. rho: float = 7850.0 #: Uniform distributed load along the span [N/m] (positive #: value, acts in ``-z``). q: float = 1.0e3 @property def published_values(self) -> tuple[PublishedValue, ...]: I = self.width * self.height**3 / 12.0 # noqa: E741 delta = self.q * self.L**4 / (8.0 * self.E * I) return ( PublishedValue( name="tip_deflection", value=delta, unit="m", source="Timoshenko 1955 §5.4", formula="delta = q L^4 / (8 E I)", # UDL bending is less sensitive to shear-lock than a # tip shear load (the bending moment tapers linearly # from root to tip), so the same EAS formulation # gets tighter than 5 % at the reference mesh. tolerance=0.05, ), )
[docs] def build_model(self, **mesh_params: Any) -> femorph_solver.Model: nx = int(mesh_params.get("nx", 40)) ny = int(mesh_params.get("ny", 3)) nz = int(mesh_params.get("nz", 3)) xs = np.linspace(0.0, self.L, nx + 1) ys = np.linspace(0.0, self.width, ny + 1) zs = np.linspace(0.0, self.height, nz + 1) grid = pv.StructuredGrid( *np.meshgrid(xs, ys, zs, indexing="ij") ).cast_to_unstructured_grid() m = femorph_solver.Model.from_grid(grid) m.assign( ELEMENTS.HEX8(integration="enhanced_strain"), material={"EX": self.E, "PRXY": self.nu, "DENS": self.rho}, ) pts = np.asarray(m.grid.points) clamped = np.where(pts[:, 0] < 1e-9)[0] m.fix(nodes=(clamped + 1).tolist(), dof="ALL") # Distribute the UDL across the top face (z = height). We # lump per-column so each structured-grid column along the # span carries its share of the line-load ``q * dx``; # within the column the force is split equally across # transverse nodes (ny+1 of them). This matches the # consistent-load nodal equivalent for linear hexes under # a uniform surface traction. top = np.where(pts[:, 2] > self.height - 1e-9)[0] # Number of columns along the span = nx+1; per-column line- # load is ``q * dx`` where ``dx = L / nx``. Trapezoid rule: # end columns carry half-weight. dx = self.L / nx for n in top: x = pts[n, 0] # End columns (x = 0 or x = L) get half-weight; # interior columns full-weight. col_weight = 0.5 if (x < 1e-9 or x > self.L - 1e-9) else 1.0 # Split the per-column line-load across the transverse # (y) nodes in the column. Trapezoid rule again: edge # y-nodes get half the weight of interior ones. Same # treatment avoids over-loading the corners. y = pts[n, 1] y_weight = 0.5 if (y < 1e-9 or y > self.width - 1e-9) else 1.0 # Normalising factors so the total force integrates # to ``q * L`` exactly. # sum(col_weights) = nx (trapezoid rule over nx+1 # nodes with edge half-weight). # sum(y_weights) = ny. fz = -(self.q * dx * col_weight * y_weight) / ny m.apply_force(int(n + 1), fz=fz) m._bench_clamped_nodes = clamped # type: ignore[attr-defined] # Tip mid-line nodes for extract (x == L, y ≈ width/2, top + bottom avg). tip_mask = pts[:, 0] > self.L - 1e-9 m._bench_tip_nodes = np.where(tip_mask)[0] # type: ignore[attr-defined] return m
[docs] def extract(self, model: femorph_solver.Model, result: Any, name: str) -> float: u = np.asarray(result.displacement).reshape(-1, 3) if name == "tip_deflection": tip_nodes = model._bench_tip_nodes # type: ignore[attr-defined] # Centroidal tip deflection in -z → report the magnitude. return float(-u[tip_nodes, 2].mean()) raise KeyError(f"unknown quantity {name!r}")