Python coprocessing example

This example is used to demonstrate how the co-processing library can be used with a python based simulation code. Note that this example requires MPI to be available on your system as well as pyMPI to initialize and finalize MPI from the python script. The executable takes in a python co-processing script and a number of time steps to be run for. Note to remember to set your system environment properly. See [[1]] for details.

Serial python driver code

<source lang="python"> import sys if len(sys.argv) != 3:

   print "command is 'python <python driver code> <script name> <number of time steps>'"
   sys.exit(1)

import paraview import paraview.vtk as vtk

def coProcess(grid, time, step, scriptname):

   import vtkCoProcessorPython # import libvtkCoProcessorPython for older PV versions
   if scriptname.endswith(".py"):
       scriptname = scriptname.rstrip(".py")
   try:
       cpscript = __import__(scriptname)
   except:
       print 'Cannot find ', scriptname, ' -- no coprocessing will be performed.'
       return

   datadescription = vtkCoProcessorPython.vtkCPDataDescription()
   datadescription.SetTimeData(time, step)
   datadescription.AddInput("input")
   cpscript.RequestDataDescription(datadescription)
   inputdescription = datadescription.GetInputDescriptionByName("input")
   if inputdescription.GetIfGridIsNecessary() == False:
       return

   inputdescription.SetGrid(grid)
   inputdescription.SetWholeExtent(grid.GetWholeExtent())
   cpscript.DoCoProcessing(datadescription)

try:

   numsteps = int(sys.argv[2])

except ValueError:

   print 'the last argument should be a number, setting the number of time steps to 10'
   numsteps = 10

for step in range(numsteps):

   # assume simulation time starts at 0
   time = step/float(numsteps)

   # create the input to the coprocessing library.  normally
   # this will come from the adaptor
   imageData = vtk.vtkImageData()
   imageData.SetOrigin(0, 0, 0)
   imageData.SetSpacing(.1, .1, .1)
   imageData.SetExtent(0, 10, 0, 12, 0, 12)
   imageData.SetWholeExtent(imageData.GetExtent())
   pointArray = vtk.vtkDoubleArray()
   pointArray.SetNumberOfTuples(imageData.GetNumberOfPoints())
   for i in range(imageData.GetNumberOfPoints()):
       pointArray.SetValue(i, i)
   pointArray.SetName("pointData")
   imageData.GetPointData().AddArray(pointArray)

   # "perform" coprocessing.  results are outputted only if
   # the passed in script says we should at time/step
   coProcess(imageData, time, step, sys.argv[1])

</source>

Sample coprocessing script

<source lang="python"> try: paraview.simple except: from paraview.simple import *

def RequestDataDescription(datadescription):

   "Callback to populate the request for current timestep"
   timestep = datadescription.GetTimeStep()
   input_name = 'input'
   if (timestep % 1 == 0) :
       datadescription.GetInputDescriptionByName(input_name).AllFieldsOn()
       datadescription.GetInputDescriptionByName(input_name).GenerateMeshOn()
   else:
       datadescription.GetInputDescriptionByName(input_name).AllFieldsOff()
       datadescription.GetInputDescriptionByName(input_name).GenerateMeshOff()

def DoCoProcessing(datadescription):

   "Callback to do co-processing for current timestep"
   cp_writers = []
   timestep = datadescription.GetTimeStep()

   grid = CreateProducer( datadescription, "input" )
   ImageWriter1 = CreateWriter( XMLPImageDataWriter, "input_grid_%t.pvti", 1, cp_writers )

   for writer in cp_writers:
       if timestep % writer.cpFrequency == 0:
           writer.FileName = writer.cpFileName.replace("%t", str(timestep))
           writer.UpdatePipeline()

   # explicitly delete the proxies -- we do it this way to avoid problems with prototypes
   tobedeleted = GetNextProxyToDelete()
   while tobedeleted != None:
       Delete(tobedeleted)
       tobedeleted = GetNextProxyToDelete()

def GetNextProxyToDelete():

   iter = servermanager.vtkSMProxyIterator()
   iter.Begin()
   while not iter.IsAtEnd():
       if iter.GetGroup().find("prototypes") != -1:
           iter.Next()
           continue
       proxy = servermanager._getPyProxy(iter.GetProxy())
       proxygroup = iter.GetGroup()
       iter.Next()
       if proxygroup != 'timekeeper' and proxy != None and proxygroup.find("pq_helper_proxies") == -1 :
           return proxy

   return None

def CreateProducer(datadescription, gridname):

   "Creates a producer proxy for the grid"
   if not datadescription.GetInputDescriptionByName(gridname):
       raise RuntimeError, "Simulation input name '%s' does not exist" % gridname
   grid = datadescription.GetInputDescriptionByName(gridname).GetGrid()
   producer = TrivialProducer()
   producer.GetClientSideObject().SetOutput(grid)
   if grid.IsA("vtkImageData") == True or grid.IsA("vtkStructuredGrid") == True or grid.IsA("vtkRectilinearGrid") == True:
       extent = datadescription.GetInputDescriptionByName(gridname).GetWholeExtent()
       producer.WholeExtent= [ extent[0], extent[1], extent[2], extent[3], extent[4], extent[5] ]
   producer.UpdatePipeline()
   return producer

def CreateWriter(proxy_ctor, filename, freq, cp_writers):

   writer = proxy_ctor()
   writer.FileName = filename
   writer.add_attribute("cpFrequency", freq)
   writer.add_attribute("cpFileName", filename)
   cp_writers.append(writer)
   return writer

</source>