.. DO NOT EDIT. .. THIS FILE WAS AUTOMATICALLY GENERATED BY SPHINX-GALLERY. .. TO MAKE CHANGES, EDIT THE SOURCE PYTHON FILE: .. "verif-manual/vm-018.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-018.py: .. _ref_vm18: Out-of-Plane Bending of a Curved Bar ------------------------------------ Problem description: - A portion of a horizontal circular ring, built-in at A, is loaded by a vertical (Z) load F applied at the end B. The ring has a solid circular cross-section of diameter d. Determine the deflection :math:`\delta` at end B, the maximum bending stress :math:`\sigma_{Bend}` , and the maximum torsional shear stress τ. Reference: - S. Timoshenko, Strength of Materials, Part I, Elementary Theory and Problems, 3rd Edition, D. Van Nostrand Co., Inc., New York, NY, 1955, pg. 412, eq. 241. Analysis type(s): - Static Analysis ``ANTYPE=0`` Element type(s): - Elastic Curved Pipe Element (PIPE18) - 3-D 3 Node Pipe Element (PIPE289) .. figure:: ../_static/vm18_setup1.png :align: center :width: 400 :alt: VM18 Curved Bar Problem Sketch .. figure:: ../_static/vm18_setup2.png :align: center :width: 400 :alt: VM18 Curved Bar Finite Element Model Material properties: - :math:`E = 30 \cdot 10^6 psi` - :math:`\mu = 0.3` Geometric properties: - :math:`r = 100 in` - :math:`d = 2 in` - :math:`\theta = 90°` Loading: - :math:`F = 50 lb` Analysis Assumptions and Modeling Notes: - Node 10 is arbitrarily located on the radius of curvature side of the element to define the plane of the elbow when PIPE18 elements are used. The wall thickness is set to half the diameter for a solid bar. Since the section has no hole in the middle, ovalization cannot occur and PIPE289 elements can be used to determine the deflection and stresses. .. GENERATED FROM PYTHON SOURCE LINES 74-102 .. code-block:: Python # sphinx_gallery_thumbnail_path = '_static/vm18_setup1.png' # Importing the `launch_mapdl` function from the `ansys.mapdl.core` module from ansys.mapdl.core import launch_mapdl import numpy as np import pandas as pd # Launch MAPDL with specified settings mapdl = launch_mapdl(loglevel="WARNING", print_com=True, remove_temp_dir_on_exit=True) # Clear the existing database mapdl.clear() # Run the FINISH command to exists normally from a processor mapdl.finish() # Set the ANSYS version mapdl.com("ANSYS MEDIA REL. 2022R2 (05/13/2022) REF. VERIF. MANUAL: REL. 2022R2") # Run the /VERIFY command mapdl.run("/VERIFY,VM18") # Set the title of the analysis mapdl.title("VM18 OUT-OF-PLANE BENDING OF A CURVED BAR") # Enter the model creation /Prep7 preprocessor 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 103-106 Define element type and real properties ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~ Use Elastic Curved Pipe element (PIPE18) and set KEYOPT(6)=2 for printing member forces. .. GENERATED FROM PYTHON SOURCE LINES 106-111 .. code-block:: Python mapdl.et(1, "PIPE18", "", "", "", "", "", 2) # Define geometry parameters (OD, wall thickness, radius) using "r" command (real constant) mapdl.r(1, 2, 1, 100) .. rst-class:: sphx-glr-script-out .. code-block:: none REAL CONSTANT SET 1 ITEMS 1 TO 6 2.0000 1.0000 100.00 0.0000 0.0000 0.0000 .. GENERATED FROM PYTHON SOURCE LINES 112-116 Define material ~~~~~~~~~~~~~~~ Set up the material and its type (a single material), Young's modulus of 30e6 and Poisson's ratio NUXY of 0.3 is specified. .. GENERATED FROM PYTHON SOURCE LINES 116-119 .. code-block:: Python mapdl.mp("EX", 1, 30e6) mapdl.mp("NUXY", 1, 0.3) .. rst-class:: sphx-glr-script-out .. code-block:: none MATERIAL 1 NUXY = 0.3000000 .. GENERATED FROM PYTHON SOURCE LINES 120-124 Define geometry ~~~~~~~~~~~~~~~ Set up the nodes and elements. This creates a mesh just like in the problem setup. .. GENERATED FROM PYTHON SOURCE LINES 124-133 .. code-block:: Python # Define nodes mapdl.n(1, 100) mapdl.n(2, "", 100) mapdl.n(10) # Define element mapdl.e(1, 2, 10) .. rst-class:: sphx-glr-script-out .. code-block:: none 1 .. GENERATED FROM PYTHON SOURCE LINES 134-138 Define boundary conditions and load ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~ Fix all dofs at node 1. Specify nodal force F = -50 lb along Z direction at node 2. Then exit prep7 processor. .. GENERATED FROM PYTHON SOURCE LINES 138-150 .. code-block:: Python mapdl.d(1, "ALL") # Define boundary conditions mapdl.f(2, "FZ", -50) # Define load # Selects all entities mapdl.allsel() # Element plot mapdl.eplot(vtk=False) # Finish preprocessing processor mapdl.finish() .. image-sg:: /verif-manual/images/sphx_glr_vm-018_001.png :alt: vm 018 :srcset: /verif-manual/images/sphx_glr_vm-018_001.png :class: sphx-glr-single-img .. rst-class:: sphx-glr-script-out .. code-block:: none ***** ROUTINE COMPLETED ***** CP = 0.000 .. GENERATED FROM PYTHON SOURCE LINES 151-154 Solve ~~~~~ Enter solution mode and solve the system. .. GENERATED FROM PYTHON SOURCE LINES 154-167 .. code-block:: Python mapdl.slashsolu() # Set the analysis type to STATIC mapdl.antype("STATIC") # Set output options mapdl.outpr("BASIC", 1) # Perform the solution mapdl.solve() # exists solution processor mapdl.finish() .. rst-class:: sphx-glr-script-out .. code-block:: none FINISH SOLUTION PROCESSING *** NOTE *** CP = 0.000 TIME= 00:00:00 Distributed parallel processing has been reactivated. ***** ROUTINE COMPLETED ***** CP = 0.000 .. GENERATED FROM PYTHON SOURCE LINES 168-171 Post-processing ~~~~~~~~~~~~~~~ Enter post-processing. Compute deflection and stress quantities. .. GENERATED FROM PYTHON SOURCE LINES 171-191 .. code-block:: Python mapdl.post1() # Set the current results set to the last set to be read from result file mapdl.set("LAST") # Get displacement results at node 2 in the Z direction def_z = mapdl.get("DEF", "NODE", 2, "U", "Z") # Create an element table for bending stresses using ETABLE command strs_ben = mapdl.etable("STRS_BEN", "NMISC", 91) # Create an element table for shear stresses using ETABLE command strs_shr = mapdl.etable("STRS_SHR", "LS", 4) # Get bending stresses (ETAB: STRS_BEN) for element 1 strss_b = mapdl.get("STRSS_B", "ELEM", 1, "ETAB", "STRS_BEN") # Get shear stresses (ETAB: STRS_SHR) for element 1 strss_t = mapdl.get("STRSS_T", "ELEM", 1, "ETAB", "STRS_SHR") .. GENERATED FROM PYTHON SOURCE LINES 192-194 Verify the results. ~~~~~~~~~~~~~~~~~~~ .. GENERATED FROM PYTHON SOURCE LINES 194-217 .. code-block:: Python # Set target values target_val = [-2.648, 6366, -3183] # Fill result values sim_res = [def_z, strss_b, strss_t] col_headers = ["TARGET", "Mechanical APDL", "RATIO"] row_headers = ["Deflection (in)", "Stress_Bend (psi)", "Shear Stress (psi)"] data = [target_val, sim_res, np.abs(target_val) / np.abs(sim_res)] title = f""" ------------------- VM18 RESULTS COMPARISON --------------------- PIPE18: ------- """ print(title) print(pd.DataFrame(np.transpose(data), row_headers, col_headers)) .. rst-class:: sphx-glr-script-out .. code-block:: none ------------------- VM18 RESULTS COMPARISON --------------------- PIPE18: ------- TARGET Mechanical APDL RATIO Deflection (in) -2.648 -2.649729 0.999348 Stress_Bend (psi) 6366.000 6366.197750 0.999969 Shear Stress (psi) -3183.000 -3183.098880 0.999969 .. GENERATED FROM PYTHON SOURCE LINES 218-220 Finish the post-processing processor. ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~ .. GENERATED FROM PYTHON SOURCE LINES 220-222 .. code-block:: Python mapdl.finish() .. rst-class:: sphx-glr-script-out .. code-block:: none EXIT THE MAPDL POST1 DATABASE PROCESSOR ***** ROUTINE COMPLETED ***** CP = 0.000 .. GENERATED FROM PYTHON SOURCE LINES 223-225 Clears the database without restarting. ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~ .. GENERATED FROM PYTHON SOURCE LINES 225-227 .. code-block:: Python mapdl.run("/CLEAR,NOSTART") .. rst-class:: sphx-glr-script-out .. code-block:: none CLEAR MAPDL DATABASE AND RESTART Ansys Mechanical Enterprise Academic Student .. GENERATED FROM PYTHON SOURCE LINES 228-230 Set a new title for the analysis ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~ .. GENERATED FROM PYTHON SOURCE LINES 230-232 .. code-block:: Python mapdl.title("VM18 OUT-OF-PLANE BENDING OF A CURVED BAR Using PIPE289 ELEMENT MODEL") .. rst-class:: sphx-glr-script-out .. code-block:: none TITLE= VM18 OUT-OF-PLANE BENDING OF A CURVED BAR Using PIPE289 ELEMENT MODEL .. GENERATED FROM PYTHON SOURCE LINES 233-235 Switches to the preprocessor (PREP7) ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~ .. GENERATED FROM PYTHON SOURCE LINES 235-237 .. code-block:: Python 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 238-241 Define element type and section properties ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~ Use 3-D 3-Node Pipe element (PIPE289) and set KEYOPT(4)= 2 Thick pipe theory. .. GENERATED FROM PYTHON SOURCE LINES 241-245 .. code-block:: Python mapdl.et(1, "PIPE289", "", "", "", 2) mapdl.sectype(1, "PIPE") # Set section type PIPE mapdl.secdata(2, 1, 16) # Set section data (OD, wall thickness) .. rst-class:: sphx-glr-script-out .. code-block:: none SECTION ID NUMBER: 1 PIPE SECTION NAME IS: PIPE SECTION DATA SUMMARY: Outside Diameter = 2.0000 Thickness = 1.0000 Area = 3.1414 Iyy = 0.78529 Torsion Constant = 1.5706 Shear Correction-yy = 0.86083 Num of Circum Cells = 16 .. GENERATED FROM PYTHON SOURCE LINES 246-250 Define material ~~~~~~~~~~~~~~~ Set up the material and its type (a single material), Young's modulus of 30e6 and Poisson's ratio NUXY of 0.3 is specified. .. GENERATED FROM PYTHON SOURCE LINES 250-253 .. code-block:: Python mapdl.mp("EX", 1, 30e6) mapdl.mp("NUXY", 1, 0.3) .. rst-class:: sphx-glr-script-out .. code-block:: none MATERIAL 1 NUXY = 0.3000000 .. GENERATED FROM PYTHON SOURCE LINES 254-258 Define geometry ~~~~~~~~~~~~~~~ Set up the nodes and elements. This creates a mesh just like in the problem setup. .. GENERATED FROM PYTHON SOURCE LINES 258-275 .. code-block:: Python mapdl.csys(1) # Set coordinate system to 1 mapdl.n(1, 100) # Define nodes # Generate additional nodes mapdl.ngen(19, 1, 1, "", "", "", 5) # Define element mapdl.e(1, 3, 2) # Generate additional elements from an existing pattern mapdl.egen(9, 2, -1) # Reset coordinate system to global mapdl.csys(0) .. rst-class:: sphx-glr-script-out .. code-block:: none ACTIVE COORDINATE SYSTEM SET TO 0 (CARTESIAN) .. GENERATED FROM PYTHON SOURCE LINES 276-280 Define boundary conditions and load ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~ Fix all dofs at node 1. Specify nodal force F = -50 lb along Z direction at node 19. Then exit prep7 processor. .. GENERATED FROM PYTHON SOURCE LINES 280-291 .. code-block:: Python mapdl.d(1, "ALL") mapdl.f(19, "FZ", -50) # Selects all entities mapdl.allsel() # Element plot mapdl.eplot(vtk=False) # exists pre-processing processor mapdl.finish() .. image-sg:: /verif-manual/images/sphx_glr_vm-018_002.png :alt: vm 018 :srcset: /verif-manual/images/sphx_glr_vm-018_002.png :class: sphx-glr-single-img .. rst-class:: sphx-glr-script-out .. code-block:: none ***** ROUTINE COMPLETED ***** CP = 0.000 .. GENERATED FROM PYTHON SOURCE LINES 292-295 Solve ~~~~~ Enter solution mode and solve the system. .. GENERATED FROM PYTHON SOURCE LINES 295-308 .. code-block:: Python mapdl.slashsolu() # Set the analysis type to STATIC mapdl.antype("STATIC") # Set output options mapdl.outpr("BASIC", 1) # Perform the solution mapdl.solve() # exists solution processor mapdl.finish() .. rst-class:: sphx-glr-script-out .. code-block:: none FINISH SOLUTION PROCESSING ***** ROUTINE COMPLETED ***** CP = 0.000 .. GENERATED FROM PYTHON SOURCE LINES 309-312 Post-processing ~~~~~~~~~~~~~~~ Enter post-processing. Compute deflection and stress quantities. .. GENERATED FROM PYTHON SOURCE LINES 312-337 .. code-block:: Python mapdl.post1() # Set the current results set to the last set mapdl.set("LAST") mapdl.graphics("POWER") # Set graphics mode to POWER mapdl.eshape(1) # Set element shape mapdl.view(1, 1, 1, 1) # Set view # Get displacement results at node 19 in the Z direction def_z = mapdl.get("DEF", "NODE", 19, "U", "Z") # Create an element table for bending stresses using ETABLE command strs_ben = mapdl.etable("STRS_BEN", "SMISC", 35) # Get bending stresses (ETAB: STRS_BEN) for element 1 using ETABLE command strss_b = mapdl.get("STRSS_B", "ELEM", 1, "ETAB", "STRS_BEN") # for graphics displays mapdl.show(option="REV") # Plot elemtal solution values for SXY component mapdl.plesol("S", "XY") # Get minimum shear stress shear_sxy = mapdl.get("SHEAR", "PLNSOL", 0, "MIN") mapdl.show("close") .. image-sg:: /verif-manual/images/sphx_glr_vm-018_003.png :alt: vm 018 :srcset: /verif-manual/images/sphx_glr_vm-018_003.png :class: sphx-glr-single-img .. GENERATED FROM PYTHON SOURCE LINES 338-340 Verify the results. ~~~~~~~~~~~~~~~~~~~ .. GENERATED FROM PYTHON SOURCE LINES 340-361 .. code-block:: Python # Set target values target_val = [-2.648, 6366, -3183] # Fill result values sim_res = [def_z, strss_b, shear_sxy] col_headers = ["TARGET", "Mechanical APDL", "RATIO"] row_headers = ["Deflection (in)", "Stress_Bend (psi)", "Shear Stress (psi)"] data = [target_val, sim_res, np.abs(target_val) / np.abs(sim_res)] title = f""" PIPE289: -------- """ print(title) print(pd.DataFrame(np.transpose(data), row_headers, col_headers)) .. rst-class:: sphx-glr-script-out .. code-block:: none PIPE289: -------- TARGET Mechanical APDL RATIO Deflection (in) -2.648 -2.649533 0.999422 Stress_Bend (psi) 6366.000 -6446.631350 0.987492 Shear Stress (psi) -3183.000 -3199.472410 0.994852 .. GENERATED FROM PYTHON SOURCE LINES 362-364 Finish the post-processing processor. ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~ .. GENERATED FROM PYTHON SOURCE LINES 364-366 .. code-block:: Python mapdl.finish() .. rst-class:: sphx-glr-script-out .. code-block:: none EXIT THE MAPDL POST1 DATABASE PROCESSOR ***** ROUTINE COMPLETED ***** CP = 0.000 .. GENERATED FROM PYTHON SOURCE LINES 367-369 Stop MAPDL. ~~~~~~~~~~~ .. GENERATED FROM PYTHON SOURCE LINES 369-370 .. code-block:: Python mapdl.exit() .. rst-class:: sphx-glr-timing **Total running time of the script:** (0 minutes 1.543 seconds) .. _sphx_glr_download_verif-manual_vm-018.py: .. only:: html .. container:: sphx-glr-footer sphx-glr-footer-example .. container:: sphx-glr-download sphx-glr-download-jupyter :download:`Download Jupyter notebook: vm-018.ipynb ` .. container:: sphx-glr-download sphx-glr-download-python :download:`Download Python source code: vm-018.py ` .. only:: html .. rst-class:: sphx-glr-signature `Gallery generated by Sphinx-Gallery `_