A few months ago, I decided it was time for me to buckle down and learn PyAnsys. We also needed a replacement for our Ansys ACT tool that output results in a 3MF file format, which became out of date. This series is the result of that effort, sharing with you all that I learned along the way.
In Part 1 of this series, we covered:
- How to install and set up PyAnsys on your machine
- Create a PyAnsys script that reads an Ansys Mechanical result file
- Pull mesh and results from the file
- Plot results to the screen using Python
- Write a VTK file of the results
If you have not run PyAnsys before, then please go back to the first post and make sure you add the PyAnsys modules to your Python install using pip as shown. Just to whet your appetite, here are some of the things we have plotted on our Stratasys J55 and J850 Polyjet printers:
Let’s Use PyAnsys to Write a File That We Can Use to 3D Print Results
Creating a VTK file lets you use other tools to convert to other formats, but when push comes to shove, we want to write 3D files directly and skip the extra step. We also need to do more than create files with equivalent stress from the first solution set in the file.
That is why, in Part 2, we are going to:
- Parameterize things
- Support multiple result types
- Distort the mesh
- Write a WRL/VRML file we can use for 3D Printing
We will step through each bit of code for each of these steps. The zip file at the bottom of this post will have it all in one file, but you can also copy and paste as we go.
1: Clean Things Up and Add Parameters
In this first section of the code, we will do pretty much the same thing we did in that section for Part 1 of these tutorials. The big difference is some new modules and more parameters to support new capabilities. As you can see we will use pyvista and numpy. You will see why later on.
Then we create new parameters for other result types, how much distortion scaling we want, and support for more than one output format. Nothing fancy here.
#### SECTION 1 #### #1.1: Get all modules loaded from ansys.dpf import post from ansys.dpf import core as dpf from ansys.dpf.core import operators as coreops from ansys.dpf.core import fields_container_factory from ansys.dpf.core.plotter import DpfPlotter # For the new capabilities, we need pyvista, and numpy along with os import pyvista as pv import numpy as np import os # 1.2: Specify variables for filenames and directory # Make sure you change this to your directory and RST file # Note that Python converst / to \ on windows so you don't have to. rstFile = "twrtest1-str.rst" rstDir = "C:/Users/eric.miller/OneDrive - PADT/Development/pyAnsys/twrtest/" outroot= "twrtest1" #1.3: Now define some parameters we will need to # get other result types, scale distortion, and write to more than one format # They ansys value to retrieve. Valid result types are: # ux, uy,uz, usum, seqv, s1,s2,s3,sx,sy,sz, sxy, sxz, syz, tmp rsttype = "usum" # A scale factor for the max deflection value on the # model as a percent of the largest value of the bounding box pcntdfl = 5 # The 3d result file format. For now we will check for vtk or wrl outtype = "wrl" # The result step/mode to grab rstnum = 1 #1.4: Build the output file name from the file root and # the result type and file type outfname = outroot + "-" + rsttype + "-" + str(rstnum) + "." + outtype #1.5: Let the user know what you have print ("==========================================================================") print ("Starting the translation of your Ansys results to a different file format") print (" Input file: ", rstFile) print (" Directory: ", rstDir) print (" Output file: ", outfname) print (" Solution step/mode: ", rstnum) print (" Result type: ", rsttype) print (" Percent deflection distortion: ", pcntdfl) print ("-------------------------------------------------------------------------")
2: Open the file, get mesh, and get the right stress value
As mentioned in Part 1, opening the file and getting the mesh is so easy, you almost feel like you left something out. Things get more complicated with a bunch of if statements looking for the type of result to read. You can see that there is a displacement, stress, and thermal object and each has a slightly different way of grabbing value. This section could use some error checking if you want to add it.
Note how when we grab the temperature, displacement, and stress values, we are now specifying a result number. It is not by load-step, sub-step. Whatever you want, put it into rstval to use in making the plot and the file.
If this section looks confusing, go back and review section 2 from the first part of this tutorial.
#### SECTION 2 #### #2.1: Open the result file and load the solution and mesh os.chdir(rstDir) path = os.getcwd() mysol = post.load_solution(rstFile) mymesh = mysol.mesh #2.2: Grab the results from the specified solution step, # then get the request result type (thermal or displacement and stress) print ("++ Getting result information from file") # Displacement, stress, and thermal for the solution set number if rsttype == "tmp": thrm1 = mysol.temperature(set=rstnum) else: dsp1 = mysol.displacement(set=rstnum) str1 = mysol.stress(set=rstnum) #2.3: Use an if statement to pull in the specific result values # Start with displacement if rsttype == "u": rstval = dsp1.vector.result_fields_container elif rsttype == "ux": rstval = dsp1.x.result_fields_container elif rsttype == "uy": rstval = dsp1.y.result_fields_container elif rsttype == "uz": rstval = dsp1.z.result_fields_container elif rsttype == "usum": rstval = dsp1.norm.result_fields_container # Now check for stresses elif rsttype == "seqv": rstval = str1.von_mises.result_fields_container elif rsttype == "s1": rstval = str1.principal_1.result_fields_container elif rsttype == "s2": rstval = str1.principal_2.result_fields_container elif rsttype == "s3": rstval = str1.principal_3.result_fields_container elif rsttype == "sx": rstval = str1.xx.result_fields_container elif rsttype == "sy": rstval = str1.yy.result_fields_container elif rsttype == "sz": rstval = str1.zz.result_fields_container elif rsttype == "sxy": rstval = str1.xy.result_fields_container elif rsttype == "sxz": rstval = str1.xz.result_fields_container elif rsttype == "syz": rstval = str1.yz.result_fields_container # Last, do temperatures elif rsttype == "tmp": rstval = thrm1.scalar.result_fields_container
3. Grab deflections and distort the mesh
Here is our first completely new section. Since we want to show results on a distorted mesh, we have to get the displacement values, multiply them by a scale factor, then move the nodes to their deflected locations. As with everything in Python, there are dozens of methods and functions to help with this. For this code, I kept things kind of simple.
First we figure out how big the model is, then get the max displacement for the result set we are working with. Knowing that value, we can now calculate a scale factor using the user’s input. We use the math operators in PyAnsys to first scale the deflection values, then to add those values to the nodal position.
#### SECTION 3 ##### # 3.1: If this is thermal, just copy the # undistored mesh to the variables we will us to plot and write if rsttype == "tmp": dflmesh = mymesh.deep_copy() newcoord = dflmesh.nodes.coordinates_field else: #3.2: Not thermal so get info needed to calculate a distorted mesh print("++ Calculating deflection distortion") # Grab the total distortion at each node usumval = dsp1.vector.result_fields_container # Get model extents. Feed the min_max operator "mymesh" extop = dpf.operators.min_max.min_max(mymesh.nodes.coordinates_field) coordmin = extop.outputs.field_min() coordmax = extop.outputs.field_max() # Look at the X, Y, Z min and max coordinate # value and find the biggest demention of the mesh # There is probably a clever python way to do this in one line dltmax =0.0 for i in range(0,3): if (coordmax.data[i]-coordmin.data[i]) > dltmax: dltmax = coordmax.data[i]-coordmin.data[i] # Get the maximum deflection value from usumval umaxop = dpf.operators.min_max.min_max(field=usumval) umax = 0.0 for i in range(0,3): if umaxop.outputs.field_max().data[i] > umax: umax = umaxop.outputs.field_max().data[i] # Calculate a scale factor that is the specified # percentage of the max deflection devided by the max size dflsf = pcntdfl/100.0 sclfact = dflsf*dltmax/umax #3.3: Scale the deflection values then distort the nodal coordinates # Get a copy of the mesh to distort dflmesh = mymesh.deep_copy() newcoord = dflmesh.nodes.coordinates_field # Scale the deflection field by the scale factor sclop = dpf.operators.math.scale(usumval,float(sclfact)) scaledusum = sclop.outputs.field() # Add the scaled deflection field to the nodal positions addop = dpf.operators.math.add(newcoord, scaledusum) tempcoord = addop.outputs.field() # Overwrite the nodal positions of dflmesh with the deflected ones newcoord.data = tempcoord.data
4. Create a plot to the screen to verify what you have
No changes here from before, just build a simple window and show it. To exit the window, close it or press q.
#### SECTION 4 #### #4.1: create a Pyansys plotter object first print ("++ Making plot (press q to exit)") plt = DpfPlotter() #4.2: add the result value to the deflected mesh and the plot plt.add_field(rstval,dflmesh,show_edges=False) plt.show_figure( show_axes=True, parallel_projection=False, background="#aaaaaa", )
Here is what I end up with for a thermal model I use (rth included in zip file)
5: Write out a VTK or WRL File
Here is where we get a bit fancy. First we use an if/elif to use the write code for either VTK or WRL.
For that WRL file, PyAnsys doesn’t have a tool to write to WRL, so we have to use pyvista. And Pyvista can’t take the PyAnsys results directly. We have to move the nodal or elemental result info from the result objects into what pyvista calls a grid as either point data or cell data. Once we do that, we can feed that grid to the pyvista plotter tools, and use export_vrml to write the file.
#### SECTION 5 #### #5.1: I writing a VTK file, then do the same thing we did for the first post # setup a pyansys vtk_export operator, add mesh, filename, # and result values, then run teh operator print ("++ Making output file") if outtype == "vtk": vtkop = coreops.serialization.vtk_export() vtkop.inputs.mesh.connect(dflmesh) vtkop.inputs.file_path.connect(outfname) vtkop.inputs.fields1.connect(rstval) vv = vtkop.run() print(vv) print ("++ Results written to: " + outfname) #5.2: If doing a WRL, use pyVista to create a VRML file for visualization # and 3D Printing # The pyvista module does not work with the Ansys mesh object, # so you have to do some fancy converstion to get the mesh into a # grid it can use # It also doesn't automatically handle element values vs nodal value, # so we have to handle that elif outtype == "wrl": # Get the result field and the deflected mesh into some variables. field = rstval meshed_region = dflmesh # We are going to need the name of the field in "field" # in order to grab it later fieldName = '_'.join(field.name.split("_")[:-1]) # Look in the field to see if the location of the results is nodal or elemental if field.location == dpf.locations.nodal: mesh_location = meshed_region.nodes elif field.location == dpf.locations.elemental: mesh_location = meshed_region.elements else: raise ValueError( "Only elemental or nodal location are supported for plotting." ) # Next we are going to use a mask to get a grid them assign resutls to # either the element or nodal data overall_data = np.full(len(mesh_location), np.nan) ind, mask = mesh_location.map_scoping(field.scoping) overall_data[ind] = field.data[mask] grid = meshed_region.grid if field.location == dpf.locations.nodal: grid.point_data[fieldName] = overall_data elif field.location == dpf.locations.elemental: grid.cell_data[fieldName] = overall_data #5.3: Now that we have a grid and values on the grid, # use the pyvista plotter operator to make the file objplt = pv.Plotter() objplt.add_mesh(grid) #close, data is not lined up with nodes objplt.export_vrml(outfname) print ("++ Results written to: " + outfname) #5.4: Close out print ("---------------------------------------------------------------------------") print ("Done")
And, here is what it looks like in Stratasys GrabCad, ready to go into our J850 Prime system.
And what it looks like in the real world:
Next Steps for Our PyAnsys Journey
You can find the Python script, result files for my test model and the housing shown above, and some WRL files I made in this zip file:
As you can see, the code as is works if you need VTK or WRL, and you just need to use it occasionally. To make things more useful, in the next part of this series, we will add a graphical user interface (GUI) and support a few more file formats.
When you are ready, you can move on to Part 3.
If this stuff gets you interested in doing more automation with PyAnsys scripting, but you don’t have the time to do it yourself, contact PADT at firstname.lastname@example.org, and we can talk about how we can help.