.. DO NOT EDIT. .. THIS FILE WAS AUTOMATICALLY GENERATED BY SPHINX-GALLERY. .. TO MAKE CHANGES, EDIT THE SOURCE PYTHON FILE: .. "gallery/elements/solid185/example_solid185.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_elements_solid185_example_solid185.py: .. _ref_solid185_example: HEX8 — uniaxial tension of a unit cube ========================================== Single HEX8 element (a 1 m unit cube) with one face clamped in the pulling direction and the opposite face loaded. The displacement field is uniform — ``u_x = σ L / E`` on the loaded face — which gives a clean check on the full 24 × 24 stiffness plus the boundary/load handling. .. GENERATED FROM PYTHON SOURCE LINES 12-22 .. code-block:: Python from __future__ import annotations import numpy as np import pyvista as pv from vtkmodules.util.vtkConstants import VTK_HEXAHEDRON import femorph_solver from femorph_solver import ELEMENTS .. GENERATED FROM PYTHON SOURCE LINES 23-27 Problem data ------------ Steel. Unit cube ⇒ cross-section area ``A = 1`` and length ``L = 1``, so applying a total force ``F`` produces nominal stress ``σ = F``. .. GENERATED FROM PYTHON SOURCE LINES 27-32 .. code-block:: Python E = 2.1e11 # Pa NU = 0.3 RHO = 7850.0 F_TOTAL = 1.0e6 # N distributed over the +x face (4 corner nodes) .. GENERATED FROM PYTHON SOURCE LINES 33-38 Build the model --------------- Eight corners in VTK_HEXAHEDRON order. Clamp the -x face (nodes 1, 4, 5, 8) in x and pin one node's y/z to prevent rigid-body rotation. Split the total +x force onto the four +x-face nodes. .. GENERATED FROM PYTHON SOURCE LINES 38-73 .. code-block:: Python coords = np.array( [ [0.0, 0.0, 0.0], # 1 [1.0, 0.0, 0.0], # 2 [1.0, 1.0, 0.0], # 3 [0.0, 1.0, 0.0], # 4 [0.0, 0.0, 1.0], # 5 [1.0, 0.0, 1.0], # 6 [1.0, 1.0, 1.0], # 7 [0.0, 1.0, 1.0], # 8 ], dtype=np.float64, ) cells = np.array([8, 0, 1, 2, 3, 4, 5, 6, 7], dtype=np.int64) cell_types = np.array([VTK_HEXAHEDRON], dtype=np.uint8) grid = pv.UnstructuredGrid(cells, cell_types, coords) m = femorph_solver.Model.from_grid(grid) m.assign(ELEMENTS.HEX8, material={"EX": E, "PRXY": NU, "DENS": RHO}) # Clamp x=0 face (nodes 1, 4, 5, 8) in UX for nn in (1, 4, 5, 8): m.fix(nodes=[nn], dof="UX") # Roller supports to kill rigid-body y/z motion without restraining Poisson for nn in (1, 2, 5, 6): m.fix(nodes=[nn], dof="UY") for nn in (1, 2, 3, 4): m.fix(nodes=[nn], dof="UZ") # Apply F_TOTAL on x=1 face, split evenly over its 4 corners F_each = F_TOTAL / 4.0 for nn in (2, 3, 6, 7): m.apply_force(nn, fx=F_each) .. GENERATED FROM PYTHON SOURCE LINES 74-80 Static solve + analytical comparison ------------------------------------ Expected: ``u_x = σ L / E`` with ``σ = F_TOTAL / A = F_TOTAL`` and ``L = 1`` → ``u_x = F_TOTAL / E`` on the +x face. Poisson contraction is ``u_y = u_z = −ν u_x`` but only the free nodes on those faces actually display it (here only UY at the +y face, UZ at the +z face). .. GENERATED FROM PYTHON SOURCE LINES 80-99 .. code-block:: Python res = m.solve() dof = m.dof_map() ux_expected = F_TOTAL / E uy_expected = -NU * ux_expected print(f"Expected u_x (+x face) = {ux_expected:.6e} m") for nn in (2, 3, 6, 7): row = np.where((dof[:, 0] == nn) & (dof[:, 1] == 0))[0][0] print(f" node {nn} UX = {res.displacement[row]:.6e}") assert np.isclose(res.displacement[row], ux_expected, rtol=1e-8) print(f"Expected u_y (+y face) = {uy_expected:.6e} m") for nn in (3, 4, 7, 8): row = np.where((dof[:, 0] == nn) & (dof[:, 1] == 1))[0][0] val = res.displacement[row] print(f" node {nn} UY = {val:.6e}") assert np.isclose(val, uy_expected, rtol=1e-6) .. rst-class:: sphx-glr-script-out .. code-block:: none Expected u_x (+x face) = 4.761905e-06 m node 2 UX = 4.761905e-06 node 3 UX = 4.761905e-06 node 6 UX = 4.761905e-06 node 7 UX = 4.761905e-06 Expected u_y (+y face) = -1.428571e-06 m node 3 UY = -1.428571e-06 node 4 UY = -1.428571e-06 node 7 UY = -1.428571e-06 node 8 UY = -1.428571e-06 .. GENERATED FROM PYTHON SOURCE LINES 100-102 Plot the deformed cube, coloured by displacement magnitude ---------------------------------------------------------- .. GENERATED FROM PYTHON SOURCE LINES 102-120 .. code-block:: Python grid = m.grid.copy() disp = np.zeros((grid.n_points, 3), dtype=np.float64) for i, nn in enumerate(grid.point_data["ansys_node_num"]): rows = np.where(dof[:, 0] == int(nn))[0] for r in rows: disp[i, int(dof[r, 1])] = res.displacement[r] grid.point_data["displacement"] = disp plotter = pv.Plotter(off_screen=True) plotter.add_mesh(grid, style="wireframe", color="gray", line_width=2) plotter.add_mesh( grid.warp_by_vector("displacement", factor=5.0e4), scalars=np.linalg.norm(disp, axis=1), show_edges=True, scalar_bar_args={"title": "disp [m]"}, ) plotter.add_axes() plotter.show() .. image-sg:: /gallery/elements/solid185/images/sphx_glr_example_solid185_001.png :alt: example solid185 :srcset: /gallery/elements/solid185/images/sphx_glr_example_solid185_001.png :class: sphx-glr-single-img .. rst-class:: sphx-glr-timing **Total running time of the script:** (0 minutes 0.207 seconds) .. _sphx_glr_download_gallery_elements_solid185_example_solid185.py: .. only:: html .. container:: sphx-glr-footer sphx-glr-footer-example .. container:: sphx-glr-download sphx-glr-download-jupyter :download:`Download Jupyter notebook: example_solid185.ipynb ` .. container:: sphx-glr-download sphx-glr-download-python :download:`Download Python source code: example_solid185.py ` .. container:: sphx-glr-download sphx-glr-download-zip :download:`Download zipped: example_solid185.zip ` .. only:: html .. rst-class:: sphx-glr-signature `Gallery generated by Sphinx-Gallery `_