.. DO NOT EDIT.
.. THIS FILE WAS AUTOMATICALLY GENERATED BY SPHINX-GALLERY.
.. TO MAKE CHANGES, EDIT THE SOURCE PYTHON FILE:
.. "examples/LinearizedElasticity/Elas8.py"
.. LINE NUMBERS ARE GIVEN BELOW.

.. only:: html

    .. note::
        :class: sphx-glr-download-link-note

        :ref:`Go to the end <sphx_glr_download_examples_LinearizedElasticity_Elas8.py>`
        to download the full example code.

.. rst-class:: sphx-glr-example-title

.. _sphx_glr_examples_LinearizedElasticity_Elas8.py:


Elas8
=====

A cantilever cylinder undergoing torsion, loaded either by a surface traction.

.. GENERATED FROM PYTHON SOURCE LINES 12-93







.. tab-set::



   .. tab-item:: Static Scene



            
     .. image-sg:: /examples/LinearizedElasticity/images/sphx_glr_Elas8_001.png
        :alt: Elas8
        :srcset: /examples/LinearizedElasticity/images/sphx_glr_Elas8_001.png
        :class: sphx-glr-single-img
     


   .. tab-item:: Interactive Scene



       .. offlineviewer:: /home/docs/checkouts/readthedocs.org/user_builds/easyfea/checkouts/v3.1.0/docs/examples/LinearizedElasticity/images/sphx_glr_Elas8_001.vtksz






.. tab-set::



   .. tab-item:: Static Scene



            
     .. image-sg:: /examples/LinearizedElasticity/images/sphx_glr_Elas8_002.png
        :alt: Elas8
        :srcset: /examples/LinearizedElasticity/images/sphx_glr_Elas8_002.png
        :class: sphx-glr-single-img
     


   .. tab-item:: Interactive Scene



       .. offlineviewer:: /home/docs/checkouts/readthedocs.org/user_builds/easyfea/checkouts/v3.1.0/docs/examples/LinearizedElasticity/images/sphx_glr_Elas8_002.vtksz






.. tab-set::



   .. tab-item:: Static Scene



            
     .. image-sg:: /examples/LinearizedElasticity/images/sphx_glr_Elas8_003.png
        :alt: Elas8
        :srcset: /examples/LinearizedElasticity/images/sphx_glr_Elas8_003.png
        :class: sphx-glr-single-img
     


   .. tab-item:: Interactive Scene



       .. offlineviewer:: /home/docs/checkouts/readthedocs.org/user_builds/easyfea/checkouts/v3.1.0/docs/examples/LinearizedElasticity/images/sphx_glr_Elas8_003.vtksz



.. rst-class:: sphx-glr-script-out

 .. code-block:: none

    
    phi analytical = 1.556930e-07 rad
    phi numerical  = 1.556930e-07 rad
    relative error = 0.000 %





|

.. code-block:: Python
   :lineno-start: 13


    import numpy as np

    from EasyFEA import Terminal, ElemType, Models, Simulations, PyVista
    from EasyFEA.Geoms import Circle

    if __name__ == "__main__":
        Terminal.Clear()

        # --------------------------------------------
        # Mesh
        # --------------------------------------------

        L = 100  # mm
        D = 30

        circle = Circle((0, 0), D, n=(1, 0, 0))
        circle.meshSize = circle.length / 20

        mesh = circle.Mesh_Extrude([], [L, 0, 0], [L // circle.meshSize], ElemType.PRISM6)

        PyVista.Plot_Mesh(mesh).show()

        # --------------------------------------------
        # Simu
        # --------------------------------------------

        model = Models.Elastic.Isotropic(3, E=210000, v=0.3)
        simu = Simulations.Elastic(mesh, model)

        nodesX0 = mesh.Nodes_Conditions(lambda x, y, z: x == 0)
        simu.add_dirichlet(nodesX0, [0] * 3, ["x", "y", "z"])

        nodesXL = mesh.Nodes_Conditions(lambda x, y, z: x == L)
        yC, zC = 0, 0
        T = 10  # N.mm
        J = np.pi * D**4 / 32  # polar second moment of area
        p = T / J  # so that the traction t = p * r integrates to a moment M_x = c
        simu.add_surfLoad(
            nodesXL,
            [
                lambda x, y, z: -p * (z - zC),
                lambda x, y, z: p * (y - yC),
            ],
            ["y", "z"],
        )

        u = simu.Solve()

        # --------------------------------------------
        # Results
        # --------------------------------------------
        # Saint-Venant torsion of a circular shaft: phi = T * L / (G * J)
        phi_ana = T * L / (model.get_mu() * J)  # rad

        # numerical twist from the tip-face nodes: a rigid rotation theta about x
        # gives uy = -theta*dz, uz = theta*dy, so theta = (dy*uz - dz*uy) / r^2
        uy = simu.Result("uy")[nodesXL]
        uz = simu.Result("uz")[nodesXL]
        dy = mesh.coord[nodesXL, 1] - yC
        dz = mesh.coord[nodesXL, 2] - zC
        r2 = dy**2 + dz**2
        phi_num = np.mean((dy * uz - dz * uy) / r2, where=r2 > 2)

        err = abs(phi_num - phi_ana) / abs(phi_ana) * 100
        Terminal.MyPrint(f"\nphi analytical = {phi_ana:.6e} rad")
        Terminal.MyPrint(f"\nphi numerical  = {phi_num:.6e} rad")
        Terminal.MyPrint(f"\nrelative error = {err:.3f} %")

        PyVista.Plot_BoundaryConditions(simu).show()

        pltr = PyVista._Plotter(shape=(3, 1))

        PyVista.Plot(simu, "uy", plotMesh=True, plotter=pltr, verticalColobar=False)
        pltr.subplot(1, 0)
        PyVista.Plot(simu, "uz", plotMesh=True, plotter=pltr, verticalColobar=False)
        pltr.subplot(2, 0)
        PyVista.Plot(
            simu, "displacement_norm", plotMesh=True, plotter=pltr, verticalColobar=False
        )
        pltr.show()


.. rst-class:: sphx-glr-timing

   **Total running time of the script:** (0 minutes 0.867 seconds)


.. _sphx_glr_download_examples_LinearizedElasticity_Elas8.py:

.. only:: html

  .. container:: sphx-glr-footer sphx-glr-footer-example

    .. container:: sphx-glr-download sphx-glr-download-jupyter

      :download:`Download Jupyter notebook: Elas8.ipynb <Elas8.ipynb>`

    .. container:: sphx-glr-download sphx-glr-download-python

      :download:`Download Python source code: Elas8.py <Elas8.py>`

    .. container:: sphx-glr-download sphx-glr-download-zip

      :download:`Download zipped: Elas8.zip <Elas8.zip>`


.. only:: html

 .. rst-class:: sphx-glr-signature

    `Gallery generated by Sphinx-Gallery <https://sphinx-gallery.github.io>`_