[vtkusers] "manually" specify extents for vtkXMLPRectilinearGridWriter

[Crni Gorac]

About 3 months ago, I posted here, asking on how to properly use
vtkXMLPRectilinearGridWriter:
    http://www.vtk.org/pipermail/vtkusers/2012-July/124770.html

Basically, I'm still facing the same problem.  Attached is short sample
program, that will create pieces of rectilinear grid, and then try to
write them in parallel if multiple instances of this process launched,
for example through "mpirun -np 2 ./foo".  Individual .vtr files seems
to be written correctly, however extents are messed up within master
.pvtr file, so for example Paraview is unable to properly display the
whole grid.

I'd appreciate any advice, in particular pointers to source code
samples, on how to properly utilize vtkXMLPRectilinearGridWriter to
write pieces of the grid into multiple files in parallel.  Or, on an
alternative approach on how to save pieces of grid from multiple MPI
processes, so that these could be later loaded into ParaView as single
entity.

Thanks.

PS. Side question, as my dataset is actually structured grid: How to
define "physical" size of the box for structured grid, so that I don't
have to resort to using rectilinear grid for its Set[XYZ]Coordinates()
methods?
-------------- next part --------------
#include <mpi.h>
#include <vtkFloatArray.h>
#include <vtkPointData.h>
#include <vtkRectilinearGrid.h>
#include <vtkXMLPRectilinearGridWriter.h>

#define X0 0
#define Y0 0
#define Z0 0
#define X1 1
#define Y1 1
#define Z1 1
#define NX 100
#define NY 100
#define NZ 100
#define FILENAME "foo.pvtr"

static float f(const float x, const float y, const float z)
{
    return x + y * y + z * z * z;
}

int main(int argc, char **argv)
{
    MPI_Init(&argc, &argv);

    int size;
    int rank;
    MPI_Comm_size(MPI_COMM_WORLD, &size);
    MPI_Comm_rank(MPI_COMM_WORLD, &rank);

    int z_lo = NZ * rank / size;
    int nz = NZ * (rank + 1) / size - z_lo;
    float dx = (float) (X1 - X0) / NX;
    float dy = (float) (Y1 - Y0) / NY;
    float dz = (float) (Z1 - Z0) / NZ;

    vtkFloatArray *x = vtkFloatArray::New();
    x->SetNumberOfValues(NY);
    for (int i = 0; i < NX; ++i)
        x->SetValue(i, X0 + dx / 2 + i * dx);

    vtkFloatArray *y = vtkFloatArray::New();
    y->SetNumberOfValues(NY);
    for (int i = 0; i < NY; ++i)
        y->SetValue(i, Y0 + dy / 2 + i * dy);

    vtkFloatArray *z = vtkFloatArray::New();
    z->SetNumberOfValues(nz);
    for (int i = 0; i < nz; ++i)
        z->SetValue(i, Z0 + z_lo * dz + dz / 2 + i * dz);
 
    vtkFloatArray *p = vtkFloatArray::New();
    p->SetNumberOfValues(nz * NY * NX);
    for (int i = 0; i < nz; ++i)
        for (int j = 0; j < NY; ++j)
            for (int k = 0; k < NX; ++k)
                p->SetValue(i * NY * NX + j * NX + k, f(x->GetValue(k), y->GetValue(j), z->GetValue(i)));

    vtkRectilinearGrid *grid = vtkRectilinearGrid::New();
    grid->SetDimensions(NX, NY, nz);
    grid->SetXCoordinates(x);
    grid->SetYCoordinates(y);
    grid->SetZCoordinates(z);
    grid->GetPointData()->SetScalars(p);

    vtkXMLPRectilinearGridWriter *writer = vtkXMLPRectilinearGridWriter::New();
    writer->SetNumberOfPieces(size);
    writer->SetStartPiece(rank);
    writer->SetEndPiece(rank);
    writer->SetInput(grid);
    writer->SetFileName(FILENAME);
    writer->Write();

    x->Delete();
    y->Delete();
    z->Delete();
    p->Delete();
    grid->Delete();
    writer->Delete();

    MPI_Finalize();
}

// Local Variables:
// c-basic-offset: 4
// End: