.. DO NOT EDIT. .. THIS FILE WAS AUTOMATICALLY GENERATED BY SPHINX-GALLERY. .. TO MAKE CHANGES, EDIT THE SOURCE PYTHON FILE: .. "gallery/elements/shell181/example_shell181.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_shell181_example_shell181.py: .. _ref_shell181_example: QUAD4_SHELL — uniaxial stretch of a flat square plate ===================================================== Single QUAD4_SHELL element (a 1 m × 1 m flat plate of thickness ``t``) in pure in-plane tension. Out-of-plane DOFs are suppressed so the problem reduces to plane-stress membrane behaviour: ``u_x = σ L / E`` and the Poisson contraction ``u_y = −ν u_x``. .. 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_QUAD import femorph_solver from femorph_solver import ELEMENTS .. GENERATED FROM PYTHON SOURCE LINES 23-26 Problem data ------------ Steel plate, 10 mm thick, pulled with a total 100 kN along +x. .. GENERATED FROM PYTHON SOURCE LINES 26-33 .. code-block:: Python E = 2.1e11 NU = 0.3 RHO = 7850.0 THK = 0.010 # 10 mm L = 1.0 F_TOTAL = 1.0e5 # N spread over the +x edge (two corner nodes) .. GENERATED FROM PYTHON SOURCE LINES 34-41 Build the model --------------- Corner ordering is VTK_QUAD: node 1 = (0,0), 2 = (1,0), 3 = (1,1), 4 = (0,1). We clamp UX on the -x edge, anchor UY on a single corner to kill the in-plane translation mode, and suppress every out-of-plane DOF (UZ, ROTX, ROTY, ROTZ) so the element behaves like a pure membrane. .. GENERATED FROM PYTHON SOURCE LINES 41-78 .. code-block:: Python coords = np.array( [ [0.0, 0.0, 0.0], [L, 0.0, 0.0], [L, L, 0.0], [0.0, L, 0.0], ], dtype=np.float64, ) cells = np.array([4, 0, 1, 2, 3], dtype=np.int64) cell_types = np.array([VTK_QUAD], dtype=np.uint8) grid = pv.UnstructuredGrid(cells, cell_types, coords) m = femorph_solver.Model.from_grid(grid) m.assign( ELEMENTS.QUAD4_SHELL, material={"EX": E, "PRXY": NU, "DENS": RHO}, real=(THK,), ) # -x edge clamped in UX m.fix(nodes=[1], dof="UX") m.fix(nodes=[4], dof="UX") # Anchor UY on one node to stop rigid-body in-plane translation m.fix(nodes=[1], dof="UY") # Pure membrane: zero all out-of-plane DOFs everywhere for nn in (1, 2, 3, 4): m.fix(nodes=[nn], dof="UZ") m.fix(nodes=[nn], dof="ROTX") m.fix(nodes=[nn], dof="ROTY") m.fix(nodes=[nn], dof="ROTZ") F_each = F_TOTAL / 2.0 m.apply_force(2, fx=F_each) m.apply_force(3, fx=F_each) .. GENERATED FROM PYTHON SOURCE LINES 79-84 Static solve + analytical comparison ------------------------------------ Stress on the +x face is ``σ = F / (width · t)``. On a 1 m width with 10 mm thickness: ``σ = 10 MPa``. Expected ``u_x = σ L / E`` and ``u_y = −ν u_x`` on the free corner nodes. .. GENERATED FROM PYTHON SOURCE LINES 84-106 .. code-block:: Python res = m.solve() dof = m.dof_map() sigma = F_TOTAL / (L * THK) ux_expected = sigma * L / E uy_expected = -NU * ux_expected print(f"In-plane stress σ = {sigma:.3e} Pa") print(f"Expected u_x = {ux_expected:.6e} m") for nn in (2, 3): row = np.where((dof[:, 0] == nn) & (dof[:, 1] == 0))[0][0] val = res.displacement[row] print(f" node {nn} UX = {val:.6e}") assert np.isclose(val, ux_expected, rtol=1e-8) # Node 3 is free in UY (node 1 is anchored, so check there) row_uy3 = np.where((dof[:, 0] == 3) & (dof[:, 1] == 1))[0][0] uy3 = res.displacement[row_uy3] print(f"Expected u_y (node 3) = {uy_expected:.6e}") print(f"Computed u_y (node 3) = {uy3:.6e}") assert np.isclose(uy3, uy_expected, rtol=1e-6) .. rst-class:: sphx-glr-script-out .. code-block:: none In-plane stress σ = 1.000e+07 Pa Expected u_x = 4.761905e-05 m node 2 UX = 4.761905e-05 node 3 UX = 4.761905e-05 Expected u_y (node 3) = -1.428571e-05 Computed u_y (node 3) = -1.428571e-05 .. GENERATED FROM PYTHON SOURCE LINES 107-111 Plot the deformed plate ----------------------- Scatter displacements onto point_data and warp the mesh. A modest magnification (~1000×) makes the Poisson contraction visible. .. GENERATED FROM PYTHON SOURCE LINES 111-131 .. 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: d_idx = int(dof[r, 1]) if d_idx < 3: disp[i, d_idx] = 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=2.0e3), scalars=np.linalg.norm(disp, axis=1), show_edges=True, scalar_bar_args={"title": "|u| [m]"}, ) plotter.add_axes() plotter.show() .. image-sg:: /gallery/elements/shell181/images/sphx_glr_example_shell181_001.png :alt: example shell181 :srcset: /gallery/elements/shell181/images/sphx_glr_example_shell181_001.png :class: sphx-glr-single-img .. rst-class:: sphx-glr-timing **Total running time of the script:** (0 minutes 0.194 seconds) .. _sphx_glr_download_gallery_elements_shell181_example_shell181.py: .. only:: html .. container:: sphx-glr-footer sphx-glr-footer-example .. container:: sphx-glr-download sphx-glr-download-jupyter :download:`Download Jupyter notebook: example_shell181.ipynb ` .. container:: sphx-glr-download sphx-glr-download-python :download:`Download Python source code: example_shell181.py ` .. container:: sphx-glr-download sphx-glr-download-zip :download:`Download zipped: example_shell181.zip ` .. only:: html .. rst-class:: sphx-glr-signature `Gallery generated by Sphinx-Gallery `_