.. DO NOT EDIT. .. THIS FILE WAS AUTOMATICALLY GENERATED BY SPHINX-GALLERY. .. TO MAKE CHANGES, EDIT THE SOURCE PYTHON FILE: .. "verif-manual/vm-021.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-021.py: .. _ref_vm21: Tie Rod with Lateral Loading ---------------------------- Problem description: - A tie rod is subjected to the action of a tensile force F and a uniform lateral load p. Determine the maximum deflection :math:`z_{max}`, the slope :math:`\theta` at the left-hand end, and the maximum bending moment :math:`M_{max}`. In addition, determine the same three quantities for the unstiffened tie rod (F = 0). Reference: - S. Timoshenko, Strength of Materials, Part II, Elementary Theory and Problems, 3rd Edition, D. Van Nostrand Co., Inc., New York, NY, 1956, pg. 42, article 6. Analysis type(s): - Static, Stress Stiffening Analysis ``ANTYPE=0`` Element type(s): - 3-D 2 node beam (BEAM188) .. image:: ../_static/vm21_setup.png :width: 400 :alt: VM21 Tie Rod Problem Sketch Material properties: - :math:`E = 30 \cdot 10^6 psi` Geometric properties: - :math:`l = 200 in` - :math:`b = h = 2.5 in` Loading: - :math:`F = 21,972.6 lb` - :math:`p = 1.79253 lb/in` Analysis Assumptions and Modeling Notes: - Due to symmetry, only one-half of the beam is modeled. The full load is applied for each iteration. The first solution represents the unstiffened case. The second solution represents the stiffened case. .. GENERATED FROM PYTHON SOURCE LINES 65-93 .. code-block:: Python # sphinx_gallery_thumbnail_path = '_static/vm21_setup.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 options mapdl = launch_mapdl(loglevel="WARNING", print_com=True, remove_temp_dir_on_exit=True) # Clear the current 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 for VM21 mapdl.run("/VERIFY,VM21") # Set the title of the analysis mapdl.title("VM21 TIE ROD WITH LATERAL LOADING NO STREES STIFFENING") # 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 94-97 Define element type and section properties ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~ Use 3-D 2-Node Beam Element and specify cubic shape function via setting Keyopt(3)=3. .. GENERATED FROM PYTHON SOURCE LINES 97-104 .. code-block:: Python mapdl.et(1, "BEAM188") mapdl.keyopt(1, 3, 3) # Set KEYOPT(3) to 3 cubic shape function mapdl.sectype(1, "BEAM", "RECT") # Specify section properties for the beam element mapdl.secdata(2.5, 2.5) # Define section data .. rst-class:: sphx-glr-script-out .. code-block:: none SECTION ID NUMBER IS: 1 BEAM SECTION TYPE IS: Rectangle BEAM SECTION NAME IS: COMPUTED BEAM SECTION DATA SUMMARY: Area = 6.2500 Iyy = 3.2552 Iyz = 0.41633E-16 Izz = 3.2552 Warping Constant = 0.23513E-01 Torsion Constant = 5.5714 Centroid Y = 0.88818E-17 Centroid Z = 0.17764E-16 Shear Center Y = 0.27711E-16 Shear Center Z =-0.14282E-15 Shear Correction-xy = 0.84211 Shear Correction-yz =-0.95989E-15 Shear Correction-xz = 0.84211 Beam Section is offset to CENTROID of cross section .. GENERATED FROM PYTHON SOURCE LINES 105-109 Define material ~~~~~~~~~~~~~~~ Set up the material and its type (a single material), Young's modulus of 30e6 and Poisson's ratio PRXY of 0.3 is specified. .. GENERATED FROM PYTHON SOURCE LINES 109-113 .. code-block:: Python mapdl.mp("EX", 1, 30e6) mapdl.mp("PRXY", "", 0.3) .. rst-class:: sphx-glr-script-out .. code-block:: none MATERIAL 1 PRXY = 0.3000000 .. GENERATED FROM PYTHON SOURCE LINES 114-118 Define geometry ~~~~~~~~~~~~~~~ Set up the nodes and elements. This creates a mesh just like in the problem setup. .. GENERATED FROM PYTHON SOURCE LINES 118-131 .. code-block:: Python mapdl.n(1) # define nodes mapdl.n(5, 100) # Generate additional nodes mapdl.fill() # Define elements mapdl.e(1, 2) # Generate additional elements from an existing pattern mapdl.egen(4, 1, 1) .. rst-class:: sphx-glr-script-out .. code-block:: none GENERATE 4 TOTAL SETS OF ELEMENTS WITH NODE INCREMENT OF 1 SET IS SELECTED ELEMENTS IN RANGE 1 TO 1 IN STEPS OF 1 MAXIMUM ELEMENT NUMBER= 4 .. GENERATED FROM PYTHON SOURCE LINES 132-138 Define boundary conditions and loading ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~ Apply a displacement boundary condition in the UY, ROTX and ROTZ directions to all nodes. Specify symmetry degree-of-freedom constraints on nodes, surface normal to X-dir (default). Apply a tensile force in X-dir, F = -21972.6 lb and a uniform lateral load, p = 1.79253 lb/in. Then exit prep7 processor. .. GENERATED FROM PYTHON SOURCE LINES 138-162 .. code-block:: Python mapdl.d("ALL", "UY", "", "", "", "", "ROTX", "ROTZ") mapdl.d(1, "UZ") # Select nodes for symmetry boundary mapdl.nsel("S", "", "", 5) mapdl.dsym("SYMM", "X") # Select all nodes mapdl.nsel("ALL") # Apply nodal force along x-direction mapdl.f(1, "FX", -21972.6) # Specifies surface loads on beam and pipe elements. mapdl.sfbeam("ALL", 1, "PRES", 1.79253) # Selects all entities mapdl.allsel() # Element plot mapdl.eplot() # Finish pre-processing processor mapdl.finish() .. tab-set:: .. tab-item:: Static Scene .. image-sg:: /verif-manual/images/sphx_glr_vm-021_001.png :alt: vm 021 :srcset: /verif-manual/images/sphx_glr_vm-021_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-021_001.vtksz .. rst-class:: sphx-glr-script-out .. code-block:: none ***** ROUTINE COMPLETED ***** CP = 0.000 .. GENERATED FROM PYTHON SOURCE LINES 163-166 Solve ~~~~~ Enter solution mode and solve the system. .. GENERATED FROM PYTHON SOURCE LINES 166-175 .. code-block:: Python mapdl.slashsolu() # Set the analysis type to STATIC mapdl.antype("STATIC") # 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 176-179 Post-processing ~~~~~~~~~~~~~~~ Enter post-processing. Compute displacement and rotation quantities. .. GENERATED FROM PYTHON SOURCE LINES 179-196 .. code-block:: Python mapdl.post1() # Select nodes for results output mapdl.nsel("S", "", "", 1, 5, 4) # Print displacement quantities in Z direction mapdl.prnsol("U", "Z") # Print rotation quantities in Y direction mapdl.prnsol("ROT", "Y") # Select all nodes mapdl.nsel("ALL") # Print results solution mapdl.prrsol() .. rst-class:: sphx-glr-script-out .. code-block:: none 'PRINT REACTION SOLUTIONS PER NODE\n *****MAPDL VERIFICATION RUN ONLY*****\n DO NOT USE RESULTS FOR PRODUCTION\n \n ***** POST1 TOTAL REACTION SOLUTION LISTING ***** \n \n LOAD STEP= 1 SUBSTEP= 1 \n TIME= 1.0000 LOAD CASE= 0 \n \n THE FOLLOWING X,Y,Z SOLUTIONS ARE IN THE GLOBAL COORDINATE SYSTEM \n \n NODE FX FY FZ MX MY MZ \n 1 0.0000 179.25 0.32697E-014 -0.69256E-014\n 2 0.0000 -0.21095E-015 -0.10747E-013\n 3 0.0000 -0.21095E-015 -0.71644E-014\n 4 0.0000 -0.21095E-015 -0.35822E-014\n 5 21973. 0.0000 -0.10548E-015 -8962.7 0.11439E-012\n\n TOTAL VALUES\n VALUE 21973. 0.0000 179.25 0.25314E-014 -8962.7 0.85973E-013' .. GENERATED FROM PYTHON SOURCE LINES 197-199 Inline functions in PyMAPDL to query node ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~ .. GENERATED FROM PYTHON SOURCE LINES 199-203 .. code-block:: Python q = mapdl.queries RGHT_END = q.node(200, 0, 0) # store node number to 'RGHT_END' with coordinate (4,0,0) LFT_END = q.node(0, 0, 0) # store node number to 'LFT_END' with coordinate (4,0,0) .. GENERATED FROM PYTHON SOURCE LINES 204-206 Retrieve nodal deflection and rotation ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~ .. GENERATED FROM PYTHON SOURCE LINES 206-216 .. code-block:: Python # Get results at node RGHT_END uz_mx_c2 = mapdl.get("UZ_MX_C2", "NODE", RGHT_END, "U", "Z") # Get results at node LFT_END slope_c2 = mapdl.get("SLOPE_C2", "NODE", LFT_END, "ROT", "Y") # exists post-processing processor 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 217-220 Post-processing ~~~~~~~~~~~~~~~ Enter /post26 time-history post-processing processor. .. GENERATED FROM PYTHON SOURCE LINES 220-233 .. code-block:: Python mapdl.post26() # Specifies the total reaction force data to be stored at nodes associated to RGHT_END mapdl.rforce(2, RGHT_END, "M", "Y") # Store results mapdl.store() # Get maximum moment at node RGHT_END m_mx_c2 = mapdl.get("M_MX_C2", "VARI", 2, "EXTREM", "VMAX") # exists post-processing processor mapdl.finish() .. rst-class:: sphx-glr-script-out .. code-block:: none ***** ROUTINE COMPLETED ***** CP = 0.000 .. GENERATED FROM PYTHON SOURCE LINES 234-236 Set a new title for the analysis ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~ .. GENERATED FROM PYTHON SOURCE LINES 236-238 .. code-block:: Python mapdl.title("VM21 TIE ROD WITH LATERAL LOADING STRESS STIFFENING PRESENT") .. rst-class:: sphx-glr-script-out .. code-block:: none TITLE= VM21 TIE ROD WITH LATERAL LOADING STRESS STIFFENING PRESENT .. GENERATED FROM PYTHON SOURCE LINES 239-242 Solve ~~~~~ Enter new solution mode and solve the nonlinear system including stress stiffening. .. GENERATED FROM PYTHON SOURCE LINES 242-257 .. code-block:: Python mapdl.slashsolu() # Set number of substeps to 5 mapdl.nsubst(5) # Activate auto time stepping mapdl.autots("ON") # Activate nonlinear geometry mapdl.nlgeom("ON") # Set a smaller convergence tolerance mapdl.cnvtol("F", "", 0.0001, "", 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 258-261 Post-processing ~~~~~~~~~~~~~~~ Enter post-processing. Compute displacement and rotation quantities. .. GENERATED FROM PYTHON SOURCE LINES 261-281 .. code-block:: Python mapdl.post1() # Select nodes for results output mapdl.nsel("S", "", "", 1, 5, 4) # Print displacement quantities in Z direction mapdl.prnsol("U", "Z") # Print rotation quantities in Y direction mapdl.prnsol("ROT", "Y") # Print results solution mapdl.prrsol() # Get results at node RGHT_END uz_mx_c1 = mapdl.get("UZ_MX_C1", "NODE", RGHT_END, "U", "Z") # Get results at node LFT_END slope_c1 = mapdl.get("SLOPE_C1", "NODE", LFT_END, "ROT", "Y") # exists post-processing processor 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 282-285 Post-processing ~~~~~~~~~~~~~~~ Enter /post26 time-history post-processing processor. .. GENERATED FROM PYTHON SOURCE LINES 285-296 .. code-block:: Python mapdl.post26() # Specifies the total reaction force data to be stored at nodes associated to RGHT_END mapdl.rforce(2, RGHT_END, "M", "Y") # Store results mapdl.store() # Get maximum moment at node RGHT_END m_mx_c1 = mapdl.get("M_MX_C1", "VARI", 2, "EXTREM", "VMAX") .. GENERATED FROM PYTHON SOURCE LINES 297-299 Verify the results. ~~~~~~~~~~~~~~~~~~~ .. GENERATED FROM PYTHON SOURCE LINES 299-344 .. code-block:: Python # Set target values target_res = [-0.19945, 0.0032352, -4580.1] target_res_strss = [-0.38241, 0.0061185, -8962.7] # Fill result values sim_res = [uz_mx_c1, slope_c1, m_mx_c1] sim_res_strss = [uz_mx_c2, slope_c2, m_mx_c2] title = f""" ------------------- VM21 RESULTS COMPARISON --------------------- F neq 0 (stiffened): -------------------- """ col_headers = ["TARGET", "Mechanical APDL", "RATIO"] row_headers = ["Z_max, in", "Slope, rad", "M_max , in-lb"] data = [target_res, sim_res, np.abs(target_res) / np.abs(sim_res)] print(title) print(pd.DataFrame(np.transpose(data), row_headers, col_headers)) title = f""" F = 0 (unstiffened): -------------------- """ row_headers = ["Z_max, in", "Slope, rad", "M_max , in-lb"] data = [ target_res_strss, sim_res_strss, np.abs(target_res_strss) / np.abs(sim_res_strss), ] print(title) print(pd.DataFrame(np.transpose(data), row_headers, col_headers)) .. rst-class:: sphx-glr-script-out .. code-block:: none ------------------- VM21 RESULTS COMPARISON --------------------- F neq 0 (stiffened): -------------------- TARGET Mechanical APDL RATIO Z_max, in -0.199450 -0.199560 0.999449 Slope, rad 0.003235 0.003235 0.999957 M_max , in-lb -4580.100000 -4579.818510 1.000061 F = 0 (unstiffened): -------------------- TARGET Mechanical APDL RATIO Z_max, in -0.382410 -0.382554 0.999623 Slope, rad 0.006118 0.006119 1.000000 M_max , in-lb -8962.700000 -8962.650000 1.000006 .. GENERATED FROM PYTHON SOURCE LINES 345-347 Finish the post-processing processor. ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~ .. GENERATED FROM PYTHON SOURCE LINES 347-349 .. code-block:: Python mapdl.finish() .. rst-class:: sphx-glr-script-out .. code-block:: none ***** ROUTINE COMPLETED ***** CP = 0.000 .. GENERATED FROM PYTHON SOURCE LINES 350-352 Stop MAPDL. ~~~~~~~~~~~ .. GENERATED FROM PYTHON SOURCE LINES 352-353 .. code-block:: Python mapdl.exit() .. rst-class:: sphx-glr-timing **Total running time of the script:** (0 minutes 0.696 seconds) .. _sphx_glr_download_verif-manual_vm-021.py: .. only:: html .. container:: sphx-glr-footer sphx-glr-footer-example .. container:: sphx-glr-download sphx-glr-download-jupyter :download:`Download Jupyter notebook: vm-021.ipynb ` .. container:: sphx-glr-download sphx-glr-download-python :download:`Download Python source code: vm-021.py ` .. only:: html .. rst-class:: sphx-glr-signature `Gallery generated by Sphinx-Gallery `_