.. DO NOT EDIT. .. THIS FILE WAS AUTOMATICALLY GENERATED BY SPHINX-GALLERY. .. TO MAKE CHANGES, EDIT THE SOURCE PYTHON FILE: .. "gallery/verification/example_verify_cantilever_udl.py" .. LINE NUMBERS ARE GIVEN BELOW. .. only:: html .. note:: :class: sphx-glr-download-link-note :ref:`Go to the end ` to download the full example code. .. rst-class:: sphx-glr-example-title .. _sphx_glr_gallery_verification_example_verify_cantilever_udl.py: Cantilever under uniformly distributed load — Euler–Bernoulli closed form ========================================================================= For a clamped-free prismatic cantilever of length :math:`L`, flexural rigidity :math:`EI`, and uniform transverse line load :math:`w` (force per unit length), the tip deflection is .. math:: \delta_\text{tip} = \frac{w L^{4}}{8 E I}, \qquad I = \frac{w_b h^{3}}{12} (Timoshenko 1955 §5.4, Cook 2002 §2.4–§2.5). The corresponding slope and curvature at the clamp are .. math:: \theta_\text{tip} = \frac{w L^{3}}{6 E I}, \qquad \kappa(0) = \frac{M(0)}{E I} = \frac{w L^{2}}{2 E I}, so the maximum bending stress at the root extreme fibre is :math:`\sigma_\text{max} = E\, \kappa(0)\, c = w L^{2} c / (2 I)`. References ---------- * Timoshenko, S. P. *Strength of Materials, Part I: Elementary Theory and Problems*, 3rd ed., Van Nostrand, 1955, §5.4. * Cook, R. D., Malkus, D. S., Plesha, M. E., Witt, R. J. (2002) *Concepts and Applications of Finite Element Analysis*, 4th ed., Wiley, §2.4 (Hermite cubics) — the cantilever-UDL example is the canonical test of consistent-load integration on the beam element. Implementation note ------------------- The cantilever is meshed as a 3D HEX8 plate at slenderness :math:`L/h = 20` with three transverse layers; the ``enhanced_strain`` formulation (Simo-Rifai 1990) is selected to skirt the shear-locking that plain bilinear hexes suffer on slender bending. The UDL is distributed across the top-surface nodes as a per-node ``F_y`` proportional to the nodal influence area. .. GENERATED FROM PYTHON SOURCE LINES 45-54 .. code-block:: Python from __future__ import annotations import numpy as np import pyvista as pv import femorph_solver from femorph_solver import ELEMENTS .. GENERATED FROM PYTHON SOURCE LINES 55-57 Geometry + material ------------------- .. GENERATED FROM PYTHON SOURCE LINES 57-78 .. code-block:: Python L = 1.0 WIDTH = 0.05 HEIGHT = 0.05 E = 2.0e11 NU = 0.30 RHO = 7850.0 # Total load per unit length [N/m]; integrated over L gives the # total force for the equivalent point-load comparison. w = 1000.0 W_total = w * L I = WIDTH * HEIGHT**3 / 12.0 # noqa: E741 # Closed-form magnitude is wL^4 / (8EI); sign convention here is # negative because the applied UDL points in -y. delta_published = -w * L**4 / (8.0 * E * I) print(f"problem: L={L} m, w={w} N/m (downward), EI={E * I:.4e} N·m²") print(f"published δ_tip = -w L^4 / (8 E I) = {delta_published:+.6e} m") .. rst-class:: sphx-glr-script-out .. code-block:: none problem: L=1.0 m, w=1000.0 N/m (downward), EI=1.0417e+05 N·m² published δ_tip = -w L^4 / (8 E I) = -1.200000e-03 m .. GENERATED FROM PYTHON SOURCE LINES 79-81 Build a HEX8 plate cantilever --------------------------------- .. GENERATED FROM PYTHON SOURCE LINES 81-98 .. code-block:: Python NX, NY, NZ = 40, 3, 3 xs = np.linspace(0.0, L, NX + 1) ys = np.linspace(0.0, WIDTH, NY + 1) zs = np.linspace(0.0, 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": E, "PRXY": NU, "DENS": RHO}, ) pts = np.asarray(m.grid.points) clamped = np.where(pts[:, 0] < 1e-9)[0] m.fix(nodes=(clamped + 1).tolist(), dof="ALL") .. GENERATED FROM PYTHON SOURCE LINES 99-107 Distribute the UDL as nodal F_y on the top surface -------------------------------------------------- A uniform line load ``w`` (force/length) on the beam axis maps to a per-node :math:`F_y` proportional to the nodal influence area on the top surface. For a regular structured grid this is constant per interior node and half at the edges. Total distributed force = ``w * L`` is invariant. .. GENERATED FROM PYTHON SOURCE LINES 107-115 .. code-block:: Python top = np.where(pts[:, 2] > HEIGHT - 1e-9)[0] # Equal weights per top-surface node — coarse but acceptable on # a structured grid for this verification. fy_per_node = -W_total / len(top) for n in top: m.apply_force(int(n + 1), fy=fy_per_node) .. GENERATED FROM PYTHON SOURCE LINES 116-118 Solve and compare ----------------- .. GENERATED FROM PYTHON SOURCE LINES 118-134 .. code-block:: Python res = m.solve() u = np.asarray(res.displacement).reshape(-1, 3) tip = np.where(pts[:, 0] > L - 1e-9)[0] delta_computed = float(u[tip, 1].mean()) rel_err = (delta_computed - delta_published) / delta_published print(f"computed δ_tip = {delta_computed:+.6e} m") print(f"relative error vs textbook = {rel_err * 100:+.2f} %") # 10 % tolerance at this mesh — HEX8 EAS recovers the # Euler-Bernoulli answer with ~5-8 % residual at L/h=20 on a # 40x3x3 grid. assert abs(rel_err) < 0.10, f"cantilever UDL deflection failed: {rel_err:.2%}" .. rst-class:: sphx-glr-script-out .. code-block:: none computed δ_tip = -1.202263e-03 m relative error vs textbook = +0.19 % .. GENERATED FROM PYTHON SOURCE LINES 135-137 Render the deflected shape -------------------------- .. GENERATED FROM PYTHON SOURCE LINES 137-148 .. code-block:: Python warped = m.grid.copy() warped.points = m.grid.points + u warped["|u|"] = np.linalg.norm(u, axis=1) plotter = pv.Plotter(off_screen=True, window_size=(640, 240)) plotter.add_mesh(m.grid, style="wireframe", color="grey", opacity=0.3) plotter.add_mesh(warped, scalars="|u|", cmap="viridis", show_edges=False) plotter.view_xy() plotter.camera.zoom(1.05) plotter.show() .. image-sg:: /gallery/verification/images/sphx_glr_example_verify_cantilever_udl_001.png :alt: example verify cantilever udl :srcset: /gallery/verification/images/sphx_glr_example_verify_cantilever_udl_001.png :class: sphx-glr-single-img .. rst-class:: sphx-glr-timing **Total running time of the script:** (0 minutes 0.428 seconds) .. _sphx_glr_download_gallery_verification_example_verify_cantilever_udl.py: .. only:: html .. container:: sphx-glr-footer sphx-glr-footer-example .. container:: sphx-glr-download sphx-glr-download-jupyter :download:`Download Jupyter notebook: example_verify_cantilever_udl.ipynb ` .. container:: sphx-glr-download sphx-glr-download-python :download:`Download Python source code: example_verify_cantilever_udl.py ` .. container:: sphx-glr-download sphx-glr-download-zip :download:`Download zipped: example_verify_cantilever_udl.zip ` .. only:: html .. rst-class:: sphx-glr-signature `Gallery generated by Sphinx-Gallery `_