Programmable Source: Vorticity and gradient from numpy data

ParaView Support
1

Hello,

I have a programmable source that creates a Delaunay2D mesh from Numpy arrays. Now I want to calculate the vorticity from the vector field but “algs.vorticity(Vector_Speed)” do not work because “Vector_Speed” is a numpy object. The same happens with the “algs.gradient(Vx)” (because “Vx” is a numpy array)

Since I really want to use the “algs.vorticity” function, I tried to convert the numpy vector (i.e., “Vector_Speed”) into a vtk object and then to wrap this element but I cannot make it work. Please, can you provide me some guidelines to solve this issue?

Thanks in advance!

PS: the programmable source must be a vtkPolyData (vtkImage is not an option for my application)

Miguel

# MODULES LOAD:
import numpy as np
import vtk
from vtk.numpy_interface import algorithms as algs
from vtk.numpy_interface import dataset_adapter as dsa
import paraview.vtk.util.numpy_support as vnp


# NUMPY INPUT DATA:
X=np.array([0,1,2,0,1,2])   # X-coordinate
Y=np.array([0,0,0,1,1,1])
Z=np.array([0,0,0,0,0,0])
Vx=np.array([0,1,2,0,1,2])   #X-velocity
Vy=np.array([0,0,0,0,0,0])
Vz=np.array([0,0,0,0,0,0])


# DELAUNAY PLANE CREATION:
points = vtk.vtkPoints()  						
 
for i in range(X.size):
	points.InsertNextPoint(X[i],Y[i],Z[i])
 
aPolyData = vtk.vtkPolyData()
aPolyData.SetPoints(points)
 
aCellArray = vtk.vtkCellArray()
 
boundary = vtk.vtkPolyData()
boundary.SetPoints(aPolyData.GetPoints())
boundary.SetPolys(aCellArray)
delaunay = vtk.vtkDelaunay2D()
if vtk.VTK_MAJOR_VERSION <= 5:
	delaunay.SetInput(aPolyData.GetOutput())
	delaunay.SetSource(boundary)
else:
	delaunay.SetInputData(aPolyData)
	delaunay.SetSourceData(boundary)

delaunay.Update()

self.GetOutput().ShallowCopy(delaunay.GetOutput())



# GRADIENT AND VORTICITY:
Vector_Speed = algs.make_vector(Vx,Vy,Vz)
#VORTICITY = algs.vorticity(Vector_Speed)   #This fuction is not working here because Vector_Speed is a numpy element


#CONVERTING THE NUMPY VECTOR INTO A VTK ARRAY
vtkdata = vnp.numpy_to_vtk(Vector_Speed.ravel(), deep=True, array_type=vtk.VTK_FLOAT)
vtkdata.SetNumberOfComponents(3)
vtkdata.SetName("Vector_Speed")


#HERE I TRY TO WRAP THE VTK ELEMENT (but fails!)
WrapData = dsa.WrapDataObject(vtkdata.GetOutput())
Vector_Speed = WrapData.PointData["Vector_Speed"]


# THIS IS WHAT I WANT TO CALCULATE
VORTICITY = algs.vorticity(Vector_Speed)
Vx_GRADIENT = algs.gradient(Vx)
2
3

Thanks for the info Mathieu !

I was wrapping the wrong vtk object. Here is the good version of the code for those that may have the same problem (Programmable source - vtk PolyData output).

Best regards

Miguel

# MODULES LOAD:
import numpy as np
import vtk
from vtk.numpy_interface import algorithms as algs
from vtk.numpy_interface import dataset_adapter as dsa


# NUMPY INPUT DATA:
X=np.array([0,1,2,0,1,2,0,1,2])   # X-coordinate
Y=np.array([0,0,0,1,1,1,2,2,2])
Z=np.array([0,0,0,0,0,0,0,0,0])
Vx=np.array([0,1,0,0,0,0,0,-1,0])   #X-velocity
Vy=np.array([0,0,0,-1,0,1,0,0,0])
Vz=np.array([0,0,0,0,0,0,0,0,0])


# DELAUNAY PLANE CREATION:
points = vtk.vtkPoints()  						
 
for i in range(X.size):
	points.InsertNextPoint(X[i],Y[i],Z[i])
 
aPolyData = vtk.vtkPolyData()
aPolyData.SetPoints(points)
 
aCellArray = vtk.vtkCellArray()
boundary = vtk.vtkPolyData()
boundary.SetPoints(aPolyData.GetPoints())
boundary.SetPolys(aCellArray)
delaunay = vtk.vtkDelaunay2D()
if vtk.VTK_MAJOR_VERSION <= 5:
	delaunay.SetInput(aPolyData.GetOutput())
	delaunay.SetSource(boundary)
else:
	delaunay.SetInputData(aPolyData)
	delaunay.SetSourceData(boundary)

delaunay.Update()
self.GetOutput().ShallowCopy(delaunay.GetOutput())

# WRAP THE VTK ELEMENT:
WrapData = dsa.WrapDataObject(delaunay.GetOutput())
WrapData.PointData.append(Vx, "Vx")
WrapData.PointData.append(Vy, "Vy")
WrapData.PointData.append(Vz, "Vz")
Vx = WrapData.PointData["Vx"]
Vy = WrapData.PointData["Vy"]
Vz = WrapData.PointData["Vz"]


# GRADIENT AND VORTICITY:
Vector_Speed = algs.make_vector(Vx,Vy,Vz)
VORTICITY = algs.vorticity(Vector_Speed)
Vx_GRADIENT = algs.gradient(Vx)



# OUTPUT DATA:
output.PointData.append(Vector_Speed, "Vector_Speed")
output.PointData.append(VORTICITY, "VORTICITY")
output.PointData.append(Vx_GRADIENT, "Vx_GRADIENT")