Resources

Transforming Stress into a User-Defined Coordinate Frame in LS-DYNA

Written by Mark Lytell | Aug 8, 2024 7:33:05 PM

Introduction

Often, solving an engineering problem using simulation tools such as LS-DYNA or Ansys Mechanical requires that tensorial results be compared to an analytical solution.  One such example is a pressurized cylinder in which the circumferential, or hoop, stress is the main stress measure of interest.  If such a simulation is conducted on a geometry that has complex end treatments, for example, one would like to use the hoop stress, say in the center of the vessel, to validate that the simulation is providing accurate results in areas of more complex geometry where analytic results do not exist.

In most modern postprocessors, stress components may be fringe plotted in either cartesian coordinates, principal directions, or even in cylindrical or spherical coordinates.  However, it may necessary to rotate the results a tensorial result into an orientation that is not available "out-of-the-box" in said postprocessors.  In this article, we will provide a verifiable example of tensor transformation using a Python 3.11 script to read stress fringe plot results for a pressurized cylinder from LS-PrePost in Cartesian coordinates and subsequently rotate the stress tensor into cylindrical coordinates to compute the hoop stress.

Case Study

The model that we will use to illustrate the concept is a d = 100 mm mean diameter x l = 150 mm long x t = 2mm thick cylinder subject to an internal pressure of p = 4 MPa solved using LS-DYNA.  A textbook image1 of the problem is shown here:  

FE Model

The material for the cylinder is the standard linear isotropic elastic Structural Steel from Workbench Engineering Data, and the mesh is composed of 2.5 mm sized shell elements, using default section properties.  The mesh sizing provides a uniform grid.

 

Loads and Boundary Conditions

The internal pressure, p, was chosen to be 4 MPa so that the hoop stress would equal 100 MPa according to the well-known formula

where r = mean cylinder radius and t = wall thickness.  Additionally, isostatic constraints (a.k.a. 3-2-1 constraints) are applied to the vertices shown in yellow in order to provide stress-free constraint.

Analysis Settings

The analysis type is set to implicit with a pseudo end-time of 0.01 s.  Implicit analysis is chosen for speedy analysis run time.

Simulation Results

Using a cylindrical coordinate system in Ansys Mechanical, we can easily postprocess the (unaveraged) hoop stress, resulting in 100 MPa, thus matching theory.

This result shows that our simulation matches theory.

Reconstructing the Hoop Stress using Python

Now that the hoop stress results have been found using Ansys Workbench LS-DYNA, we will reconstruct the same hoop stress results by processing the Cartesian stress components that are output from LS-PrePost.

Exporting Cartesian Stress Components

To obtain the data necessary to rotate the stress tensor, all six components must be exported from LS-PrePost.  The process to export each stress component will be illustrated for the X-component and must be repeated for the Y-, Z-, XY-, YZ- and ZX- components:

  1. Load the binary D3Plot file as usual and set the results time to the end state.
  2. Plot the X-stress in a fringe plot. 
  3. Output the results using the settings shown below, naming the file "xstress", for example.
  4. Hit Done and repeat for the other 5 stress components.
  5. Make sure that all six of the files are in the same directory and the python script.

Output File Format

The format of the output files is in LS-DYNA keyword format which is a fixed-width file format.  The Python script parses each section of the file taking advantage of the fixed-width structure.

*ELEMENT_SHELL Section

The element connectivity table is within the *ELEMENT_SHELL section of the file, a small except of which is shown here:

The column definitions and widths are as follows:

Column Data Width
1 Element ID 8
2 Part ID 8
3 Node 1 8
4 Node 2 8
5 Node 3 8
6 Node 4 8

*NODE Section

Similarly, the nodal coordinates are given in the *NODE section, an excerpt of which is shown here:

The column definitions and widths are as follows:

Column Data Width
1 Node ID 8
2 X-Coordinate 16
3 Y-Coordinate 16
4 Z-Coordinate 16

 

$SHELL_ELEMENT_RESULTS Section

The $SHELL_ELEMENT_RESULTS section provides the output result at every integration point of the element.  In this case we are using shell elements with a single integration point.  The excerpt of this section is shown here:

The column definitions and widths are as follows: 

Column Data Width
1 Element ID 8
2 Element Result 16

 

$SHELL_ELEMENT_CENTROID Section

The $SHELL_ELEMENT_CENTROID section provides the centroid and volume of each element in the deformed state.  The excerpt of this section is shown here:

The column definitions and widths are as follows:

Column Data Width
1 Element ID 10
2 X-Coordinate 16
3 Y-Coordinate 16
4 Z-Coordinate 16
5 Volume 16

 

Python Algorithm

The following pseudocode explains the algorithm used to transform the stress:

  1. Read each of the six stress files and parse the sections to obtain the element table, nodal coordinates, element results (unaveraged stress) and element centroidal coordinates (in the deformed configuration).
  2. Merge the element table with the nodal coordinates and element results.
  3. Form the elemental result (stress) tensor, S.
  4. Calculate the normal unit vector to each element at its centroid, using the cross product from two element edges according to the right-hand rule.  We call this the element Z-direction. 
    1. For the particular shell elements used in this example, the two edges were formed from node 1 to 4 (1->4), and from node 1 to 2 (1->2).  Thus 1->4 x 1->2 provides the outward element normal, followed by a normalization to make it a unit vector.
  5. Since we are using cylindrical coordinates, we form the element X-direction unit vector by crossing the global Z-axis with the element Z-direction vector.
    1. This elemental X-direction vector is not orthogonal to the element normal vector, so we orthonormalize it by using Gram-Schmidt orthogonalization.
  6. Having the elemental X- and Z-unit vectors, we take Z x X to get the elemental Y-direction unit vector, thus completing our orthonormal set.
  7. Form the transformation matrix Q using the X-, Y-, and Z- elemental unit vectors as columns.
  8. Rotate the stress tensor using the tensor transformation rule, Sloc = Q'SQ, where ()' indicates transpose.
  9. Write the results to an output .csv file.

Python Computed Hoop Stress

An excerpt from the results of the Python script is shown here, where the image on the left displays the element stress in Cartesian coordinates; and the image on the right shows the corresponding hoop stress, S11_loc, resulting from the transformation and the components of the transformation matrix Q:

   

 

From the image on the rights, we see that the computed hoop stress, S11_loc, matches very well to the output from Ansys Mechanical.

Conclusions

Although the hoop stress is easily obtained from within Workbench Ansys LS-DYNA or LS-PrePost (as the 1st principal stress), the steps in this article show how to accomplish the rotation of tensorial results for more general orientations by transforming results in Cartesian coordinates.  The procedure outlined above can also be used in a larger automation framework using the tools of PyAnsys.

Going Further

  1. Modify the python script to look at stress (or strain) in other orientations or apply to your own model to obtain stress (or strain) results of interest.
  2. Modify the code for use with solid elements or shell elements with multiple integration points.
  3. Utilize this script along with PyAnsys to automate the process.
Downloadable Items

Ansys 2024 R2 Workbench Archive, python code and results

References

1.  Pressure vessel image taken from https://engineersfield.com/stress-and-strain/