Coprocessing example

From KitwarePublic
Revision as of 15:28, 13 September 2011 by Andy.bauer (Talk | contribs)

Jump to: navigation, search

This example is used to demonstrate how the co-processing library can be used with a C++ based simulation code. In the ParaView/CoProcessing/Adaptors/FortranAdaptors directory there is code useful for integrating C or Fortran based simulation codes withthe co-processing library. Note that this example requires MPI to be available on your system. The executable takes in a python coprocessing script and a number of time steps to be run for. Note to remember to set your system environment properly. See [[1]] for details.


CoProcessingExample.cxx

#include "vtkCPDataDescription.h"
#include "vtkCPInputDataDescription.h"
#include "vtkCPProcessor.h"
#include "vtkCPPythonScriptPipeline.h"
#include "vtkElevationFilter.h"
#include "vtkPolyData.h"
#include "vtkSmartPointer.h"
#include "vtkSphereSource.h"
#include "vtkXMLUnstructuredGridReader.h"
 
#include <mpi.h>
#include <string>
 
class DataGenerator {
public:
  DataGenerator()
    {
    this->Sphere = vtkSmartPointer<vtkSphereSource>::New();
    this->Sphere->SetThetaResolution(30);
    this->Sphere->SetPhiResolution(30);
    int procId;
    MPI_Comm_rank(MPI_COMM_WORLD, &procId);
    this->Sphere->SetCenter(procId*4.0, 0, 0);
    this->Elevation = vtkSmartPointer<vtkElevationFilter>::New();
    this->Elevation->SetInputConnection(this->Sphere->GetOutputPort());
    this->Index = 0;
    }
 
  vtkSmartPointer<vtkPolyData> GetNext()
    {
    double radius = fabs(sin(0.1 * this->Index));
    this->Index++;
    this->Sphere->SetRadius(1.0 + radius);
    this->Elevation->Update();
    vtkSmartPointer<vtkPolyData> ret = vtkSmartPointer<vtkPolyData>::New();
    ret->DeepCopy(this->Elevation->GetOutput());
    return ret;
    }
 
protected:
  int Index;
  vtkSmartPointer<vtkSphereSource> Sphere;
  vtkSmartPointer<vtkElevationFilter> Elevation;
};
 
int main(int argc, char* argv[])
{
  if (argc < 3)
    {
    printf("Usage: %s <cp python file> <number of steps>\n", argv[0]);
    return 1;
    }
  // we assume that this is done in parallel
  MPI_Init(&argc, &argv);
 
  std::string cpPythonFile = argv[1];
  int nSteps = atoi(argv[2]);
 
  vtkCPProcessor* processor = vtkCPProcessor::New();
  processor->Initialize();
  vtkCPPythonScriptPipeline* pipeline = vtkCPPythonScriptPipeline::New();
 
  // read the coprocessing python file
  if(pipeline->Initialize(cpPythonFile.c_str()) == 0)
    {
    cout << "Problem reading the python script.\n";
    return 1;
    }
 
  processor->AddPipeline(pipeline);
  pipeline->Delete();
 
  if (nSteps == 0)
    {
    return 0;
    }
 
  // create a data source, typically this will come from the adaptor
  // but here we use generator to create it ourselves
  DataGenerator generator;
 
  // do coprocessing
  double tStart = 0.0;
  double tEnd = 1.0;
  double stepSize = (tEnd - tStart)/nSteps;
 
  vtkCPDataDescription* dataDesc = vtkCPDataDescription::New();
  dataDesc->AddInput("input");
 
  for (int i = 0; i < nSteps; ++i)
    {
    double currentTime = tStart + stepSize*i;
    // set the current time and time step
    dataDesc->SetTimeData(currentTime, i);
 
    // check if the script says we should do coprocessing now
    if(processor->RequestDataDescription(dataDesc) != 0)
      {
      // we are going to do coprocessing so use generator to
      // create our grid at this timestep and provide it to
      // the coprocessing library
      vtkSmartPointer<vtkDataObject> dataObject =
        generator.GetNext();
 
      dataDesc->GetInputDescriptionByName("input")->SetGrid(dataObject);
      processor->CoProcess(dataDesc);
      }
    }
 
  dataDesc->Delete();
  processor->Finalize();
  processor->Delete();
 
  MPI_Finalize();
 
  return 0;
}

CMakeLists.txt

cmake_minimum_required(VERSION 2.6)
 
PROJECT(CoProcessingExample)
 
FIND_PACKAGE(ParaView REQUIRED)
INCLUDE(${PARAVIEW_USE_FILE})
 
INCLUDE_DIRECTORIES(${VTK_MPI_INCLUDE_DIR})
IF(NOT VTK_MPI_INCLUDE_DIR)
  MESSAGE(SEND_ERROR "Must build ParaView with MPI.")
ENDIF()
 
ADD_EXECUTABLE(CoProcessingExample CoProcessingExample.cxx)
TARGET_LINK_LIBRARIES(CoProcessingExample vtkCoProcessor)

Python Scripts

The first python script below is used to just output the actual results of the example. This would correspond to a simulation run with a coarse grid in order to set up coprocessing runs for larger grids where outputting the entire simulation results can be computationally prohibitive.

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()
 
    Sphere1 = CreateProducer( datadescription, "input" )
    ParallelPolyDataWriter1 = CreateWriter( XMLPPolyDataWriter, "input_grid_%t.pvtp", 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)
    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

This second script is still rather simple and only performs a cut on the input from the simulation code. It demonstrates though how desired results can be obtained while performing coprocessing at specified time steps.

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 % 5 == 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()
 
    filename_0_pvtp = CreateProducer( datadescription, "input" )
 
    Clip2 = Clip( guiName="Clip2", InsideOut=0, UseValueAsOffset=0, Scalars=['POINTS', 'Elevation'], Value=0.0, ClipType="Plane" )
    Clip2.ClipType.Normal = [0.0, 1.0, 0.0]
    Clip2.ClipType.Origin = [1.9999999105930328, 0.0, 0.0]
    Clip2.ClipType.Offset = 0.0
 
    ParallelUnstructuredGridWriter2 = CreateWriter( XMLPUnstructuredGridWriter, "Cut_%t.pvtu", 5, 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)
  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