.. DO NOT EDIT. .. THIS FILE WAS AUTOMATICALLY GENERATED BY SPHINX-GALLERY. .. TO MAKE CHANGES, EDIT THE SOURCE PYTHON FILE: .. "verif-manual/vm-025.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-025.py: .. _ref_vm25: Stresses in a Long Cylinder --------------------------- Problem description: - A long thick-walled cylinder is initially subjected to an internal pressure p. Determine the radial displacement :math:`\delta_r` at the inner surface, the radial stress :math:`\sigma_r` , and tangential stress :math:`\sigma_t` , at the inner and outer surfaces and at the middle wall thickness. Internal pressure is then removed and the cylinder is subjected to a rotation ω about its center line. Determine the radial :math:`\sigma_r` and tangential :math:`\sigma_t` stresses at the inner wall and at an interior point located at r = Xi. Reference: - S. Timoshenko, Strength of Materials, Part II, Elementary Theory and Problems, 3rd Edition, D. Van Nostrand Co., Inc., New York, NY, 1956, pg. 213, problem 1 and pg. 213, article 42. Analysis type(s): - Static Analysis ``ANTYPE=0`` Element type(s): - 2-D 8-Node Structural Solid Elements (PLANE183) .. image:: ../_static/vm25_setup.png :width: 400 :alt: VM25 Long Cylinder Problem Sketch Material properties: - :math:`E = 30 \cdot 10^6 psi` - :math:`\mu = 0.3` - :math:`\rho = 0.00073 lb-sec^2/in^4` Geometric properties: - :math:`a = 4 inches` - :math:`b = 8 inches` - :math:`X_i = 5.43 inches` Loading: - :math:`p = 30,000 psi` - :math:`\Omega = 1000 rad/sec` Analysis Assumptions and Modeling Notes: - The axial length is arbitrarily selected. Elements are oriented such that surface stresses may be obtained at the inner and outer cylinder surfaces. POST1 is used to display linearized stresses through the thickness of the cylinder when it is subjected to an internal pressure. .. GENERATED FROM PYTHON SOURCE LINES 72-103 .. code-block:: Python # sphinx_gallery_thumbnail_path = '_static/vm25_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 # Launching MAPDL with specified settings mapdl = launch_mapdl(loglevel="WARNING", print_com=True, remove_temp_dir_on_exit=True) # Clearing the MAPDL 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 VM25 mapdl.run("/VERIFY,VM25") # Set the title of the analysis mapdl.title("VM25 Stresses in a Long Cylinder") # Enter the model creation /Prep7 preprocessor mapdl.prep7() # Deactivate automatic (smart) element sizing mapdl.smrtsize(sizlvl="OFF") .. rst-class:: sphx-glr-script-out .. code-block:: none DEACTIVATE SMART SIZING .. GENERATED FROM PYTHON SOURCE LINES 104-108 Define element type and properties ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~ Use 2-D 8-Node or 6-Node Structural Solid and specify Axisymmetric element behavior via setting Keyopt(3)=1. .. GENERATED FROM PYTHON SOURCE LINES 108-110 .. code-block:: Python mapdl.et(1, "PLANE183", "", "", 1) .. rst-class:: sphx-glr-script-out .. code-block:: none 1 .. GENERATED FROM PYTHON SOURCE LINES 111-115 Define material ~~~~~~~~~~~~~~~ Set up the material and its type (a single material), Young's modulus of 30e6, Density, rho = 0.00073 and Poisson's ratio NUXY of 0.3 is specified. .. GENERATED FROM PYTHON SOURCE LINES 115-119 .. code-block:: Python mapdl.mp("EX", 1, 30e6) mapdl.mp("DENS", 1, 0.00073) 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-154 .. code-block:: Python # Defining node 1 with coordinates (4) mapdl.n(1, 4) # Defining node 2 with coordinates (4+4/14) mapdl.n(2, "4+4/14") # Defining node 3 with coordinates (4+4/14, 0.5) mapdl.n(3, "4+4/14", 0.5) # Defining node 4 with coordinates (4, 0.5) mapdl.n(4, 4, 0.5) # Defining node 5 with coordinates (4+4/28) mapdl.n(5, "4+4/28") # Defining node 6 with coordinates (4+4/14, 0.25) mapdl.n(6, "4+4/14", 0.25) # Defining node 7 with coordinates (4+4/28, 0.5) mapdl.n(7, "4+4/28", 0.5) # Defining node 8 with coordinates (4, 0.25) mapdl.n(8, 4, 0.25) # Creating element 1 with nodes 2, 3, 4, 5, 6, 7, and 8 mapdl.e(1, 2, 3, 4, 5, 6, 7, 8) # Generating additional nodes along element 14 using the specified # parameters (from an existing pattern) mapdl.egen(14, 8, 1, "", "", "", "", "", "", "", "4/14") # Merging nodes based on their coordinates mapdl.nummrg("NODE") # Generate a mesh using EGEN command mapdl.egen(2, 111, 1, 14, "", "", "", "", "", "", "", 0.5) # Merge the nodes that share the same coordinates mapdl.nummrg("NODE") .. rst-class:: sphx-glr-script-out .. code-block:: none MERGE COINCIDENT NODES WITHIN TOLERANCE OF 0.10000E-03 NODE 4 USED FOR NODES 112 NODE 7 USED FOR NODES 116 NODE 3 USED FOR NODES 113 NODE 15 USED FOR NODES 124 NODE 11 USED FOR NODES 121 NODE 23 USED FOR NODES 132 NODE 19 USED FOR NODES 129 NODE 31 USED FOR NODES 140 NODE 27 USED FOR NODES 137 NODE 39 USED FOR NODES 148 NODE 35 USED FOR NODES 145 NODE 47 USED FOR NODES 156 NODE 43 USED FOR NODES 153 NODE 55 USED FOR NODES 164 NODE 51 USED FOR NODES 161 NODE 63 USED FOR NODES 172 NODE 59 USED FOR NODES 169 NODE 71 USED FOR NODES 180 NODE 67 USED FOR NODES 177 NODE 79 USED FOR NODES 188 NODE 75 USED FOR NODES 185 NODE 87 USED FOR NODES 196 NODE 83 USED FOR NODES 193 NODE 95 USED FOR NODES 204 NODE 91 USED FOR NODES 201 NODE 103 USED FOR NODES 212 NODE 99 USED FOR NODES 209 NODE 111 USED FOR NODES 220 NODE 107 USED FOR NODES 217 A total of 29 nodes were merged at 29 locations. .. GENERATED FROM PYTHON SOURCE LINES 155-161 Define coupling and boundary conditions ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~ Apply a displacement boundary condition in the vertical direction (UY) to all nodes. Couple the axial displacements at the unconstrained Y-dir. Apply internal pressure of 30000 psi is applied on nodes. Also, apply dummy pressure for surface printout. Then exit prep7 processor. .. GENERATED FROM PYTHON SOURCE LINES 161-192 .. code-block:: Python mapdl.nsel("S", "LOC", "Y", 0) # Select nodes located on the Y-axis # Apply a displacement boundary condition mapdl.d("ALL", "UY") mapdl.nsel("S", "LOC", "Y", 1) # Select nodes located on the positive Y-axis # Couple the axial displacements at the unconstrained Y-dir mapdl.cp(1, "UY", "ALL") mapdl.nsel( "S", "LOC", "X", 4 ) # Select nodes located on the X-axis at a specific coordinate # Apply internal pressure on nodes mapdl.sf("", "PRES", 30000) mapdl.nsel( "S", "LOC", "X", 8 ) # Select nodes located on the X-axis at a different coordinate # Apply dummy pressure for surface printout mapdl.sf("", "PRES", 1e-10) # Selects all entities mapdl.allsel() # Element plot mapdl.eplot(background="w") # Finish the pre-processing processor mapdl.finish() # Save the finite element model mapdl.save("MODEL") .. tab-set:: .. tab-item:: Static Scene .. image-sg:: /verif-manual/images/sphx_glr_vm-025_001.png :alt: vm 025 :srcset: /verif-manual/images/sphx_glr_vm-025_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-025_001.vtksz .. rst-class:: sphx-glr-script-out .. code-block:: none ALL CURRENT MAPDL DATA WRITTEN TO FILE NAME= FOR POSSIBLE RESUME FROM THIS POINT .. GENERATED FROM PYTHON SOURCE LINES 193-196 Solve ~~~~~ Enter solution mode and solve the system. .. GENERATED FROM PYTHON SOURCE LINES 196-207 .. code-block:: Python mapdl.slashsolu() # Set the analysis type to STATIC mapdl.antype("STATIC") # Output results for all nodes mapdl.outpr("", "ALL") # Perform the solution for the load step 1, which is the internal pressure mapdl.solve() # exists solution processor for load case 1 mapdl.finish() .. rst-class:: sphx-glr-script-out .. code-block:: none FINISH SOLUTION PROCESSING ***** ROUTINE COMPLETED ***** CP = 0.000 .. GENERATED FROM PYTHON SOURCE LINES 208-211 Post-processing ~~~~~~~~~~~~~~~ Enter post-processing. Compute displacement and stress components. .. GENERATED FROM PYTHON SOURCE LINES 211-216 .. code-block:: Python mapdl.post1() # Set the load step 1 and substep to 1 mapdl.set(1, 1) .. rst-class:: sphx-glr-script-out .. code-block:: none USE LOAD STEP 1 SUBSTEP 1 FOR LOAD CASE 0 SET COMMAND GOT LOAD STEP= 1 SUBSTEP= 1 CUMULATIVE ITERATION= 1 TIME/FREQUENCY= 1.0000 TITLE= VM25 Stresses in a Long Cylinder .. GENERATED FROM PYTHON SOURCE LINES 217-219 Inline functions in PyMAPDL to query node ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~ .. GENERATED FROM PYTHON SOURCE LINES 219-224 .. code-block:: Python q = mapdl.queries LFT_NODE = q.node(4, 0, 0) # store node number to 'LFT_NODE' with coordinate (4,0,0) MID_NODE = q.node(6, 0, 0) # store node number to 'MID_NODE' with coordinate (6,0,0) RT_NODE = q.node(8, 0, 0) # store node number to 'RT_NODE' with coordinate (8,0,0) .. GENERATED FROM PYTHON SOURCE LINES 225-227 Retrieve nodal deflection and section stresses ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~ .. GENERATED FROM PYTHON SOURCE LINES 227-252 .. code-block:: Python # Retrieves the displacement "DEF_4" of nodes associated to # "LFT_NODE" in the X direction def_4 = mapdl.get("DEF_4", "NODE", LFT_NODE, "U", "X") # Retrieves the stress and store it "parm" of nodes associated to a # variables 'LFT_NODE','MID_NODE' and 'RT_NODE' rst_4_c1 = mapdl.get("RST_4_C1", "NODE", LFT_NODE, "S", "X") rst_6_c1 = mapdl.get("RST_6_C1", "NODE", MID_NODE, "S", "X") rst_8_c1 = mapdl.get("RST_8_C1", "NODE", RT_NODE, "S", "X") tst_4_c1 = mapdl.get("TST_4_C1", "NODE", LFT_NODE, "S", "Z") tst_6_c1 = mapdl.get("TST_6_C1", "NODE", MID_NODE, "S", "Z") tst_8_c1 = mapdl.get("TST_8_C1", "NODE", RT_NODE, "S", "Z") # Print the nodal stress solution (COMP means all stress components) mapdl.prnsol("S", "COMP") # Define a path with the name 'STRESS' and ID 2, no limits specified mapdl.path("STRESS", 2, "", 48) mapdl.ppath(1, LFT_NODE) # Define the path points using the variable 'LFT_NODE' mapdl.ppath(2, RT_NODE) # Define the path points using the variable 'RT_NODE' mapdl.plsect("S", "Z", -1) # Display the SZ stresses in a sectional plot mapdl.plsect("S", "X", -1) # Display the SX stresses in a sectional plot mapdl.prsect(-1) # Print linearized stresses .. rst-class:: sphx-glr-script-out .. code-block:: none 'PRINT LINEARIZED STRESS THROUGH A SECTION DEFINED BY PATH= STRESS DSYS= 0\n *****MAPDL VERIFICATION RUN ONLY*****\n DO NOT USE RESULTS FOR PRODUCTION\n\n ***** POST1 LINEARIZED STRESS LISTING *****\n INSIDE NODE = 1 OUTSIDE NODE = 106\n\n LOAD STEP 1 SUBSTEP= 1\n TIME= 1.0000 LOAD CASE= 0\n\n ** AXISYMMETRIC OPTION ** RHO = 0.47392E+12\n THE FOLLOWING X,Y,Z STRESSES ARE IN SECTION COORDINATES. \n\n ** MEMBRANE **\n SX SY SZ SXY SYZ SXZ\n -0.1000E+05 0.1255E-08 0.3000E+05 0.8005E-10 0.000 0.000 \n S1 S2 S3 SINT SEQV\n 0.3000E+05 0.000 -0.1000E+05 0.4000E+05 0.3606E+05\n\n ** BENDING ** I=INSIDE C=CENTER O=OUTSIDE\n SX SY SZ SXY SYZ SXZ\n I -0.1991E+05 0.1083E-08 0.1364E+05 0.000 0.000 0.000 \n C -4948. 0.1083E-09 0.1919E-07 0.000 0.000 0.000 \n O 0.1001E+05 -0.8666E-09 -0.1364E+05 0.000 0.000 0.000 \n S1 S2 S3 SINT SEQV\n I 0.1364E+05 0.000 -0.1991E+05 0.3355E+05 0.2922E+05\n C 0.000 0.000 -4948. 4948. 4948. \n O 0.1001E+05 0.000 -0.1364E+05 0.2365E+05 0.2056E+05\n\n ** MEMBRANE PLUS BENDING ** I=INSIDE C=CENTER O=OUTSIDE\n SX SY SZ SXY SYZ SXZ\n I -0.2991E+05 0.2338E-08 0.4364E+05 0.8005E-10 0.000 0.000 \n C -0.1495E+05 0.1363E-08 0.3000E+05 0.8005E-10 0.000 0.000 \n O 6.734 0.3879E-09 0.1636E+05 0.8005E-10 0.000 0.000 \n S1 S2 S3 SINT SEQV\n I 0.4364E+05 0.000 -0.2991E+05 0.7355E+05 0.6407E+05\n C 0.3000E+05 0.000 -0.1495E+05 0.4495E+05 0.3965E+05\n O 0.1636E+05 6.734 0.000 0.1636E+05 0.1636E+05\n\n ** PEAK ** I=INSIDE C=CENTER O=OUTSIDE\n SX SY SZ SXY SYZ SXZ\n I -0.3638E-11 0.2863E-08 6264. 0.1697E-07 0.000 0.000 \n C 7193. -0.9225E-09 -2245. 0.4677E-09 0.000 0.000 \n O 0.7061E-12 0.2176E-08 3632. 0.2536E-09 0.000 0.000 \n S1 S2 S3 SINT SEQV\n I 6264. 0.000 0.000 6264. 6264. \n C 7193. 0.000 -2245. 9438. 8540. \n O 3632. 0.000 0.000 3632. 3632. \n\n ** TOTAL ** I=INSIDE C=CENTER O=OUTSIDE \n SX SY SZ SXY SYZ SXZ\n I -0.2991E+05 0.5201E-08 0.4991E+05 0.1705E-07 0.000 0.000 \n C -7758. 0.4404E-09 0.2776E+05 0.5477E-09 0.000 0.000 \n O 6.734 0.2564E-08 0.1999E+05 0.3336E-09 0.000 0.000 \n S1 S2 S3 SINT SEQV TEMP\n I 0.4991E+05 0.000 -0.2991E+05 0.7982E+05 0.6984E+05 0.000 \n C 0.2776E+05 0.000 -7758. 0.3552E+05 0.3234E+05\n O 0.1999E+05 6.734 0.000 0.1999E+05 0.1999E+05 0.000' .. GENERATED FROM PYTHON SOURCE LINES 253-255 Verify the results. ~~~~~~~~~~~~~~~~~~~ .. GENERATED FROM PYTHON SOURCE LINES 255-307 .. code-block:: Python # Set target values target_def = 0.0078666 target_strss = [-30000, -7778, 0] target_tst_strss = [50000, 27778, 20000] # Fill result values res_def = def_4 res_strss = [rst_4_c1, rst_6_c1, rst_8_c1] res_tst_strss = [tst_4_c1, tst_6_c1, tst_8_c1] title = f""" ------------------- VM25 RESULTS COMPARISON --------------------- RESULTS FOR CASE p = 30,000 psi: -------------------------------- """ col_headers = ["TARGET", "Mechanical APDL", "RATIO"] row_headers = ["Displacement, in (r = 4 in)"] data = [ [target_def, res_def, abs(target_def / res_def)], ] print(title) print(pd.DataFrame(data, row_headers, col_headers)) # Radial stress results comparison row_headers = [ "Stress_r, psi (r = 4 in)", "Stress_r, psi (r = 6 in)", "Stress_r, psi (r = 8 in)", ] data = [target_strss, res_strss, np.abs(target_strss) / np.abs(res_strss)] print(pd.DataFrame(np.transpose(data), row_headers, col_headers)) # Tangential stress results comparison row_headers = [ "Stress_t, psi (r = 4 in)", "Stress_t, psi (r = 6 in)", "Stress_t, psi (r = 8 in)", ] data = [ target_tst_strss, res_tst_strss, np.abs(target_tst_strss) / np.abs(res_tst_strss), ] print(pd.DataFrame(np.transpose(data), row_headers, col_headers)) .. rst-class:: sphx-glr-script-out .. code-block:: none ------------------- VM25 RESULTS COMPARISON --------------------- RESULTS FOR CASE p = 30,000 psi: -------------------------------- TARGET Mechanical APDL RATIO Displacement, in (r = 4 in) 0.007867 0.007867 0.999992 TARGET Mechanical APDL RATIO Stress_r, psi (r = 4 in) -30000.0 -29908.046900 1.003075 Stress_r, psi (r = 6 in) -7778.0 -7757.541500 1.002637 Stress_r, psi (r = 8 in) 0.0 6.734007 0.000000 TARGET Mechanical APDL RATIO Stress_t, psi (r = 4 in) 50000.0 49908.0469 1.001842 Stress_t, psi (r = 6 in) 27778.0 27757.5410 1.000737 Stress_t, psi (r = 8 in) 20000.0 19993.2656 1.000337 .. GENERATED FROM PYTHON SOURCE LINES 308-310 Finish the post-processing processor. ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~ .. GENERATED FROM PYTHON SOURCE LINES 310-312 .. 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 313-315 Set a new title for the analysis ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~ .. GENERATED FROM PYTHON SOURCE LINES 315-320 .. code-block:: Python mapdl.title("VM25 Stresses in a Long Cylinder - Rotation About Axis") # Resume the Finite Element (FE) "MODEL" save previously mapdl.resume("MODEL") .. rst-class:: sphx-glr-script-out .. code-block:: none RESUME ANSYS DATA FROM FILE NAME= *** MAPDL GLOBAL STATUS *** TITLE = VM25 Stresses in a Long Cylinder 1 ELEM TYPES DEFINED MAX ELEM TYPE NUMBER = 1 28 ELEMENTS DEFINED MAX ELEMENT NUMBER = 28 117 NODES DEFINED MAX NODE NUMBER = 222 1 MATERIALS DEFINED MAX MATERIAL NUMBER = 1 0 REAL CONSTS DEFINED MAX REAL CONST NUMBER = 0 0 SECTIONS DEFINED MAX SECTION NUMBER = 0 0 COORD SYS DEFINED MAX COORD SYS NUMBER = 0 1 COUPLING EQNS DEFINED MAX COUPLING EQN NUMBER = 1 ACTIVE COORDINATE SYSTEM = 0 (CARTESIAN) NUMBER OF DEFINED NODAL CONSTRAINTS = 29 NUMBER OF DEFINED NODAL LOADS = 0 NUMBER OF DEFINED ELEM SURFACE LOADS = 4 NUMBER OF DEFINED ELEM BODY LOADS = 0 NUMBER OF DEFINED NODE BODY FORCES = 0 INITIAL JOBNAME = CURRENT JOBNAME = .. GENERATED FROM PYTHON SOURCE LINES 321-324 Solve ~~~~~ Enter solution mode and solve the system. .. GENERATED FROM PYTHON SOURCE LINES 324-329 .. code-block:: Python mapdl.slashsolu() # Print all results mapdl.outpr("", "ALL") .. rst-class:: sphx-glr-script-out .. code-block:: none PRINT BASI ITEMS WITH A FREQUENCY OF ALL FOR ALL APPLICABLE ENTITIES .. GENERATED FROM PYTHON SOURCE LINES 330-335 Define boundary conditions ~~~~~~~~~~~~~~~~~~~~~~~~~~ Apply a displacement boundary condition in the vertical direction (UY) to all nodes. Rotate the cylinder with an angular velocity of 1000 rad/sec. Also, apply dummy pressure for surface printout. Then exit solution processor. .. GENERATED FROM PYTHON SOURCE LINES 335-356 .. code-block:: Python mapdl.nsel( "S", "LOC", "Y", 0 ) # Select nodes located at Y=0 to prevent rigid body motion mapdl.nsel("R", "LOC", "X", 4) # Select nodes located at X=4 # Displace all nodes in the Y-direction mapdl.d("ALL", "UY") mapdl.nsel("S", "LOC", "X", 4) # Select nodes located at X=4 # Apply a small pressure to allow stress printout mapdl.sf("", "PRES", 1e-10) mapdl.nsel("ALL") # Select all nodes # Rotate the cylinder with an angular velocity of 1000 RAD/SEC mapdl.omega("", 1000) # Solve the problem in load step 2 - centrifugal loading 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 357-360 Post-processing ~~~~~~~~~~~~~~~ Enter post-processing. Compute displacement and stress components. .. GENERATED FROM PYTHON SOURCE LINES 360-362 .. code-block:: Python mapdl.post1() .. rst-class:: sphx-glr-script-out .. code-block:: none *****MAPDL VERIFICATION RUN ONLY***** DO NOT USE RESULTS FOR PRODUCTION ***** MAPDL RESULTS INTERPRETATION (POST1) ***** USE LAST SUBSTEP ON RESULT FILE FOR LOAD CASE 0 SET COMMAND GOT LOAD STEP= 1 SUBSTEP= 1 CUMULATIVE ITERATION= 1 TIME/FREQUENCY= 1.0000 TITLE= VM25 Stresses in a Long Cylinder .. GENERATED FROM PYTHON SOURCE LINES 363-365 Inline functions in PyMAPDL to query node ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~ .. GENERATED FROM PYTHON SOURCE LINES 365-371 .. code-block:: Python q = mapdl.queries LFT_NODE = q.node(4, 0, 0) # store node number to 'LFT_NODE' with coordinate (4,0,0) XI_NODE = q.node( 5.43, 0, 0 ) # store node number to 'XI_NODE' with coordinate (5.43,0,0) .. GENERATED FROM PYTHON SOURCE LINES 372-374 Retrieve nodal deflection and section stresses ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~ .. GENERATED FROM PYTHON SOURCE LINES 374-380 .. code-block:: Python rst_4_c2 = mapdl.get("RST_4_C2", "NODE", LFT_NODE, "S", "X") tst_4_c2 = mapdl.get("TST_4_C2", "NODE", LFT_NODE, "S", "Z") rst_x_c2 = mapdl.get("RST_X_C2", "NODE", XI_NODE, "S", "X") tst_x_c2 = mapdl.get("TST_X_C2", "NODE", XI_NODE, "S", "Z") .. GENERATED FROM PYTHON SOURCE LINES 381-383 Verify the results. ~~~~~~~~~~~~~~~~~~~ .. GENERATED FROM PYTHON SOURCE LINES 383-419 .. code-block:: Python # Set target values target_strss = [0, 4753] target_tst_strss = [40588, 29436] # Fill result values res_strss = [rst_4_c2, rst_x_c2] res_tst_strss = [tst_4_c2, tst_x_c2] title = f""" RESULTS FOR CASE Rotation = 1000 rad/sec: ----------------------------------------- """ col_headers = ["TARGET", "Mechanical APDL", "RATIO"] # Radial stress results comparison row_headers = ["Stress_r, psi (r = 4 in)", "Stress_r, psi (r = 5.43 in)"] data = [target_strss, res_strss, np.abs(target_strss) / np.abs(res_strss)] print(title) print(pd.DataFrame(np.transpose(data), row_headers, col_headers)) # Tangential stress results comparison row_headers = ["Stress_t, psi (r = 4 in)", "Stress_t, psi (r = 5.43 in)"] data = [ target_tst_strss, res_tst_strss, np.abs(target_tst_strss) / np.abs(res_tst_strss), ] print(pd.DataFrame(np.transpose(data), row_headers, col_headers)) .. rst-class:: sphx-glr-script-out .. code-block:: none RESULTS FOR CASE Rotation = 1000 rad/sec: ----------------------------------------- TARGET Mechanical APDL RATIO Stress_r, psi (r = 4 in) 0.0 50.366142 0.000000 Stress_r, psi (r = 5.43 in) 4753.0 4957.113530 0.958824 TARGET Mechanical APDL RATIO Stress_t, psi (r = 4 in) 40588.0 41671.1641 0.974007 Stress_t, psi (r = 5.43 in) 29436.0 29740.9404 0.989747 .. GENERATED FROM PYTHON SOURCE LINES 420-422 Finish the post-processing processor. ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~ .. GENERATED FROM PYTHON SOURCE LINES 422-424 .. 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 425-427 Stop MAPDL. ~~~~~~~~~~~ .. GENERATED FROM PYTHON SOURCE LINES 427-428 .. code-block:: Python mapdl.exit() .. rst-class:: sphx-glr-timing **Total running time of the script:** (0 minutes 0.808 seconds) .. _sphx_glr_download_verif-manual_vm-025.py: .. only:: html .. container:: sphx-glr-footer sphx-glr-footer-example .. container:: sphx-glr-download sphx-glr-download-jupyter :download:`Download Jupyter notebook: vm-025.ipynb ` .. container:: sphx-glr-download sphx-glr-download-python :download:`Download Python source code: vm-025.py ` .. only:: html .. rst-class:: sphx-glr-signature `Gallery generated by Sphinx-Gallery `_