"""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