.. DO NOT EDIT. .. THIS FILE WAS AUTOMATICALLY GENERATED BY SPHINX-GALLERY. .. TO MAKE CHANGES, EDIT THE SOURCE PYTHON FILE: .. "verif-manual/vm-002-beam_stresses_and_deflections.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_verif-manual_vm-002-beam_stresses_and_deflections.py: .. _ref_vm2_example: Beam stresses and deflections ----------------------------- Problem description: - A standard 30 inch WF beam, with a cross-sectional area :math:`A`, is supported as shown below and loaded on the overhangs by a uniformly distributed load :math:`w`. Determine the maximum bending stress, :math:`\sigma_max`, in the middle portion of the beam and the deflection, :math:`\delta`, at the middle of the beam. Reference: - S. Timoshenko, Strength of Material, Part I, Elementary Theory and Problems, 3rd Edition, D. Van Nostrand Co., Inc., New York, NY, 1955, pg. 98, problem 4. Analysis type(s): - Static Analysis ``ANTYPE=0`` Element type(s): - 3-D 2 Node Beam (BEAM188) .. image:: ../_static/vm2_setup.png :width: 400 :alt: VM2 Problem Sketch Material properties: - :math:`E = 30 \cdot 10^6 psi` Geometric properties: - :math:`a = 120 in` - :math:`l = 240 in` - :math:`h = 30 in` - :math:`A = 50.65 in^2` - :math:`I_z = 7892 in^4` Loading: - :math:`w = (10000/12) lb/in` Analytical equations: - :math:`M` is the bending moment for the middle portion of the beam: :math:`M = 10000 \cdot 10 \cdot 60 = 6 \cdot 10^6 lb \cdot in` - Determination of the maximum stress in the middle portion of the beam is :math:`\sigma_max = \frac{M h}{2 I_z}` - The deflection, :math:`\delta`, at the middle of the beam can be defined by the formulas of the transversally loaded beam: :math:`\delta = 0.182 in` .. GENERATED FROM PYTHON SOURCE LINES 80-82 .. code-block:: Python # sphinx_gallery_thumbnail_path = '_static/vm2_setup.png' .. GENERATED FROM PYTHON SOURCE LINES 83-85 Start MAPDL ~~~~~~~~~~~ .. GENERATED FROM PYTHON SOURCE LINES 85-97 .. code-block:: Python from ansys.mapdl.core import launch_mapdl # Start mapdl and clear it. mapdl = launch_mapdl() mapdl.clear() # Enter verification example mode and the pre-processing routine. mapdl.verify() mapdl.prep7() .. rst-class:: sphx-glr-script-out .. code-block:: none *****MAPDL VERIFICATION RUN ONLY***** DO NOT USE RESULTS FOR PRODUCTION ***** MAPDL ANALYSIS DEFINITION (PREP7) ***** .. GENERATED FROM PYTHON SOURCE LINES 98-101 Define element type ~~~~~~~~~~~~~~~~~~~ Set up the element type (a beam-type). .. GENERATED FROM PYTHON SOURCE LINES 101-122 .. code-block:: Python # Type of analysis: static. mapdl.antype("STATIC") # Element type: BEAM188. mapdl.et(1, "BEAM188") # Special Features are defined by keyoptions of beam element: # KEYOPT(3) # Shape functions along the length: # Cubic mapdl.keyopt(1, 3, 3) # Cubic shape function # KEYOPT(9) # Output control for values extrapolated to the element # and section nodes: # Same as KEYOPT(9) = 1 plus stresses and strains at all section nodes mapdl.keyopt(1, 9, 3, mute=True) .. GENERATED FROM PYTHON SOURCE LINES 123-126 Define material ~~~~~~~~~~~~~~~ Set up the material. .. GENERATED FROM PYTHON SOURCE LINES 126-132 .. code-block:: Python mapdl.mp("EX", 1, 30e6) mapdl.mp("PRXY", 1, 0.3) print(mapdl.mplist()) .. rst-class:: sphx-glr-script-out .. code-block:: none LIST MATERIALS 1 TO 1 BY 1 PROPERTY= ALL MATERIAL NUMBER 1 TEMP EX 0.3000000E+08 TEMP PRXY 0.3000000 .. GENERATED FROM PYTHON SOURCE LINES 133-136 Define section ~~~~~~~~~~~~~~ Set up the cross-section properties for a beam element. .. GENERATED FROM PYTHON SOURCE LINES 136-144 .. code-block:: Python w_f = 1.048394965 w_w = 0.6856481 sec_num = 1 mapdl.sectype(sec_num, "BEAM", "I", "ISection") mapdl.secdata(15, 15, 28 + (2 * w_f), w_f, w_f, w_w) .. rst-class:: sphx-glr-script-out .. code-block:: none SECTION ID NUMBER IS: 1 BEAM SECTION TYPE IS: I Section BEAM SECTION NAME IS: ISection COMPUTED BEAM SECTION DATA SUMMARY: Area = 50.650 Iyy = 7892.0 Iyz =-0.78160E-12 Izz = 590.47 Warping Constant = 0.12403E+06 Torsion Constant = 14.962 Centroid Y =-0.11661E-14 Centroid Z = 15.048 Shear Center Y = 0.51867E-13 Shear Center Z = 15.048 Shear Correction-xy = 0.54626 Shear Correction-yz =-0.75912E-14 Shear Correction-xz = 0.38629 Beam Section is offset to CENTROID of cross section .. GENERATED FROM PYTHON SOURCE LINES 145-149 Define geometry ~~~~~~~~~~~~~~~ Set up the nodes and elements. Create nodes then create elements between nodes. .. GENERATED FROM PYTHON SOURCE LINES 149-160 .. code-block:: Python # Define nodes for node_num in range(1, 6): mapdl.n(node_num, (node_num - 1) * 120, 0, 0) # Define one node for the orientation of the beam cross-section. orient_node = mapdl.n(6, 60, 1) # Print the list of the created nodes. print(mapdl.nlist()) .. rst-class:: sphx-glr-script-out .. code-block:: none LIST ALL SELECTED NODES. DSYS= 0 *****MAPDL VERIFICATION RUN ONLY***** DO NOT USE RESULTS FOR PRODUCTION NODE X Y Z THXY THYZ THZX 1 0.0000 0.0000 0.0000 0.00 0.00 0.00 2 120.00 0.0000 0.0000 0.00 0.00 0.00 3 240.00 0.0000 0.0000 0.00 0.00 0.00 4 360.00 0.0000 0.0000 0.00 0.00 0.00 5 480.00 0.0000 0.0000 0.00 0.00 0.00 6 60.000 1.0000 0.0000 0.00 0.00 0.00 .. GENERATED FROM PYTHON SOURCE LINES 161-162 Define elements .. GENERATED FROM PYTHON SOURCE LINES 162-173 .. code-block:: Python for elem_num in range(1, 5): mapdl.e(elem_num, elem_num + 1, orient_node) # Print the list of the created elements. print(mapdl.elist()) # Display elements with their nodes numbers. mapdl.eplot(show_node_numbering=True, line_width=5, cpos="xy", font_size=40) .. tab-set:: .. tab-item:: Static Scene .. image-sg:: /verif-manual/images/sphx_glr_vm-002-beam_stresses_and_deflections_001.png :alt: vm 002 beam stresses and deflections :srcset: /verif-manual/images/sphx_glr_vm-002-beam_stresses_and_deflections_001.png :class: sphx-glr-single-img .. tab-item:: Interactive Scene .. offlineviewer:: /home/runner/work/pymapdl-examples/pymapdl-examples/doc/source/verif-manual/images/sphx_glr_vm-002-beam_stresses_and_deflections_001.vtksz .. rst-class:: sphx-glr-script-out .. code-block:: none LIST ALL SELECTED ELEMENTS. (LIST NODES) *****MAPDL VERIFICATION RUN ONLY***** DO NOT USE RESULTS FOR PRODUCTION ELEM MAT TYP REL ESY SEC NODES 1 1 1 1 0 1 1 2 6 2 1 1 1 0 1 2 3 6 3 1 1 1 0 1 3 4 6 4 1 1 1 0 1 4 5 6 .. GENERATED FROM PYTHON SOURCE LINES 174-177 Define boundary conditions ~~~~~~~~~~~~~~~~~~~~~~~~~~ Application of boundary conditions (BC). .. GENERATED FROM PYTHON SOURCE LINES 177-189 .. code-block:: Python # BC for the beams seats mapdl.d(2, "UX", lab2="UY") mapdl.d(4, "UY") # BC for all nodes of the beam mapdl.nsel("S", "LOC", "Y", 0) mapdl.d("ALL", "UZ") mapdl.d("ALL", "ROTX") mapdl.d("ALL", "ROTY") mapdl.nsel("ALL") .. rst-class:: sphx-glr-script-out .. code-block:: none array([1, 2, 3, 4, 5, 6], dtype=int32) .. GENERATED FROM PYTHON SOURCE LINES 190-194 Define distributed loads ~~~~~~~~~~~~~~~~~~~~~~~~ Apply a distributed force of :math:`w = (10000/12) lb/in` in the y-direction. .. GENERATED FROM PYTHON SOURCE LINES 194-204 .. code-block:: Python # Parametrization of the distributed load. w = 10000 / 12 # Application of the surface load to the beam element. mapdl.sfbeam(1, 1, "PRES", w) mapdl.sfbeam(4, 1, "PRES", w) mapdl.finish() .. rst-class:: sphx-glr-script-out .. code-block:: none ***** ROUTINE COMPLETED ***** CP = 0.000 .. GENERATED FROM PYTHON SOURCE LINES 205-208 Solve ~~~~~ Enter solution mode and solve the system. Print the solver output. .. GENERATED FROM PYTHON SOURCE LINES 208-215 .. code-block:: Python mapdl.run("/SOLU") out = mapdl.solve() mapdl.finish() print(out) .. rst-class:: sphx-glr-script-out .. code-block:: none ***** MAPDL SOLVE COMMAND ***** *** NOTE *** CP = 0.000 TIME= 00:00:00 There is no title defined for this analysis. *** SELECTION OF ELEMENT TECHNOLOGIES FOR APPLICABLE ELEMENTS *** ---GIVE SUGGESTIONS ONLY--- ELEMENT TYPE 1 IS BEAM188 . KEYOPT(1)=1 IS SUGGESTED FOR NON-CIRCULAR CROSS SECTIONS AND KEYOPT(3)=2 IS ALWAYS SUGGESTED. ELEMENT TYPE 1 IS BEAM188 . KEYOPT(15) IS ALREADY SET AS SUGGESTED. *****MAPDL VERIFICATION RUN ONLY***** DO NOT USE RESULTS FOR PRODUCTION S O L U T I O N O P T I O N S PROBLEM DIMENSIONALITY. . . . . . . . . . . . .3-D DEGREES OF FREEDOM. . . . . . UX UY UZ ROTX ROTY ROTZ ANALYSIS TYPE . . . . . . . . . . . . . . . . .STATIC (STEADY-STATE) GLOBALLY ASSEMBLED MATRIX . . . . . . . . . . .SYMMETRIC *** NOTE *** CP = 0.000 TIME= 00:00:00 Present time 0 is less than or equal to the previous time. Time will default to 1. *** NOTE *** CP = 0.000 TIME= 00:00:00 The conditions for direct assembly have been met. No .emat or .erot files will be produced. *** NOTE *** CP = 0.000 TIME= 00:00:00 Internal nodes from 7 to 14 are created. 8 internal nodes are used for quadratic and/or cubic options of BEAM188, PIPE288, and/or SHELL208. D I S T R I B U T E D D O M A I N D E C O M P O S E R ...Number of elements: 4 ...Number of nodes: 14 ...Decompose to 0 CPU domains ...Element load balance ratio = 0.000 L O A D S T E P O P T I O N S LOAD STEP NUMBER. . . . . . . . . . . . . . . . 1 TIME AT END OF THE LOAD STEP. . . . . . . . . . 1.0000 NUMBER OF SUBSTEPS. . . . . . . . . . . . . . . 1 STEP CHANGE BOUNDARY CONDITIONS . . . . . . . . NO PRINT OUTPUT CONTROLS . . . . . . . . . . . . .NO PRINTOUT DATABASE OUTPUT CONTROLS. . . . . . . . . . . .ALL DATA WRITTEN FOR THE LAST SUBSTEP *** NOTE *** CP = 0.000 TIME= 00:00:00 Predictor is ON by default for structural elements with rotational degrees of freedom. Use the PRED,OFF command to turn the predictor OFF if it adversely affects the convergence. Range of element maximum matrix coefficients in global coordinates Maximum = 2.999405619E+10 at element 0. Minimum = 2.999405619E+10 at element 0. *** ELEMENT MATRIX FORMULATION TIMES TYPE NUMBER ENAME TOTAL CP AVE CP 1 4 BEAM188 0.000 0.000000 Time at end of element matrix formulation CP = 0. DISTRIBUTED SPARSE MATRIX DIRECT SOLVER. Number of equations = 60, Maximum wavefront = 0 Memory available (MB) = 0.0 , Memory required (MB) = 0.0 Distributed sparse solver maximum pivot= 0 at node 0 . Distributed sparse solver minimum pivot= 0 at node 0 . Distributed sparse solver minimum pivot in absolute value= 0 at node 0 . *** ELEMENT RESULT CALCULATION TIMES TYPE NUMBER ENAME TOTAL CP AVE CP 1 4 BEAM188 0.000 0.000000 *** NODAL LOAD CALCULATION TIMES TYPE NUMBER ENAME TOTAL CP AVE CP 1 4 BEAM188 0.000 0.000000 *** LOAD STEP 1 SUBSTEP 1 COMPLETED. CUM ITER = 1 *** TIME = 1.00000 TIME INC = 1.00000 NEW TRIANG MATRIX .. GENERATED FROM PYTHON SOURCE LINES 216-221 Post-processing ~~~~~~~~~~~~~~~ Enter post-processing. To get the stress and deflection results from the middle node and cross-section of the beam we can use :meth:`Mapdl.get_value `. .. GENERATED FROM PYTHON SOURCE LINES 221-233 .. code-block:: Python # Enter the post-processing routine and select the first load step. mapdl.post1() mapdl.set(1) # Get the maximum stress at the middle of the beam. s_eqv_max = mapdl.get_value("secr", 2, "s", "eqv", "max") # Get the deflection at the middle of the beam. mid_node_uy = mapdl.get_value(entity="NODE", entnum=3, item1="u", it1num="y") .. GENERATED FROM PYTHON SOURCE LINES 234-239 Check results ~~~~~~~~~~~~~ Now that we have the results we can compare the nodal displacement and stress experienced by middle node of the beam to the known stresses -11,400 psi and 0.182 inches of the deflection. .. GENERATED FROM PYTHON SOURCE LINES 239-260 .. code-block:: Python # Results obtained by hand-calculations. stress_target = 11400.0 deflection_target = 0.182 # Calculate the deviation. stress_ratio = s_eqv_max / stress_target deflection_ratio = mid_node_uy / deflection_target # Print output results. output = f""" ----------------------------- VM3 RESULTS COMPARISON ----------------------------- | TARGET | Mechanical APDL | RATIO | ---------------------------------------------------------------------------------- Stress{stress_target:18.3f} {s_eqv_max:16.3f} {stress_ratio:14.3f} Deflection{deflection_target:14.3f} {mid_node_uy:16.3f} {deflection_ratio:14.3f} ---------------------------------------------------------------------------------- """ print(output) .. rst-class:: sphx-glr-script-out .. code-block:: none ----------------------------- VM3 RESULTS COMPARISON ----------------------------- | TARGET | Mechanical APDL | RATIO | ---------------------------------------------------------------------------------- Stress 11400.000 11440.746 1.004 Deflection 0.182 0.182 1.003 ---------------------------------------------------------------------------------- .. GENERATED FROM PYTHON SOURCE LINES 261-262 Stop MAPDL. .. GENERATED FROM PYTHON SOURCE LINES 262-263 .. code-block:: Python mapdl.exit() .. rst-class:: sphx-glr-timing **Total running time of the script:** (0 minutes 0.829 seconds) .. _sphx_glr_download_verif-manual_vm-002-beam_stresses_and_deflections.py: .. only:: html .. container:: sphx-glr-footer sphx-glr-footer-example .. container:: sphx-glr-download sphx-glr-download-jupyter :download:`Download Jupyter notebook: vm-002-beam_stresses_and_deflections.ipynb ` .. container:: sphx-glr-download sphx-glr-download-python :download:`Download Python source code: vm-002-beam_stresses_and_deflections.py ` .. only:: html .. rst-class:: sphx-glr-signature `Gallery generated by Sphinx-Gallery `_