.. DO NOT EDIT. .. THIS FILE WAS AUTOMATICALLY GENERATED BY SPHINX-GALLERY. .. TO MAKE CHANGES, EDIT THE SOURCE PYTHON FILE: .. "examples/ex02_plate-with-hole.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_examples_ex02_plate-with-hole.py: Plate with a Hole ----------------- .. topic:: Plane stress linear analysis. * create and mesh a plate with a hole * define a solid body with a linear-elastic plane stress material * create an external pressure load * export and plot stress results A plate with length :math:`2L`, height :math:`2h` and a hole with radius :math:`r` is subjected to a uniaxial tension :math:`p=-100` MPa. What is being looked for is the von Mises stress distribution and the concentration of normal stress :math:`\sigma_{11}` over the hole. .. image:: ../../examples/ex02_plate-with-hole_sketch.svg :width: 400px .. GENERATED FROM PYTHON SOURCE LINES 24-30 .. code-block:: Python h = 1 L = 2 r = 0.3 import numpy as np .. GENERATED FROM PYTHON SOURCE LINES 32-36 Let's create a meshed plate with a hole out of quad cells. Only a quarter model of the plate is considered. The mesh generation is carried out by filling the area between the edge of the hole and the top line. Then, this section is duplicated, mirrored and merged with another rectangle. .. GENERATED FROM PYTHON SOURCE LINES 36-51 .. code-block:: Python import felupe as fem phi = np.linspace(1, 0.5, 21) * np.pi / 2 line = fem.mesh.Line(n=21) curve = line.copy(points=r * np.vstack([np.cos(phi), np.sin(phi)]).T) top = line.copy(points=np.vstack([np.linspace(0, h, 21), np.linspace(h, h, 21)]).T) face = curve.fill_between(top, n=np.linspace(0, 1, 21) ** 1.3 * 2 - 1) rect = fem.mesh.Rectangle(a=(h, 0), b=(L, h), n=21) mesh = fem.mesh.concatenate([face, face.mirror(normal=[-1, 1, 0]), rect]) mesh = mesh.sweep(decimals=5) mesh.plot().show() .. tab-set:: .. tab-item:: Static Scene .. image-sg:: /examples/images/sphx_glr_ex02_plate-with-hole_001.png :alt: ex02 plate with hole :srcset: /examples/images/sphx_glr_ex02_plate-with-hole_001.png :class: sphx-glr-single-img .. tab-item:: Interactive Scene .. offlineviewer:: /home/docs/checkouts/readthedocs.org/user_builds/felupe/checkouts/latest/docs/examples/images/sphx_glr_ex02_plate-with-hole_001.vtksz .. GENERATED FROM PYTHON SOURCE LINES 52-55 A numeric quad-region created on the mesh in combination with a vector-valued displacement field represents the plate. The Boundary conditions for the symmetry planes are generated on the displacement field. .. GENERATED FROM PYTHON SOURCE LINES 55-61 .. code-block:: Python region = fem.RegionQuad(mesh) displacement = fem.Field(region, dim=2) field = fem.FieldContainer([displacement]) boundaries = fem.dof.symmetry(displacement) .. GENERATED FROM PYTHON SOURCE LINES 62-65 The material behaviour is defined through a built-in isotropic linear-elastic material formulation for plane stress problems. A solid body applies the linear-elastic material formulation on the displacement field. .. GENERATED FROM PYTHON SOURCE LINES 65-68 .. code-block:: Python umat = fem.LinearElasticPlaneStress(E=210000, nu=0.3) solid = fem.SolidBody(umat, field) .. GENERATED FROM PYTHON SOURCE LINES 69-72 The external uniaxial tension is applied by a pressure load on the right end at :math:`x=L`. Therefore, a boundary region in combination with a field has to be created at :math:`x=L`. .. GENERATED FROM PYTHON SOURCE LINES 72-77 .. code-block:: Python region_boundary = fem.RegionQuadBoundary(mesh, mask=mesh.points[:, 0] == L) field_boundary = fem.FieldContainer([fem.Field(region_boundary, dim=2)]) load = fem.SolidBodyPressure(field_boundary, pressure=-100) .. GENERATED FROM PYTHON SOURCE LINES 78-88 The simulation model is now ready to be solved. The equivalent von Mises Cauchy stress will be plotted. For the two-dimensional case it is calculated by Eq. :eq:`svM`. Stresses, located at quadrature-points of cells, are shifted to and averaged at mesh- points. .. math:: :label: svM \sigma_{vM} = \sqrt{\sigma_{11}^2 + \sigma_{22}^2 - \sigma_{11} \ \sigma_{22} + 3 \ \sigma_{12}^2} .. GENERATED FROM PYTHON SOURCE LINES 88-93 .. code-block:: Python step = fem.Step(items=[solid, load], boundaries=boundaries) job = fem.Job(steps=[step]).evaluate() solid.plot("Equivalent of Cauchy Stress", show_edges=False, project=fem.topoints).show() .. tab-set:: .. tab-item:: Static Scene .. image-sg:: /examples/images/sphx_glr_ex02_plate-with-hole_002.png :alt: ex02 plate with hole :srcset: /examples/images/sphx_glr_ex02_plate-with-hole_002.png :class: sphx-glr-single-img .. tab-item:: Interactive Scene .. offlineviewer:: /home/docs/checkouts/readthedocs.org/user_builds/felupe/checkouts/latest/docs/examples/images/sphx_glr_ex02_plate-with-hole_002.vtksz .. rst-class:: sphx-glr-script-out .. code-block:: none /home/docs/checkouts/readthedocs.org/user_builds/felupe/envs/latest/lib/python3.12/site-packages/felupe/mechanics/_solidbody.py:345: UserWarning: Cauchy stress tensor can't be evaluated on a 2d-Field. Falling-back to the Kirchhoff stress tensor. warnings.warn( .. GENERATED FROM PYTHON SOURCE LINES 94-96 The normal stress distribution :math:`\sigma_{11}` over the hole at :math:`x=0` is plotted with matplotlib. .. GENERATED FROM PYTHON SOURCE LINES 96-106 .. code-block:: Python import matplotlib.pyplot as plt plt.plot( fem.tools.project(solid.evaluate.cauchy_stress(), region)[:, 0, 0][mesh.x == 0], mesh.points[:, 1][mesh.x == 0] / h, "o-", ) plt.xlabel(r"$\sigma_{11}(x=0, y)$ in MPa $\longrightarrow$") plt.ylabel(r"$y/h$ $\longrightarrow$") .. image-sg:: /examples/images/sphx_glr_ex02_plate-with-hole_003.png :alt: ex02 plate with hole :srcset: /examples/images/sphx_glr_ex02_plate-with-hole_003.png :class: sphx-glr-single-img .. rst-class:: sphx-glr-script-out .. code-block:: none /home/docs/checkouts/readthedocs.org/user_builds/felupe/envs/latest/lib/python3.12/site-packages/felupe/mechanics/_solidbody.py:345: UserWarning: Cauchy stress tensor can't be evaluated on a 2d-Field. Falling-back to the Kirchhoff stress tensor. warnings.warn( .. rst-class:: sphx-glr-timing **Total running time of the script:** (0 minutes 1.850 seconds) .. _sphx_glr_download_examples_ex02_plate-with-hole.py: .. only:: html .. container:: sphx-glr-footer sphx-glr-footer-example .. container:: sphx-glr-download sphx-glr-download-jupyter :download:`Download Jupyter notebook: ex02_plate-with-hole.ipynb ` .. container:: sphx-glr-download sphx-glr-download-python :download:`Download Python source code: ex02_plate-with-hole.py ` .. only:: html .. rst-class:: sphx-glr-signature `Gallery generated by Sphinx-Gallery `_