Hi all,

This topic has come up a number of times in recent conversations, so I thought I would share the general workflow for importing VTK data directly into python (numpy arrays).

Below is a snippet of python code that should largely work out of the box for 3D cartesian models, where all one has to do is provide the full path and name of the .pvtu in a string. Please note that lines in bold are comments, and you would need to add a # symbol before them in python. I excluded them as they produced odd formatting in the post.

Please let me know if you have any comments, questions, or suggestions.

Cheers,

John

Load modules

`import numpy np; import vtk as vtk; from vtk.util import numpy_support`

Load vtu data (pvtu directs to vtu files)

`reader = vtk.vtkXMLPUnstructuredGridReader()`

`reader.SetFileName('full_path_and_name_of_pvtu_file')`

`reader.Update()`

Example for how to set the file name:

`reader.SetFileName('/home/naliboff/solution-00000.pvtu')`

Get the coordinates of nodes in the mesh

`nodes_vtk_array= reader.GetOutput().GetPoints().GetData()`

Convert nodal vtk data to a numpy array

`nodes_numpy_array = vtk.util.numpy_support.vtk_to_numpy(nodes_vtk_array)`

Extract x, y and z coordinates from numpy array

x,y,z= nodes_numpy_array[:,0] , nodes_numpy_array[:,1] , nodes_numpy_array[:,2]

Determine the number of scalar fields contained in the .pvtu file

`number_of_fields = reader.GetOutput().GetPointData().GetNumberOfArrays()`

Determine the name of each field and place it in an array.

`field_names = []`

`for i in range(number_of_fields):`

`field_names.append(reader.GetOutput().GetPointData().GetArrayName(i))`

Determine the index of the field strain_rate

`idx = field_names.index("strain_rate")`

Example: Extract values of strain_rate

`field_vtk_array = reader.GetOutput().GetPointData().GetArray(idx)`

`strain_rate = numpy_support.vtk_to_numpy(field_vtk_array)`